      SUBROUTINE READ(DLAT,DLONG,DATE,HT,LQ,JK,*)
C     SUBROUTINE READ.  READ TITLE FOR JOB. FOR EACH IONOGRAM, READ ID
C     (NO. OF POINTS,  YR,DAY #,GMT OF DATA START); PROCESSING SWITCHES:
C        MAP=1 IF SATELLITE POSITION TO BE INTERPOLATED FROM INPUT WORLD
C             MAP DATA; MAP=0 IF THIS DATA IS INPUT DIRECTLY.
C        ISAT = SATELLITE ID: 1=AL-1; 2=AL-2; 3=ISIS-1; 4=ISIS-2.
C             IF ZERO OR BLANK, VARIABLE POSITION ANALYSIS NOT PERFORMED
C        ICODE IDENTIFIES CONDITIONS APPLIED IN EACH ANALYSIS (SEE MAIN)
C        IORD=1 IF O-TRACE VIRTUAL HEIGHTS ARE TO BE COMPUTED FROM
C             X-TRACE N(H) PROFILE. IF ZERO, OMIT THIS CALCULATION.
C             IP=1 IF PLOTTING DESIRED. IP=0 OR BLANK OTHERWISE.
C        METH IDENTIFIES TYPE OF LAMINATION USED IN ANALYSIS (1=LINEAR
C             IN LOG(N); 2(OR BLANK)=PARABOLIC IN LOG(N)).
C     AND SCALED H'(F) DATA.
C     Added JK argument to make compiler happy with ENTRY below.  Original
C     code was dubious.
      IMPLICIT REAL*8(A-H,O-Z)
      REAL*8 APPF(8)/2.,2.165,2.39,2.75,3.,5.,5.52,6./, REALF(8)/2.,2.1,
     R  2.3,2.7,3.,5.,5.5,6./, RAD/57.29577951D0/
      CHARACTER*8 SATID(5)/'AL-1:   ','AL-2:   ','ISIS-1: '
     R  ,'ISIS-2: ','        '/
      REAL*4 DLAT,DLONG,DATE,HT,BN,BE,BV,BT
      LOGICAL*4 TEST
      LOGICAL*1 STRING(80)
      COMMON /CALC/ AA(50),BB(50),F(50),HP(50),FH(50),TH(70),DEN(70)
      COMMON /MISC/ DP,FH1,ASUM,BSUM,FA,CPH,X1,X2,GMT1
      COMMON /TIME/IH,IM,IS
      COMMON /IO/ IYR,IDY,ITI,LMT,TITLE(40),DIP,FHS,MON,IDAY
      COMMON /SWJTCH/ ISAT,LOCSAT,IFXS,NCASES,ICODE(4),N,METHOD,IORD,IP
      DIMENSION ITLO(3),DDEC(27),DET(27),DNO(27),DD(27,4),EQM(27,4)
      SAVE
      DATA ISKP/0/, K/0/
C     OPEN(UNIT=3,NAME='NHINPUT.DAT',TYPE='OLD',FORM='FORMATTED')
      IF (ISKP .NE. 0) GO TO 100
C     LIST ALL USER INPUT DATA CARDS.
      I1=0
      I=0
      J=1
      WRITE(6,22)J
 1400 READ(3,23,END=1300) STRING
      I1=I1+1
      IF(I1.LT.29) GO TO 1400
      I=I+1
      IF(MOD(I,57).NE.0) GO TO 1500
      J=J+1
      WRITE(6,22)J
 1500 WRITE(6,24) STRING
      GO TO 1400
 1300 REWIND 3
C     SET UP AND TEST FIELDG. INITIALIZE PARAMETERS FOR THIS JOB.
      LQ = 0
      AA(1)=0.0
      BB(1)=0.0
      BB(2)=0.0
      BB(3)=0.
      CALL FIELDG(-60.,-180.,100.,1960.,10,LQ,BN,BE,BV,BT)
      WRITE (6,11) BN,BE,BV,BT
C     INPUT: DECLINATION FROM NAUTICAL ALMANAC, EQUATION OF TIME FROM EP
C     DATA AT NOON GREENWICH, IN DECIMAL UNITS, FOR 4 SUCCESSIVE YEARS.
      DO 25 J=1,4
  25  READ(3,18)  (DNO(I),DD(I,J),EQM(I,J),I=1,27)
      DO 20 I=1,27
   20 DNO(I)=DNO(I)+0.5
      RADI=1.0D0/RAD
      STM=1.0D0/60.0D0
      HTD=1.0D0/24.0D0
      CALL INREF
      READ (3,1) TITLE
      ISKP=1
  100 READ(3,2,END=200) N,IYR,IDY,IH,IM,IS,MAP,ISAT,ICODE,IORD,IP,METH
      DATE=FLOAT(IYR)+FLOAT(IDY)/365.0+1900.0
      ITI=IH*10000 + IM*100 + IS
      GMT1=DFLOAT(IH)+DFLOAT(IM)/60.0+DFLOAT(IS)/3600.0
      IF (MAP .LT. 1) GO TO 300
      TIME=GMT1
      CALL ORBINT(0,TIME,DLAT,DLONG,HT)
      IF(MAP.EQ.2) READ(3,3) FHS
      GO TO 400
  300 READ(3,3) HT,DIP,FHS,DLAT,DLONG
  400 READ(3,3) (HP(I),F(I),I=1,N)
      IF (ISAT .NE. 3) GO TO 800
C     CORRECT FREQUENCIES SCALED FOR ISIS-1 IONOGRAMS.
      DO 900 I=1,N
      FREQ=F(I)
      IF((FREQ.GT.2..AND.FREQ.LT.3.) .OR. (FREQ.GT.5..AND.FREQ.LT.6.))
     I  F(I)=TAB(APPF,REALF,8,1,FREQ)
  900 CONTINUE
C     CALCULATE LOCAL MEAN TIME
  800 TLO=GMT1+DLONG/15.D0+24.D0
      DO 700 I=1,3
      ITLO(I)=TLO
  700 TLO=(TLO-ITLO(I))*60.0
      LMT=ITLO(1)*1.E4+ITLO(2)*1.0E2+ITLO(3)
      IF (LMT .GT. 240000) LMT=LMT-240000
C     COMPUTE MONTH AND DAY FROM DAY OF YEAR. COMPUTE SOLAR ZENITH ANGLE
C     USING DECLINATION AND EQN OF TIME FOR CORRESPONDING YEAR IN PERIOD
C     1968 TO 1971.
      CALL DAYMON(IYR,IDY,MON,IDAY)
      I=K
      K=MOD(IYR,4)+1
      IF(K.EQ.I) GO TO 1200
      DO 1100 I=1,27
      DDEC(I)=DD(I,K)
 1100 DET(I)=EQM(I,K)
 1200 CPH=DLAT*RADI
      DDATE=GMT1*HTD+DFLOAT(IDY)
      ET=TAB(DNO,DET,27,4,DDATE)*STM
      HA=(DABS(180.0D0-DLONG-15.0D0*(GMT1-ET)))*RADI
      DEC=TAB(DNO,DDEC,27,4,DDATE)*RADI
      COSZ=DSIN(DEC)*DSIN(CPH)+DCOS(DEC)*DCOS(CPH)*DCOS(HA)
      ZENI=DACOS(COSZ)*RAD
      HA=HA*RAD
      DEC=DEC*RAD
      ET=ET*60.0D0
C     COMPUTE DIP AND FH AT ORIGINAL SATELLITE POSITION
      CALL FIELDG(DLAT,DLONG,HT,DATE,10,LQ,BN,BE,BV,BT)
      DIP=ATAN2(BV,SQRT(BT**2-BV**2))*RAD
C     CPH=COS(PROPAGATION ANGLE). VERTICAL PROPAGATION IS ASSUMED.
      CPH=DCOS((90.0-DABS(DIP))/RAD)
      FH(1)=BT*2.8000E-05
C     EXT. RAY EXIT FREQ. MUST BE GREATER THAN FHS
      FXS=FH(1)+0.001
      IF (F(1).GE.FXS) GO TO 710
      F(1)=FXS
      WRITE (6,17)
  710 CONTINUE
      FH1=FH(1)
      IF(MAP.EQ.2) F(1)=F(1)-FHS+FH1
      FHS = FH(1)
      TH(1) = HT
C     DETERMINE NCASES, THE NUMBER OF WAYS IN WHICH THIS IONOGRAM IS TO
C     BE PROCESSED.
      IF (ICODE(1) .EQ. 0) GO TO 1000
C     EXCLUDE VARIABLE SATELLITE ANALYSIS IF ORBIT DATA OR SATELLITE ID
C     NOT SPECIFIED, AND CALL FOR FIXED POSITION ANALYSIS IF NOT
C     PREVIOUSLY PERFORMED.
      TEST=MAP.LT.1 .OR. ISAT.EQ.0
      DO 600 I=1,4
      J=I
      IF (TEST .AND. ICODE(I).GT.2) ICODE(I)=ICODE(I)-2
      IF (ICODE(I) .EQ. 0 .OR. (TEST .AND. J.GT.2)) GO TO 500
  600 CONTINUE
      J=5
      GO TO 500
 1000 ICODE(1)=1
      J=2
  500 NCASES=J-1
      J=ICODE(1)
C     INITIALIZE METHOD FOR TYPE OF LAMINATION DESIRED. WRITE OUT INPUT.
      ENTRY PRINT(JK)
C     J=JK
C     See if ENTRY was called or we dropped down from subroutine (jk=-1).
      if( jk .ge. 0 ) j = jk
      METHOD = 2
      IF(METH.EQ.1)METHOD=1
      WRITE (6,4) TITLE
      IF(METHOD.EQ.1)WRITE(6,19)
      IF(METHOD.EQ.2)WRITE(6,21)
      WRITE(6,5) IDY,MON,IDAY,IYR,ITI,LMT,DEC,ET,ZENI
      WRITE (6,6) ITI,HT,DLAT,DLONG,DIP,FHS
      I2=ISAT
      IF(ISAT.EQ.0) I2=5
      WRITE(6,16) SATID(I2)
      IF (J .EQ. 1) WRITE (6,12)
      IF (J .EQ. 2) WRITE (6,13)
      IF (J .EQ. 3) WRITE (6,14)
      IF (J .EQ. 4) WRITE (6,15)
      IF (ISAT .EQ. 3) WRITE (6,9)
      WRITE (6,7)
      WRITE (6,8) (F(I),HP(I),I=1,N)
      RETURN
  200 RETURN 1
    1 FORMAT (20A4)
    2 FORMAT(3I10,4X,3I2,2I10,6X,4I1,I5,I2,I3)
    3 FORMAT (8F10.4)
    4 FORMAT(1x,20A4/1H ,20A4/' GEOMAGNETIC FIELD MODEL GSFC 9/65, NMAX 
     F=10 (HENDRICKS AND CAIN, 1966).')
    5 FORMAT(' INPUT DATA'/'  DAY',I4,' (',I2,'/',I2,') OF 19',I2,' AT',
     F  I7,' GMT  LOCAL MEAN TIME',I7/' DECLINATION',F7.2,' DEG  EQN OF 
     FTIME',F7.2,' MIN  ZENITH ANGLE',F8.2,' DEG')
    6 FORMAT(' SATELLITE POSITION AT ',I6,' GMT'/'  HT=',F6.1,'  LAT=',
     F       F7.2,'  LONG=',F7.2,'  DIP=',F7.2,'  FHS=',F5.3)
    7 FORMAT(1x,'X TRACE'/1H ,5(3X,'F',6X,'HP',2X))
    8 FORMAT(5(1X,F6.3,1X,F6.1))
    9 FORMAT(1H ,8X,'CORRECTED SCALED FREQUENCIES FOLLOW.')
   11 FORMAT(1H ,'SEPTEMBER 65 COEFFICIENTS'/1H ,2X,'LAT',5X,'LONG',5X,
     F       'HT',4X,'DATE',6X,'BN',11X,'BE',11X,'BV',11X,'BT'/1H ,
     F      '-60.00',2X,'-180.00',2X,'100.0',2X,'1960.0',4F13.4)
   12 FORMAT(1H ,8X,'FIXED SATELLITE POSITION AND KNOWN FXS ASSUMED.')
   13 FORMAT(1H ,8X,'FIXED SATELLITE POSITION AND UNKNOWN FXS ASSUMED.')
   14 FORMAT(1H ,8X,'VARIABLE SATELLITE POSITION AND KNOWN FXS ASSUMED.'
     F       )
   15 FORMAT(1H ,8X,'VARIABLE SATELLITE POSITION AND UNKNOWN FXS ASSUMED
     F.')
   16 FORMAT(' ',A8)
   17 FORMAT(1H ,'FREQUENCY READ TOO LOW ON FIRST DATA POINT.')
   18 FORMAT(4(F4.0,2F8.5))
   19 FORMAT(' LAMINATIONS LINEAR IN LOG(N)')
   21 FORMAT(' LAMINATIONS PARABOLIC IN LOG(N)')
   22 FORMAT(//'      ***** LISTING OF ALL INPUT DATA CARDS: PAGE',I3,
     F  ' *****'//)
   23 FORMAT(80A1)
   24 FORMAT(' ',80A1)
      END
