C =============================================================================
C PBO driver for hab_los.f
C
C                             Lyman beta         Lyman gamma
C   LINE_LABEL                    2                  3
C   WAVELN                    1025.72             972.537
C   FLUORB                       0.882              0.839
C   FLUORESCENCE WAVELN       6562.725 (3)       4861.287 (3)
C      & WEIGHTS              6562.771 (1)       4861.298 (1)
C
C required input definitions:
C   FLUX_C   . . . . solar "line center" flux
C   LINE_LABEL . . . index identifying parent resonance line
C   WAVELN . . . . . parent resonance wavelength
C   FLUORB . . . . . branching ratio for fluorescence of interest
C   MAXLOS . . . . . maximum array size
C   ILOS . . . . . . number of lines of sight (ILOS </= MAXLOS)
C   SZA  . . . . . . array of observer solar zenith angles
C   ZNTH . . . . . . array of observation zenith angles
C                      (0 <\= ZNTH < 90 deg)
C   AZI  . . . . . . array of observation azimuth angles with respect to
C                      solar plane (interior angles:  0 </= AZI </= 180 deg)
C   LPROF  . . . . . line profile precision
C
C main output definitions:
C   ACER . . . . . . array of LOS radiances (R) for unit solar flux
C                      ACER(1,#) single scattering component
C                      ACER(2,#) multiple scattering component
C                      ACER(3,#) sum of components
C   PRFL . . . . . . corresponding line profiles of length LPROF for each ACER,
C                      giving one wing of emission profile for each component,
C                      normalized to total intensity at each line center
C   SDPNT  . . . . . LPROF "speed" quadrature point array,
C   WSPD   . . . . . LPROF "speed" quadrature weight array
C
C secondary output definitions:
C   LABEL_OUT  . . . name of fluorescence feature
C   WVLN_OUT . . . . wavelength of fluorescence feature (needed in translating
C                    spectral displacement to wavelength)
C   DOPPLER_SPEED  . Doppler speed at reference exobase temperature
C   EXOBASE_T  . . . reference exobase temperature
C   HSHAD  . . . . . array of LOS shadow heights (km) (useful for plotting)
C
C =============================================================================

C PRELIMINARY
      IMPLICIT  REAL*8 (A-H,O-Z)
      PARAMETER (MAXLOS = 198)                       !! increase as needed
      PARAMETER (LPROF  =  16)                       !! OPTIONS: 4 8 16
      DIMENSION SZA(MAXLOS), ZNTH(MAXLOS), AZI(MAXLOS)
      DIMENSION HSHAD(MAXLOS), ACER(3,MAXLOS), PRFL(3,0:LPROF,MAXLOS)
      DIMENSION SPDPT(LPROF), WSPD(LPROF)
      DIMENSION SHAD(maxlos), ADUM(9)

      CHARACTER*20 LABEL_OUT

CC --- i added this here -em
C COLLECT BASIC PARAMETERS
      OPEN(01,NAME='inputs_los.dat',STATUS='OLD')
      READ(01,*) ILOS, LINE_LABEL, WAVELN, FLUORB, FLUX_C
      WRITE(06,*) ILOS, LINE_LABEL, WAVELN, FLUORB, FLUX_C
      FLUORB = 1 - FLUORB
CC      READ(01,*) ILOS
      DO 101 LK = 1,ILOS
        READ(01,*) sequence, inttime, timeobs, targaz, targzd, sunaz, 
     &  sunzd, diffaz, shadalt, zdbin, gallat, gallon, 
     &  vlsr, igeo, dopw, error_igeo, error_dopw
CC
C        WRITE(06,*) sequence
C        WRITE(06,*) inttime
C        WRITE(06,*) timeobs
C        WRITE(06,*) targaz
C        WRITE(06,*) sunaz
C        WRITE(06,*) zdbin
C        WRITE(06,*) gallat
C        WRITE(06,*) gallon
C        WRITE(06,*) vlsr
C        WRITE(06,*) igeo
C        WRITE(06,*) dopw
C        WRITE(06,*) error_igeo
C        WRITE(06,*) error_dopw
CC
        SZA(LK)  = sunzd
        ZNTH(LK) = targzd
        AZI(LK)  = diffaz
        SHAD(LK) = shadalt

        WRITE(06,*) 'sequence', sequence
        WRITE(06,*) 'inttime', inttime
C        WRITE(06,*) 'sza', sunzd
C        WRITE(06,*) SZA(LK)
C        WRITE(06,*) 'znth', targzd
C        WRITE(06,*) ZNTH(LK)
C        WRITE(06,*) 'azi', diffaz
C        WRITE(06,*) AZI(LK)
C        WRITE(06,*) 'shad', shadalt
C        WRITE(06,*) SHAD(LK)
 
C        READ(01,*) ADUM
C        READ(01,*) asda, atargzd, adiffaz, ashadalt
C        SZA(LK)  = asda + 90.0
C        ZNTH(LK) = atargzd
C        AZI(LK)  = adiffaz
C        SHAD(LK) = ashadalt
C        READ(01,*) ADUM
 101  CONTINUE
      CLOSE(01)

CCC COLLECT BASIC PARAMETERS
CC      OPEN(01, NAME='inputs_los.dat', STATUS='OLD')
CC      READ(01,*) LINE_LABEL, WAVELN, FLUORB, FLUX_C
CC      READ(01,*) ILOS
CC      DO LK = 1,ILOS
CC        READ(01,*) asda, atargzd, adiffaz
CC        SZA(LK)  = asda + 90.0
CC        ZNTH(LK) = atargzd
CC        AZI(LK)  = adiffaz
CC      END DO
CC      CLOSE(01)

C open am inputs here to get amd_exo... -em, 11-29-00
      OPEN(02,NAME='saveam.rt',STATUS='OLD')
      READ(02,*) AMLINE_LABEL
      READ(02,*) AMGLAT,AMGLONG
      READ(02,*) AMIDAY, AMIYEAR
      READ(02,*) AMUT
      WRITE(06,*) 'AMUT', AMUT
      READ(02,*) AMIAHP,AMAP1,AMAP2,AMAP3,AMAP4,AMAP5,AMAP6,AMAP7
      READ(02,*) AMF107,AMF107A
      READ(02,*) AMMOD_MSIS,AMD_EXO,AMFLUX,AMD_MAX
      READ(02,*) AMIGEOM,AMSATT,AMSATD
      CLOSE(02)

C open pm inputs here to get amd_exo... -em, 11-29-00
      OPEN(03,NAME='savepm.rt',STATUS='OLD')
      READ(03,*) PMLINE_LABEL
      READ(03,*) PMGLAT,PMGLONG
      READ(03,*) AMIDAY, AMIYEAR
      READ(03,*) PMUT
      WRITE(06,*) 'PMUT', PMUT
      READ(03,*) PMIAHP,PMAP1,PMAP2,PMAP3,PMAP4,PMAP5,PMAP6,PMAP7
      READ(03,*) PMF107,PMF107A
      READ(03,*) PMMOD_MSIS,PMD_EXO,PMFLUX,PMD_MAX
      READ(03,*) PMIGEOM,PMSATT,PMSATD
      CLOSE(03)
 
      WRITE(06,*) ILOS, LINE_LABEL, WAVELN, FLUORB, FLUX_C

      IF (LINE_LABEL.EQ.2) THEN
        LABEL_OUT = 'H ALPHA'
        WVLN_OUT  = 6562.7
      ELSE IF (LINE_LABEL.EQ.3) THEN
        LABEL_OUT = 'H BETA '
        WVLN_OUT  = 4861.3
      ELSE
        WRITE(06,'(A)') 'INPUT LINE_LABEL NOT DEFINED'
        GO TO 999
      END IF

C -----------------------------------------------------------------------------
C the big call

      CALL HAB_LOS(LINE_LABEL,MAXLOS,LPROF,WAVELN,FLUORB,
     &             SPDPT,WSPD,ILOS,SZA,ZNTH,AZI,HSHAD,ACER,PRFL,
     &             EXOBASE_T,DOPPLER_SPEED)

C -----------------------------------------------------------------------------
C OUTPUT
C example:  here only the total apparent column emission rates (ACER)
C           and corresponding line profiles are output

C OPEN OUTPUT FILE
      OPEN(10, FILE='hab_los.dat', STATUS='UNKNOWN')

C 10   FORMAT(1X,2I4,2X,1P1E10.2,0P1F9.2,1P1E10.2)
 10   FORMAT(1X,2I4,2X,1P1E10.2,0P1F9.2,1P1E10.2,1P1E10.2,1P1E10.2)
 11   FORMAT(1X,0P6F11.3)
 12   FORMAT(1X,1P6E11.3)

      WRITE(10,'(A)') LABEL_OUT
      WRITE(10,11) WAVELN, FLUORB, WVLN_OUT
C      WRITE(10,10) ILOS, LPROF, FLUX_C, EXOBASE_T, DOPPLER_SPEED
      WRITE(10,10)ILOS,LPROF,FLUX_C,EXOBASE_T,DOPPLER_SPEED,AMD_EXO,
     &PMD_EXO
      WRITE(10,12) SPDPT

      DO II = 1,ILOS
        ANORM = PRFL(3,1,II)
        DO LP = 1,LPROF
          PRFL(1,LP,II) = PRFL(1,LP,II) / ANORM
          PRFL(2,LP,II) = PRFL(2,LP,II) / ANORM
          PRFL(3,LP,II) = PRFL(3,LP,II) / ANORM
        END DO
      END DO

      DO II = 1,ILOS
        BLMR = FLUX_C * ACER(3,II)
        WRITE(10,11) SZA(II), ZNTH(II), AZI(II), HSHAD(II), BLMR
        WRITE(10,12) (PRFL(1,LP,II), LP=1,LPROF)
        WRITE(10,12) (PRFL(2,LP,II), LP=1,LPROF)
        WRITE(10,12) (PRFL(3,LP,II), LP=1,LPROF)
      END DO

      CLOSE(10)

C -----------------------------------------------------------------------------

 999  CONTINUE

      STOP
      END

C =============================================================================