PROGRAM REZEST_FDC
      IMPLICIT NONE
      INCLUDE "lun.par"
      CHARACTER*3 DET
      REAL*4 B, EPSILON, N_RL, N_RL_FRONT, LAMBDA, L, N_M, N_M_POLAR, M,
     X     P, R_M, STEREO
      REAL*4 DK_MS, DK_RES, DK
      REAL*4 DPHI_MS, DPHI_RES, DPHI_K, DPHI
      REAL*4 DTHETA_MS, DTHETA_RES, DTHETA
      REAL*4 PT, K, DP_OVER_P

      B = 1.74
      EPSILON = 0.00015
      STEREO = 6.0*3.14159/180.0
      N_RL = 0.0317
      N_RL_FRONT = 0.0212
      LAMBDA = 40*3.14159/180.0
      L = 0.58600 - 0.16895
      R_M = (0.58600 + 0.16895)/2.0
      N_M = 23
      N_M_POLAR = 10
      DET = 'CDC'

      M = 0.139

      open(UNIT=LUNOUT, FILE='cdc_pi.txt', STATUS='UNKNOWN')
      open(UNIT=LUNOUTA, FILE='cdc_pi_angle.txt', STATUS='UNKNOWN')
      open(UNIT=LUNOUTB, FILE='cdc_pi_polar.txt', STATUS='UNKNOWN')
      P = 0.01
      DO WHILE (P < 2.0)
         CALL REZEST_CURVE(DET, LAMBDA, L, EPSILON, N_M, P, M, N_RL, B,
     X        DK_MS, DK_RES, DK)
         CALL REZEST_AZIMUTH(DET, LAMBDA, L, EPSILON, N_M, P, M,
     X        N_RL, N_RL_FRONT, B, R_M, DK, DPHI_MS, DPHI_RES, DPHI_K,
     X        DPHI)
         CALL REZEST_POLAR(DET, LAMBDA, L, EPSILON, STEREO,
     X        N_M_POLAR, P, M, N_RL, N_RL_FRONT, DTHETA_MS,
     X        DTHETA_RES, DTHETA)
         PT = P*COS(LAMBDA)
         K = 0.3*B/PT
         DP_OVER_P = DK/K
         WRITE (LUNOUT,*) "p ", P, " dp/p ", DP_OVER_P
         WRITE (LUNOUTA,*) "p ", P, " dphi ", DPHI
         WRITE (LUNOUTB,*) "p ", P, " dtheta ", DTHETA
         P = P*1.02
      END DO
      CLOSE(LUNOUT)
      CLOSE(LUNOUTA)
      CLOSE(LUNOUTB)

      M = 0.938

      OPEN(UNIT=LUNOUT, FILE='cdc_p.txt', STATUS='UNKNOWN')
      OPEN(UNIT=LUNOUTA, FILE='cdc_p_angle.txt', STATUS='UNKNOWN')
      OPEN(UNIT=LUNOUTB, FILE='cdc_p_polar.txt', STATUS='UNKNOWN')
      P = 0.01
      DO WHILE (p < 2.0)
         CALL REZEST_CURVE(DET, LAMBDA, L, EPSILON, N_M, P, M, N_RL, B,
     X        DK_MS, DK_RES, DK)
         CALL REZEST_AZIMUTH(DET, LAMBDA, L, EPSILON, N_M, P, M,
     X        N_RL, N_RL_FRONT, B, R_M, DK, DPHI_MS, DPHI_RES, DPHI_K,
     X        DPHI)
         CALL REZEST_POLAR(DET, LAMBDA, L, EPSILON, STEREO,
     X        N_M_POLAR, P, M, N_RL, N_RL_FRONT, DTHETA_MS,
     X        DTHETA_RES, DTHETA)
         PT = P*COS(LAMBDA)
         K = 0.3*B/PT
         DP_OVER_P = DK/K
         WRITE (LUNOUT,*) "p ", P, " dp/p ", DP_OVER_P
         WRITE (LUNOUTA,*) "p ", P, " dphi ", DPHI
         WRITE (LUNOUTB,*) "p ", P, " dtheta ", DTHETA
         P = P*1.02
      END DO
      CLOSE(LUNOUT)
      CLOSE(LUNOUTA)
      CLOSE(LUNOUTB)

      STOP
      END