C path: /home/rc1/aer_lblrtm/src/SCCS/s.xmerge.f
C revision: 5.2
C created: 05/27/98 11:33:24
C presently: 05/27/98 12:01:38
C
C ----------------------------------------------------------------
C
SUBROUTINE XMERGE (NPTS,LFILE,MFILE,JPATHL) 8,11
C
IMPLICIT REAL*8 (V)
C
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
C XMERGE CALL ABSMRG,EMINIT,RADMRG
C
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYER,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /HVERSN/ HVRLBL,HVRCNT,HVRFFT,HVRATM,HVRLOW,HVRNCG,
* HVROPR,HVRPST,HVRPLT,HVRTST,HVRUTL,HVRXMR
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYHDR,YI1,YID(10),LSTWDF
C
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /XME/ V1R4,V2R4,DVR4,NPTR4,BOUND4,R4(4996)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
EQUIVALENCE (FSCDID(1),IHIRAC) , (FSCDID(2),ILBLF4),
* (FSCDID(3),IXSCNT) , (FSCDID(4),IAERSL),
* (FSCDID(5),IEMIT) , (FSCDID(7),IPLOT),
* (FSCDID(8),IPATHL) , (FSCDID(9),JRAD),
* (FSCDID(11),IMRG) , (FSCDID(16),LAYR1),
* (FSCDID(17),NLAYHD)
C
CHARACTER*8 HVRLBL,HVRCNT,HVRFFT,HVRATM,HVRLOW,HVRNCG,HVROPR,
* HVRPLT,HVRPST,HVRTST,HVRUTL,HVRXMR
C
C ASSIGN SCCS VERSION NUMBER TO MODULE
C
HVRXMR = '5.2'
C
IOD = 0
C
IEXFIL = 20
IAFIL = 14
C
C --------------------------------------------------------------
C
C Special case of merging optical depths from multiple files to
C multiple files, for use when calculating radiance derivatives.
C The subroutine call tree is as follows:
C
C LBLSUB called XLAYER when IHIRAC+IATM+IMRG > 0
C XLAYER called XMERGE when IMRG = 10
C XMERGE calls ABSMRG when IMRG = 10
C
IF (IMRG.EQ.10) THEN
CALL ABSMRG
(NPTS,LFILE,MFILE,JPATHL)
CALL ENDFIL (MFILE)
RETURN
ENDIF
C
C --------------------------------------------------------------
C
C WHEN IAERSL EQUALS 2 CALL ADARSL TO ADD ABSORPTION AND SCATTERING
C TO COMMON BLOCKS FOR USE IN A TRANSMITTANCE RUN
C
C --------------------------------------------------------------
C
C IEMIT = 0 => Optical depth calculation only
C
IF (IEMIT.EQ.0) THEN
IF (LAYER.EQ.1) THEN
CALL ABINIT (NPTS,MFILE,JPATHL)
ELSE
WRITE (IPR,900) XID,(YID(M),M=1,2)
CALL ABSMRG (NPTS,LFILE,MFILE,JPATHL)
ENDIF
ELSE
C
C --------------------------------------------------------------
C
C IEMIT > 0 TO REACH THIS STATEMENT
C
WRITE (IPR,900) XID,(YID(M),M=1,2)
IF (LAYER.EQ.1.AND.IAERSL.GE.1) REWIND IEXFIL
IF (IAERSL.GE.1) CALL GETEXT (IEXFIL,LAYER,IEMIT)
C
TBND = 0.
C
C -----------------------------------------------------------
C
C IEMIT = 1 => Radiance and Transmittance calculated
C
IF (IEMIT.EQ.1) THEN
IF (IMRG.NE.36.AND.IMRG.NE.46) THEN
IF (LAYER.EQ.1) THEN
IF (IPATHL.EQ.1) TBND = TMPBND
CALL EMINIT (NPTS,MFILE,JPATHL,TBND)
ELSE
IF (IPATHL.EQ.3.AND.LAYER.EQ.LH2) TBND = TMPBND
CALL RADMRG (NPTS,LFILE,MFILE,JPATHL,TBND)
ENDIF
ELSE
IF (LAYER.EQ.1) THEN
TBND = TMPBND
CALL FLINIT (NPTS,MFILE,JPATHL,TBND)
ELSE
CALL FLUXUP (NPTS,LFILE,MFILE,JPATHL,TBND)
ENDIF
ENDIF
ENDIF
C
C -----------------------------------------------------------
C
C IEMIT = 3 => Radiance, Transmittance and Radiance
C derivative calculated
C
IF (IEMIT.EQ.3) THEN
IF (LAYER.EQ.1) THEN
IF (IPATHL.EQ.1) TBND = TMPBND
CALL EMADL1 (NPTS,MFILE,JPATHL,TBND)
ELSE
IF (IPATHL.EQ.3.AND.LAYER.EQ.LH2) TBND = TMPBND
CALL EMADMG (NPTS,LFILE,MFILE,JPATHL,TBND)
ENDIF
ENDIF
ENDIF
C
C --------------------------------------------------------------
C
RETURN
C
900 FORMAT (///,1X,10A8,2X,2(1X,A8,1X))
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE XMERGI (NPTS,LFILE,MFILE,JPATHL) 3,7
C
IMPLICIT REAL*8 (V)
C
PARAMETER (MXFSC=200, MXLAY=MXFSC+3,MXZMD=3400,
* MXPDIM=MXLAY+MXZMD,IM2=MXPDIM-2,MXMOL=35,MXTRAC=22)
C
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
C XMERGI CALL ABINIT,ABSINT,ABSOUT
C
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYER,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /MSACCT/ IOD,IDIR,ITOP,ISURF,MSPTS,MSPANL(MXLAY),
* MSPNL1(MXLAY),
* MSLAY1,ISFILE,JSFILE,KSFILE,LSFILE,MSFILE,IEFILE,
* JEFILE,KEFILE
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYHDR,YI1,YID(10),LSTWDF
C
COMMON /XMI/ V1R4,V2R4,DVR4,NPTR4,BOUND4,R4(4815)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
C
EQUIVALENCE (FSCDID(1),IHIRAC) , (FSCDID(2),ILBLF4),
* (FSCDID(3),IXSCNT) , (FSCDID(4),IAERSL),
* (FSCDID(5),IEMIT) , (FSCDID(7),IPLOT),
* (FSCDID(8),IPATHL) , (FSCDID(9),JRAD),
* (FSCDID(11),IMRG) , (FSCDID(16),LAYR1),
* (FSCDID(17),NLAYHD)
C
IF (IEMIT.EQ.1) GO TO 10
IF (LAYER.EQ.LH1.AND.IANT.NE.-1) THEN
CALL ABINIT (NPTS,MFILE,JPATHL)
ELSE
C
WRITE (IPR,900) XID,(YID(M),M=1,2)
CALL ABSINT (NPTS,LFILE,MFILE,JPATHL)
ENDIF
C
GO TO 20
C
C IEMIT = 1 TO REACH THIS STATEMENT
C
10 CONTINUE
WRITE (IPR,900) XID,(YID(M),M=1,2)
IF (LAYER.EQ.1.AND.IAERSL.GE.1) REWIND IEXFIL
IF (IAERSL.GE.1) CALL GETEXT (IEXFIL,LAYER,IEMIT)
C
TBND = 0.
C
IF (IMRG.NE.35.AND.IMRG.NE.45) THEN
IF (LAYER.EQ.LH1.AND.IANT.NE.-1) THEN
IF (JPATHL.EQ.1.AND.LAYER.EQ.1) TBND = TMPBND
CALL EMINIT (NPTS,MFILE,JPATHL,TBND)
ELSE
IF (JPATHL.EQ.3.AND.LAYER.EQ.LH2) TBND = TMPBND
CALL RADINT (NPTS,LFILE,MFILE,JPATHL,TBND)
ENDIF
ELSE
IF (LAYER.EQ.LH1.AND.IANT.NE.-1) THEN
TBND = TMPBND
CALL FLINIT (NPTS,MFILE,JPATHL,TBND)
ELSE
CALL FLUXDN (NPTS,LFILE,MFILE,JPATHL,TBND)
ENDIF
ENDIF
C
20 CONTINUE
C
RETURN
C
900 FORMAT (///,1X,10A8,2X,2(1X,A8,1X))
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE ABINIT (NPTS,MFILE,JPATHL) 2,7
C
IMPLICIT REAL*8 (V)
C
COMMON ODLAY(-2:2407)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYDUM,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /ABSHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /ABSPNL/ V1P,V2P,DVP,NLIM,NSHFT,NPNTS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
DIMENSION XFILHD(2),PNLHDR(2)
DIMENSION ODLAYR(2)
C
CHARACTER*40 CEXT,CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P)
EQUIVALENCE (ODLAY(1),ODLAYR(1)) , (FSCDID(4),IAERSL),
* (FSCDID(5),IEMIT) , (FSCDID(7),IPLOT),
* (FSCDID(8),IPATHL) , (FSCDID(16),LAYR1)
C
C
C ***********************************************************
C ****** THIS SUBROUTINE INITALIZES MERGE FOR OPTICAL ******
C ****** DEPTHS ******
C ***********************************************************
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
C
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
C FOR AEROSOL RUNS, MOVE EXTID INTO YID
C
IF (IAERSL.GT.0) THEN
WRITE (CEXT,'(10A4)') EXTID
WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CYID(19:40) = CEXT(19:40)
READ (CYID,'(5A8)') (YID(I),I=3,7)
ENDIF
C
IEMIT = 0
FACT = 1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
DO 10 MOL = 1, NMOL
WK(MOL) = WK(MOL)*FACT
10 CONTINUE
WBROAD = WBROAD*FACT
LAYR1 = LAYER
WRITE (IPR,905) LAYR1,KFILE,MFILE
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
C
20 CONTINUE
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) GO TO 40
CALL BUFIN (KFILE,KEOF,ODLAYR(1),NLIM)
C
IF (FACT.EQ.2.) THEN
DO 30 I = 1, NLIM
ODLAYR(I) = ODLAYR(I)+ODLAYR(I)
30 CONTINUE
ENDIF
C
CALL ABSOUT (V1P,V2P,DVP,NLIM,1,MFILE,NPTS,ODLAYR,NPANLS)
GO TO 20
C
40 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,910) TIME1,TIM
C
RETURN
C
900 FORMAT ('0 THE TIME AT THE START OF ABINIT IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,/,'0 FILE ',I5,
* ' INITIALIZED ONTO FILE',I5)
910 FORMAT ('0 THE TIME AT THE END OF ABINIT IS ',F12.3/F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE ')
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE OPNMRG(LFILE,PATH1,L1,FM1,PATH2,L2,FM2,MFILE, 4
* PATH3,FM3)
C
C This subroutine opens file for merging in ABSMRG, when IMRG = 10
C
LOGICAL OP
CHARACTER*57 FILE1,FILE2,FILE3
CHARACTER*55 PATH1,PATH2,PATH3
CHARACTER*11 CFORM
CHARACTER*10 FM1,FM2,FM3
C
COMMON /MANE/ P0,TEMP0,NLAYRS,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYRFX,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
C 123456789-123456789-123456789-123456789-123456789-1234567
DATA FILE1 /
* ' '/,
* FILE2 /
* ' '/,
* FILE3 /
* ' '/
DATA CFORM / 'UNFORMATTED' /
C
INQUIRE (UNIT=LFILE,OPENED=OP)
IF (OP) CLOSE (LFILE)
WRITE(FILE1,FM1) PATH1,L1
OPEN(UNIT=LFILE,FILE=FILE1,FORM=CFORM,STATUS='OLD')
C
IF (L2.NE.L1) THEN
INQUIRE (UNIT=KFILE,OPENED=OP)
IF (OP) CLOSE (KFILE)
WRITE(FILE2,FM2) PATH2,L2
OPEN(UNIT=KFILE,FILE=FILE2,FORM=CFORM,STATUS='OLD')
ELSE
FILE2 = 'NO FILE USED'
ENDIF
C
INQUIRE (UNIT=MFILE,OPENED=OP)
IF (OP) CLOSE (MFILE)
WRITE(FILE3,FM3) PATH3,L2
OPEN(UNIT=MFILE,FILE=FILE3,FORM=CFORM,STATUS='UNKNOWN')
C
C Write procedure
C
WRITE(IPR,900) FILE1,FILE2,FILE3
C
RETURN
C
900 FORMAT (' Merged file: ',A57,/,
* ' With file: ',A57,/,
* ' To obtain file: ',A57,/)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE ABSMRG (NPTS,LFILE,MFILE,JPATHL) 2,19
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE ABSMRG MERGES ABSORPTION VALUES FROM KFILE INTO LFILE
C
COMMON R1(2410),OLDR1(2410)
COMMON /MANE/ P0,TEMP0,NLAYRS,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYRFX,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /ABSHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /ABSPNL/ V1P,V2P,DVP,NLIM,NSHFT,NPNTS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
DIMENSION SAVOR1(5),A1(10),A2(10),A3(10),A4(10)
DIMENSION WKSAV(35)
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
C
CHARACTER*40 CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO) , (FSCDID(4),IAERSL),
* (FSCDID(5),IEMIT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
CALL CPUTIM (TIME)
NPANLS = 0
IF (NOPR.EQ.0) WRITE (IPR,900) TIME
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
DVL = DV
LAY1SV = LAYR1
PL = PAVE
TL = TAVE
WTOTL = 0.
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
DVK = DV
LAYR1 = LAY1SV
FACT = 1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
IF (DVL.EQ.DVK) ITYPE = 0
IF (DVL.GT.DVK) ITYPE = DVK/(DVL-DVK)+0.5
IF (DVL.LT.DVK) ITYPE = -INT(DVL/(DVK-DVL)+0.5)
SAVOR1(4) = 0.
C
C ITYPE .LT. 0 IF DV(K-1) IS LESS THAN DV(K)
C
IF (ITYPE.LT.0) STOP ' ABSMRG: ITYPE LT 0 '
ITYPE = IABS(ITYPE)
WTOTK = 0.
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+FACT*WK(MOL)
WK(MOL) = FACT*WK(MOL)+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+FACT*WBROAD
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
WBROAD = FACT*WBROAD+WN2SAV
SECANT = 0.
IF (NOPR.EQ.0) WRITE (IPR,905) LAYR1,LAYER,KFILE,LFILE,MFILE
IEMIT = 0
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
DO 30 K = 1, 5
SAVOR1(K) = 0.
30 CONTINUE
ATYPE = ITYPE
AP = 1.0/(ATYPE+1.0)
IF (ITYPE.NE.0) GO TO 80
C
C 1/1 RATIO ONLY
C
40 CONTINUE
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 50
CALL BUFIN (LFILE,LEOF,OLDR1(1),NLIMO)
50 CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) GO TO 250
CALL BUFIN (KFILE,KEOF,R1(1),NLIM)
C
IF (FACT.EQ.1) THEN
DO 60 KOD = 1, NLIM
R1(KOD) = R1(KOD)+OLDR1(KOD)
60 CONTINUE
C
ELSE
DO 70 KOD = 1, NLIM
R1(KOD) = R1(KOD)+R1(KOD)+OLDR1(KOD)
70 CONTINUE
ENDIF
C
CALL ABSOUT (V1P,V2P,DVP,NLIM,1,MFILE,NPTS,R1,NPANLS)
C
GO TO 40
C
C ALL RATIOS EXCEPT 1/1
C
80 LL = ITYPE+1
DO 90 JPG = 1, ITYPE
APG = JPG
P = 1.0-(AP*APG)
C
C THE FOLLOWING ARE THE CONSTANTS FOR THE LAGRANGE 4 POINT
C INTERPOLATION.
C
A1(JPG) = -P*(P-1.0)*(P-2.0)/6.0
A2(JPG) = (P**2-1.0)*(P-2.0)*0.5
A3(JPG) = -P*(P+1.0)*(P-2.0)*0.5
A4(JPG) = P*(P**2-1.0)/6.0
90 CONTINUE
C
C ******** BEGINNING OF LOOP THAT DOES INTERPOLATION *********
C
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
CALL BUFIN (LFILE,LEOF,OLDR1(1),NLIMO)
MAXLF = NLIMO
NVS = 1
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
CALL BUFIN (KFILE,KEOF,R1(1),NLIM)
IF (KEOF.LE.0) GO TO 250
II = 1
DIF = DVP*0.01
IF (ABS(V1PO-V1P).LT.DIF) GO TO 120
C
C V1P < V1PO LASER OPTION
C
100 NVS = NVS+1
V1PN = V1PO+DVPO*(NVS-1)
IF (ABS(V1PN-V1P).LT.DIF) GO TO 120
IF (V1PN.LT.V1P) GO TO 100
C
110 II = II+1
V1PP = V1P+DVP*(II-1)
IF (ABS(V1PN-V1PP).LT.DIF) GO TO 120
IF (V1PP.LT.V1PN) GO TO 110
C
GO TO 100
120 R1(II) = FACT*R1(II)+OLDR1(NVS)
V1PN = V1PO+DVPO*(NVS-1)
V1PP = V1P+DVP*(II-1)
II = II+1
130 JJ = 1
C
DO 240 JPG = 1, LL
IF (JPG.EQ.LL) GO TO 140
IF (NVS.EQ.1) GO TO 150
GO TO 170
140 IF (FACT.EQ.1.) THEN
R1(II) = R1(II)+OLDR1(NVS)
ELSE
R1(II) = R1(II)+R1(II)+OLDR1(NVS)
ENDIF
V1PN = V1PO+DVPO*(NVS-1)
V1PP = V1P+DVP*(II-1)
GO TO 190
150 IF (SAVOR1(4).EQ.0.0) GO TO 160
OLDR1Y = SAVOR1(4)
GO TO 180
160 OLDR1Y = OLDR1(1)
GO TO 180
170 OLDR1Y = OLDR1(NVS-1)
180 OLDR1I = A1(JJ)*OLDR1Y+A2(JJ)*OLDR1(NVS)+A3(JJ)*OLDR1(NVS+1)+
* A4(JJ)*OLDR1(NVS+2)
IF (FACT.EQ.1.) THEN
R1(II) = R1(II)+OLDR1I
ELSE
R1(II) = R1(II)+R1(II)+OLDR1I
ENDIF
190 NVS = NVS+1
IF (NVS.LE.MAXLF-2) GO TO 200
SAVOR1(1) = OLDR1(NVS-1)
SAVOR1(2) = OLDR1(NVS)
SAVOR1(3) = OLDR1(NVS+1)
SAVOR1(4) = OLDR1(NVS-2)
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 210
MAXLF = NLIMO+3
CALL BUFIN (LFILE,LEOF,OLDR1(4),NLIMO)
OLDR1(1) = SAVOR1(1)
OLDR1(2) = SAVOR1(2)
OLDR1(3) = SAVOR1(3)
NVS = 2
200 II = II+1
JJ = JJ+1
IF (II.GT.NLIM) GO TO 230
GO TO 240
210 II = II+1
AVRG = (SAVOR1(3)+SAVOR1(2))*0.5
220 R1(II) = FACT*R1(II)+AVRG
II = II+1
IF (II.LE.NLIM) GO TO 220
C
C WRITE OUTPUT FILE
C
230 CALL ABSOUT (V1P,V2P,DVP,NLIM,1,MFILE,NPTS,R1,NPANLS)
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) GO TO 250
CALL BUFIN (KFILE,KEOF,R1(1),NLIM)
II = 1
240 CONTINUE
NVS = NVS-1
GO TO 130
250 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
IF (NOPR.EQ.0) WRITE (IPR,910) TIME1,TIM
RETURN
C
900 FORMAT ('0 THE TIME AT THE START OF ABSMRG IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5)
910 FORMAT (' THE TIME AT THE END OF ABSMRG IS',F12.3/F12.3,
* ' SECS. WERE REQUIRED FOR THIS ADDITION')
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE ABSINT (NPTS,LFILE,MFILE,JPATHL) 1,17
C
IMPLICIT REAL*8 (V)
C
COMMON NEWOD(2410),ODLAY(-2:2407)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYDUM,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /ABSHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /ABSPNL/ V1P,V2P,DVP,NLIM,NSHFT,NPNTS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
DIMENSION A1(0:100),A2(0:100),A3(0:100),A4(0:100)
DIMENSION OLDOD(2),ODLAYR(2)
DIMENSION WKSAV(35)
C
CHARACTER*40 CYID
REAL NEWOD
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO)
EQUIVALENCE (NEWOD(1),OLDOD(1)) , (ODLAY(1),ODLAYR(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
C ***********************************************************
C ****** THIS SUBROUTINE DOES LAYER MERGE FOR OPTICAL ******
C ****** DEPTHS USING FOUR POINT GENERAL INTERPOLATION ******
C ***********************************************************
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
C
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
DVL = DV
LAY1SV = LAYR1
PL = PAVE
TL = TAVE
WTOTL = 0.
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
DVK = DV
LAYR1 = LAY1SV
FACT = 1.
C
C IF(IPATHL.EQ.2 .AND. IANT.EQ.0) FACT=2.
C
IF (IPATHL.EQ.2.AND.IANT.EQ.0) STOP ' ABSINT: FACT=2. '
ATYPE = 9.999E09
IF (DVK.EQ.DVL) ATYPE = 0.
IF (DVL.GT.DVK) ATYPE = DVK/(DVL-DVK)+0.5
IF (DVL.LT.DVK) ATYPE = -DVL/(DVK-DVL)-0.5
IF (ATYPE.GT.0) STOP ' ABSINT; ATYPE GT 0 '
WTOTK = 0.
WRITE (IPR,905) LAYR1,LAYER,KFILE,LFILE,MFILE,ATYPE
IEMIT = 0
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+WK(MOL)*FACT
WK(MOL) = WK(MOL)*FACT+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+WBROAD*FACT
WBROAD = WBROAD*FACT+WN2SAV
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
SECANT = 0.
DV = DVL
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
C
IF (ATYPE.EQ.0.) THEN
C
C 1/1 RATIO ONLY
C
30 CONTINUE
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) GO TO 90
CALL BUFIN (KFILE,KEOF,ODLAYR(1),NLIM)
C
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
CALL BUFIN (LFILE,LEOF,OLDOD(1),NLIMO)
C
DO 40 I = 1, NLIM
NEWOD(I) = ODLAYR(I)+OLDOD(I)
40 CONTINUE
CALL ABSOUT (V1PO,V2PO,DVPO,NLIMO,1,MFILE,NPTS,NEWOD,NPANLS)
GO TO 30
C
ENDIF
C
C ALL RATIOS EXCEPT 1/1
C
DO 50 JP = 0, 100
APG = JP
P = 0.01*APG
C
C THE FOLLOW ARE THE CONSTANTS FOR THE LAGRANGE 4 POINT
C INTERPOLATION
C
A1(JP) = -P*(P-1.0)*(P-2.0)/6.0
A2(JP) = (P**2-1.0)*(P-2.0)*0.5
A3(JP) = -P*(P+1.0)*(P-2.0)*0.5
A4(JP) = P*(P**2-1.0)/6.0
50 CONTINUE
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) GO TO 90
CALL BUFIN (KFILE,KEOF,ODLAYR(1),NLIM)
C
ODLAY(-2) = ODLAY(1)
ODLAY(-1) = ODLAY(1)
ODLAY(0) = ODLAY(1)
C
RATDV = DVL/DVK
C
60 CONTINUE
C
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 90
CALL BUFIN (LFILE,LEOF,OLDOD(1),NLIMO)
C
C FJJ IS OFFSET BY 2. FOR ROUNDING PURPOSES
C
FJ1DIF = (V1PO-V1P)/DVP+1.+2.
C
C ***** BEGINNING OF LOOP THAT DOES MERGE *****
C
DO 80 II = 1, NLIMO
C
70 CONTINUE
C
FJJ = FJ1DIF+RATDV*FLOAT(II-1)
JJ = IFIX(FJJ)-2
C
IF (JJ+2.GT.NLIM) THEN
ODLAY(-2) = ODLAY(NLIM-2)
ODLAY(-1) = ODLAY(NLIM-1)
ODLAY(0) = ODLAY(NLIM)
V1PST = V1P
V2PST = V2P
NLIMST = NLIM
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
C
IF (KEOF.LE.0) THEN
V1P = V1PST
DVP = DVK
V2P = V2PST+2.*DVP
NLIM = NLIMST+2
ODLAY(NLIM-1) = ODLAY(NLIM-2)
ODLAY(NLIM) = ODLAY(NLIM-2)
ELSE
CALL BUFIN (KFILE,KEOF,ODLAYR(1),NLIM)
ENDIF
C
FJ1DIF = (V1PO-V1P)/DVP+1.+2.
GO TO 70
ENDIF
C
C JP = (FJJ-FLOAT(JJ))*100. + 0.5 - 200.
C
JP = (FJJ-FLOAT(JJ))*100.-199.5
IF (JP.GT.100) THEN
WRITE (IPR,910) JP,JJ,NLIM
STOP
ENDIF
C
C INTERPOLATE THE OLD TRANSMISSION
C
ODLAYI = A1(JP)*ODLAY(JJ-1)+A2(JP)*ODLAY(JJ)+
* A3(JP)*ODLAY(JJ+1)+A4(JP)*ODLAY(JJ+2)
IF (ODLAYI.LT.0.) ODLAYI = 0.
C
NEWOD(II) = ODLAYI+OLDOD(II)
C
80 CONTINUE
C
CALL ABSOUT (V1PO,V2PO,DVPO,NLIMO,1,MFILE,NPTS,NEWOD,NPANLS)
C
GO TO 60
C
90 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,915) TIME1,TIM
C
RETURN
C
900 FORMAT ('0 THE TIME AT THE START OF ABSINT IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5,' WITH XTYPE=',
* G15.5)
910 FORMAT ('0 JP, JJ, NLIM ',3I6)
915 FORMAT ('0 THE TIME AT THE END OF ABSINT IS ',F12.3/F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE ')
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE ABSOUT (V1PO,V2PO,DVPO,NLIMO,JLO,MFILE,NPTS,R1,NPANLS) 6,2
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE ABSOUT OUPUTS THE MERGED RESULT (R1) ONTO MFILE
C
COMMON /ABSPNI/ V1P,V2P,DVP,NLIM
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
DIMENSION PNLHDR(2),R1(*)
C
EQUIVALENCE (PNLHDR(1),V1P)
C
V1P = V1PO
V2P = V2PO
DVP = DVPO
NLIM = NLIMO
C
NPANLS = NPANLS+1
CALL BUFOUT (MFILE,PNLHDR(1),NPHDRF)
CALL BUFOUT (MFILE,R1(JLO),NLIM)
IF (NPTS.LE.0) GO TO 20
IF (NPANLS.EQ.1) WRITE (IPR,900)
WRITE (IPR,905)
NNPTS = NPTS
IF (NPTS.GT.(NLIM/2)+1) NNPTS = (NLIM/2)+1
JHILIM = JLO+NLIM-NNPTS
DO 10 I = 1, NNPTS
J = JLO+I-1
K = JHILIM+I-1
VJ = V1P+FLOAT(J-JLO)*DVP
VK = V1P+FLOAT(K-JLO)*DVP
WRITE (IPR,910) J,VJ,R1(J),K,VK,R1(K)
10 CONTINUE
20 CONTINUE
C
RETURN
C
900 FORMAT ('0 ','LOCATION WAVENUMBER',2X,'OPT DPTH',27X,
* 'LOCATION WAVENUMBER',2X,'OPT DPTH')
905 FORMAT (' ')
910 FORMAT (I8,2X,F12.6,1P,E15.7,0P,20X,I8,2X,F12.6,1P,E15.7)
C
END
C
C ----------------------------------------------------------------
C
FUNCTION BBFN (VI,DVI,V2I,XKT,VINEW,BBDEL,BBLAST) 156
C
IMPLICIT REAL*8 (V)
C
C FUNCTION BBFN CALCULATES BLACK BODY FN FOR WAVENUMBER VALUE VI
C AND CALCULATES THE WAVENUMBER VALUE (VINEW) FOR NEXT BBFN CALC.
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 23 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
DATA FACTOR / 0.003 /
C
XVI = VI
XVIOKT = XVI/XKT
EXPNEG = EXP(-XVIOKT)
GNU2 = XVI*XVI
BG2 = XVIOKT*XVIOKT
C
C Initialize BBLAST for BBLAST negative
C
IF (BBLAST.LT.0.) THEN
IF (XKT.GT.0.0) THEN
IF (XVIOKT.LE.0.01) THEN
BBLAST = RADCN1*(XVI**2)*XKT/(1.+0.5*XVIOKT)
ELSEIF (XVIOKT.LE.80.0) THEN
BBLAST = RADCN1*(XVI**3)/(EXP(XVIOKT)-1.)
ELSE
BBLAST = 0.
ENDIF
ELSE
BBLAST = 0.
ENDIF
ENDIF
C
C SET BBFN EQUAL TO BLACK BODY FUNCTION AT VI
C
C BBLAST IS BBFN(VI) FOR EACH SUBSEQUENT CALL
C
BBFN = BBLAST
C
INTVLS = 1
DELTAV2 = V2I - VI
IF (XKT.GT.0.0) THEN
C
IF (XVIOKT.LE.0.01) THEN
IF (VINEW.GE.0.0) THEN
XDELT = (GNU2 * (4.+4.*XVIOKT + BG2))/
* (10.*BG2 - 24.*XVIOKT + 8.)
DELTAV = SQRT(ABS(FACTOR*XDELT))
IF (DELTAV .GT. DELTAV2) DELTAV = DELTAV2
INTVLS = (DELTAV)/DVI
INTVLS = MAX(INTVLS,1)
VINEW = VI+DVI*FLOAT(INTVLS)
ELSE
VINEW = ABS(VINEW)
INTVLS = (VINEW-VI)/DVI
INTVLS = MAX(INTVLS,1)
ENDIF
XVINEW = VINEW
C
BBNEXT = RADCN1*(XVINEW**2)*XKT/(1.+0.5*XVINEW/XKT)
ELSEIF (XVIOKT.LE.80.0) THEN
IF (VINEW.GE.0.0) THEN
FRONT = XVIOKT/(1.-EXPNEG)
BOX = 3.0 - FRONT
DELT2C = (1./GNU2)*(2.*BOX-FRONT*(1.+BOX-FRONT*EXPNEG))
DELTAV = SQRT(ABS(FACTOR/DELT2C))
IF (DELTAV .GT. DELTAV2) DELTAV = DELTAV2
INTVLS = (DELTAV)/DVI
INTVLS = MAX(INTVLS,1)
VINEW = VI+DVI*FLOAT(INTVLS)
ELSE
VINEW = ABS(VINEW)
C
C The following IF test added for cases where
C XVIOKT > 80 on one call, and XVIOKT < 80 on
C the next call (numerical artifact causing
C the change over the XVIOKT = 80 boundary)
C
IF (VINEW.EQ.6.0E+05) THEN
BBNEXT=0.
BBDEL = (BBNEXT-BBFN)/FLOAT(INTVLS)
BBLAST = BBNEXT
RETURN
ENDIF
INTVLS = (VINEW-VI)/DVI
INTVLS = MAX(INTVLS,1)
ENDIF
XVINEW = VINEW
C
BBNEXT = RADCN1*(XVINEW**3)/(EXP(XVINEW/XKT)-1.)
ELSE
BBNEXT = 0.
VINEW = 6.0E+5
ENDIF
ELSE
BBNEXT = 0.
VINEW = 6.0E+5
ENDIF
C
BBDEL = (BBNEXT-BBFN)/FLOAT(INTVLS)
BBLAST = BBNEXT
C
RETURN
C
END
C
C ---------------------------------------------------------------
C
FUNCTION BBAD(BBVAL,VI,DVI,V2I,XKT,VDNEW,BBADDL,BBADOL) 4
C
IMPLICIT REAL*8 (V)
C
C FUNCTION BBAD calculates the derivative of the black body fn
C analytically for wavenumber value VI
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C IMPLEMENTATION: P.D. Brown
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
DATA FACTOR / 0.003 /
C
XVI = VI
XVIOKT = XVI/XKT
EXPNEG = EXP(-XVIOKT)
GNU2 = XVI*XVI
BG2 = XVIOKT*XVIOKT
TKELV = XKT*RADCN2
C
C If first call, initialize BBADOL
C
C For linear approximation,
C
C BBFN{prime} = -BBFN/(t*(1 + 0.5*(planck*clight*gnu/boltz*t)))
C
C Using the exact function,
C
C BBFN{prime} = BBFN*(planck*clight*gnu/(boltz*t*t))*
C 1/[1-exp(-planck*clight*gnu/(boltz*t))]
C where
C planck*clight/boltz = RADCN2
C boltz*t/(planck*clight) = XKT
C planck*clight*gnu/boltz*t = XVIOKT
C
C and we can solve easily for t:
C t = XKT*RADCN2.
C
IF (BBADOL.LT.0.) THEN
IF (XKT.GT.0.0) THEN
IF (XVIOKT.LE.0.01) THEN
BBADOL = BBVAL/(-TKELV*(1.+0.5*XVIOKT))
ELSEIF (XVIOKT.LE.80.0) THEN
BBADOL = BBVAL*XVIOKT/(TKELV*(1-EXPNEG))
ELSE
BBADOL = 0.
ENDIF
ELSE
BBADOL = 0.
ENDIF
ENDIF
C
C Set BBAD equal to black body function derivative
C
C BBADOL is BBAD(VI) for each subsequent call
C
BBAD = BBADOL
C
INTVLS = 1
DELTAV2 = V2I - VI
IF (XKT.GT.0.0) THEN
C
IF (XVIOKT.LE.0.01) THEN
IF (VDNEW.GE.0.0) THEN
XDELT = (GNU2 * (4.+4.*XVIOKT + BG2))/
* (10.*BG2 - 24.*XVIOKT + 8.)
DELTAV = SQRT(ABS(FACTOR*XDELT))
IF (DELTAV .GT. DELTAV2) DELTAV = DELTAV2
INTVLS = (DELTAV)/DVI
INTVLS = MAX(INTVLS,1)
VDNEW = VI+DVI*FLOAT(INTVLS)
ELSE
VDNEW = ABS(VDNEW)
INTVLS = (VDNEW-VI)/DVI
INTVLS = MAX(INTVLS,1)
ENDIF
XVDNEW = VDNEW
BBAD = BBVAL/(-TKELV*(1.+0.5*XVDNEW/XKT))
ELSEIF (XVIOKT.LE.80.0) THEN
IF (VDNEW.GE.0.0) THEN
FRONT = XVIOKT/(1.-EXPNEG)
BOX = 3.- FRONT
DELT2C = (1./GNU2)*(2.*BOX-FRONT*(1.+BOX-FRONT*EXPNEG))
DELTAV = SQRT(ABS(FACTOR/DELT2C))
IF (DELTAV .GT. DELTAV2) DELTAV = DELTAV2
INTVLS = (DELTAV)/DVI
INTVLS = MAX(INTVLS,1)
VDNEW = VI+DVI*FLOAT(INTVLS)
ELSE
VDNEW = ABS(VDNEW)
C
C The following IF test added for cases where
C XVIOKT > 80 on one call, and XVIOKT < 80 on
C the next call (numerical artifact causing
C the change over the XVIOKT = 80 boundary)
C
IF (VDNEW.EQ.6.0E+05) THEN
BBAD=0.
VDNEW = VDNEW-DVI+0.00001
BBADDL = (BBAD-BBADOL)/FLOAT(INTVLS)
BBADOL = BBAD
RETURN
ENDIF
INTVLS = (VDNEW-VI)/DVI
INTVLS = MAX(INTVLS,1)
ENDIF
XVDNEW = VDNEW
BBAD = BBVAL*(XVDNEW/XKT)/(TKELV*(1-EXP(-XVDNEW/XKT)))
ELSE
BBAD = 0.
VDNEW = 6.0E+5
ENDIF
ELSE
BBAD = 0.
VDNEW = 6.0E+5
ENDIF
C
BBADDL = (BBAD-BBADOL)/FLOAT(INTVLS)
C
VDNEW = VDNEW-DVI+0.00001
BBADOL = BBAD
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
FUNCTION EMISFN (VI,DVI,VINEM,EMDEL,EMLAST) 22,1
C
IMPLICIT REAL*8 (V)
C
C FUNCTION EMISFN CALCULATES BOUNDARY EMISSIVITY FOR WAVE NUMBER
C VALUE CORRESPONDING TO VI AND VINEM, AND THEN CALCULATES THE
C LINEAR CHANGE BETWEEN THE EMISSIVITY VALUES AT VI AND VINEM
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 23 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C ----------------------------------------------------------------
C Parameter and common block for direct input of emission function
C values
C
PARAMETER (NMAXCO=4040)
COMMON /EMSFIN/ V1EMIS,V2EMIS,DVEMIS,NLIMEM,ZEMIS(NMAXCO)
C ----------------------------------------------------------------
C
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
EQUIVALENCE (BNDEMI(1),A) , (BNDEMI(2),B) , (BNDEMI(3),C)
C
DATA FACTOR / 0.001 /
C
C ***************************************************
C Check for A < 0. If so, use input values read in from file
C "EMISSION"
C ***************************************************
C
IF (A.LT.0.) THEN
C
C Determine elements of EMISSION function to use with
C input frequency
C
NELMNT = INT((VI-V1EMIS)/DVEMIS)
C
C Test for bounds on EMISSION function
C
IF ((NELMNT.LE.0).OR.(NELMNT.GE.NLIMEM)) THEN
WRITE(*,*) 'Frequency range of calculation exceeded',
* ' emissivity input.'
WRITE(*,*) ' VI = ',VI,' V1EMIS = ',V1EMIS,' V2EMIS = ',
* V2EMIS
STOP 'ERROR IN EMISFN'
ENDIF
C
C Interpolate to obtain appropriate EMISSION value
C
V1A = V1EMIS+DVEMIS*NELMNT
V1B = V1EMIS+DVEMIS*(NELMNT+1)
CALL LINTCO(V1A,ZEMIS(NELMNT),V1B,ZEMIS(NELMNT+1),VI,ZINT,
* ZDEL)
EMISFN = ZINT
VINEM = V1B
EMDEL = ZDEL*DVI
EMLAST = ZEMIS(NELMNT+1)
RETURN
C
ENDIF
C
C ***************************************************
C The following uses a quadratic formula for emission
C ***************************************************
C
C CHECK FOR CONSTANT E (INDEPENDENT OF VI)
C IF CONSTANT RETURN LARGE VALUE FOR VINEM
C
IF (B.EQ.0..AND.C.EQ.0.) THEN
EMISFN = A
VINEM = 9.99E+9
EMDEL = 0.0
EMLAST = EMISFN
RETURN
ENDIF
C
XVI = VI
IF (EMLAST.LT.0.) THEN
EMLAST = A+B*XVI+C*XVI*XVI
ENDIF
C
C SET EMISFN EQUAL TO EMISSIVITY AT VI
C
C EMLAST IS EMISFN(VI) FOR EACH SUBSEQUENT CALL
C
EMISFN = EMLAST
C
IF (VINEM.GE.0.0) THEN
XVNEXT = XVI+FACTOR/ABS((B+2.*C*XVI))
XVNEXT = MIN(XVNEXT,(XVI+DVI*2400))
INTVLS = (XVNEXT-XVI)/DVI
INTVLS = MAX(INTVLS,1)
XVNEXT = XVI+DVI*FLOAT(INTVLS)
ELSE
XVNEXT = ABS(VINEM)
INTVLS = (XVNEXT-XVI)/DVI
INTVLS = MAX(INTVLS,1)
ENDIF
C
EMNEXT = A+B*XVNEXT+C*XVNEXT*XVNEXT
C
EMDEL = (EMNEXT-EMISFN)/FLOAT(INTVLS)
C
VINEM = XVNEXT
EMLAST = EMNEXT
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
FUNCTION REFLFN (VI,DVI,VINRF,RFDEL,RFLAST) 10,1
C
IMPLICIT REAL*8 (V)
C
C FUNCTION REFLFN CALCULATES BOUNDARY REFLECTIVITY FOR WAVE NUMBER
C VALUE CORRESPONDING TO VI AND VINRF, AND THEN CALCULATES THE
C LINEAR CHANGE BETWEEN THE REFLECTIVITY VALUES AT VI AND VINRF
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 23 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C ----------------------------------------------------------------
C Parameter and common block for direct input of reflection
C function values
C
PARAMETER (NMAXCO=4040)
COMMON /RFLTIN/ V1RFLT,V2RFLT,DVRFLT,NLIMRF,ZRFLT(NMAXCO)
C ----------------------------------------------------------------
C
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
EQUIVALENCE (BNDRFL(1),A) , (BNDRFL(2),B) , (BNDRFL(3),C)
C
DATA FACTOR / 0.001 /
C
C ***************************************************
C Check for A < 0. If so, use input values read in from file
C "REFLECTION"
C ***************************************************
C
IF (A.LT.0.) THEN
C
C Determine elements of REFLECTION function to use with
C input frequency
C
NELMNT = INT((VI-V1RFLT)/DVRFLT)
C
C Test for bounds on REFLECTION function
C
IF ((NELMNT.LE.0).OR.(NELMNT.GE.NLIMRF)) THEN
WRITE(*,*) 'Frequency range of calculation exceeded',
* ' reflectivity input.'
WRITE(*,*) ' VI = ',VI,' V1RFLT = ',V1RFLT,' V2RFLT = ',
* V2RFLT
STOP 'ERROR IN REFLFN'
ENDIF
C
C Interpolate to obtain appropriate reflection value
C
V1A = V1RFLT+DVRFLT*NELMNT
V1B = V1RFLT+DVRFLT*(NELMNT+1)
CALL LINTCO(V1A,ZRFLT(NELMNT),V1B,ZRFLT(NELMNT+1),VI,ZINT,
* ZDEL)
REFLFN = ZINT
VINRF = V1B
RFDEL = ZDEL*DVI
RFLAST = ZRFLT(NELMNT+1)
RETURN
C
ENDIF
C
C ***************************************************
C The following uses a quadratic formula for emission
C ***************************************************
C
C CHECK FOR CONSTANT R (INDEPENDENT OF VI)
C IF CONSTANT RETURN LARGE VALUE FOR VINRF
C
IF (B.EQ.0..AND.C.EQ.0.) THEN
REFLFN = A
VINRF = 9.99E+9
RFDEL = 0.0
RFLAST = REFLFN
RETURN
ENDIF
C
XVI = VI
IF (RFLAST.LT.0.) THEN
RFLAST = A+B*XVI+C*XVI*XVI
ENDIF
C
C SET REFLFN EQUAL TO REFLECTIVITY AT VI
C
C RFLAST IS REFLFN(VI) FOR EACH SUBSEQUENT CALL
C
REFLFN = RFLAST
C
IF (VINRF.GE.0.0) THEN
XVNEXT = XVI+FACTOR/ABS((B+2.*C*XVI))
XVNEXT = MIN(XVNEXT,(XVI+DVI*2400))
INTVLS = (XVNEXT-XVI)/DVI
INTVLS = MAX(INTVLS,1)
XVNEXT = XVI+DVI*FLOAT(INTVLS)
ELSE
XVNEXT = ABS(VINRF)
INTVLS = (XVNEXT-XVI)/DVI
INTVLS = MAX(INTVLS,1)
ENDIF
C
RFNEXT = A+B*XVNEXT+C*XVNEXT*XVNEXT
C
RFDEL = (RFNEXT-REFLFN)/FLOAT(INTVLS)
C
VINRF = XVNEXT
RFLAST = RFNEXT
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE EMIN (V1P,V2P,DVP,NLIM,KFILE,EM,EMB,TR,KEOF,NPANLS) 4,70
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE EMIN INPUTS OPTICAL DEPTH VALUES FROM KFILE AND
C CALCULATES SOURCE FUNCTION FOR THE LAYER.
C THIS VERSION WORKS FOR AEROSOLS AND NLTE.
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KDUMY,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /LAMCHN/ ONEPL,ONEMI,EXPMIN,ARGMIN
COMMON /BUFPNL/ V1PBF,V2PBF,DVPBF,NLIMBF
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
DIMENSION PNLHDR(2),EM(*),EMB(*),TR(*)
C
EQUIVALENCE (FSCDID(1),IHIRAC) , (FSCDID(2),ILBLF4)
EQUIVALENCE (PNLHDR(1),V1PBF)
EQUIVALENCE (FSCDID(4),IAERSL)
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) RETURN
CALL BUFIN (KFILE,KEOF,TR(1),NLIMBF)
C
C TR CONTAINS THE OPTICAL DEPTHS AT THIS STAGE
C
IF (IHIRAC.EQ.4) CALL BUFIN (KFILE,KEOF,EM(1),NLIMBF)
C
C EM CONTAINS THE OPTICAL DEPTH CORRECTIONS FOR NLTE AT THIS STAGE
C
IF (NPANLS.LT.1.AND.IAERSL.EQ.0) WRITE (IPR,900)
IF (NPANLS.LT.1.AND.IAERSL.NE.0) WRITE (IPR,905)
C
EXT = 0.
ADEL = 0.
RADFN0 = 0.
RDEL = 0.
BB = 0.
BBDEL = 0.
BBA = 0.
BBDLA = 0.
BBB = 0.
BBDLB = 0.
C
V1P = V1PBF
V2P = V2PBF
DVP = DVPBF
NLIM = NLIMBF
VI = V1P-DVP
VIDV = VI
VIBB = VI
VAER = VI
VDUM = VI
BBLAST = -1.
BBLXTA = -2.
BBLXTB = -3.
RDLAST = -1.
BBDUM = -4.
RDDUM = -1.
NLIM1 = 0
NLIM2 = 0
C
AA = 0.278
C
IF (IAERSL.EQ.0) THEN
IAFBB = -1
ELSE
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBDUM)
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDDUM)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
IAFBB = 0
IF (VITST.LT.VAER.AND.VITST.LT.VIBB) IAFBB = 1
IF (VAER.LT.VITST.AND.VAER.LT.VIBB) IAFBB = 2
ENDIF
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER ONLY ONCE
C
C - WITH XKTA=0 THIS ALGORITHM REVERTS TO THE ORIGINAL
C
IF (XKTB.LE.0.) THEN
C
C - THIS SECTION TREATS THE LTE CASE
C
IF (IHIRAC.NE.4) THEN
C
10 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 20 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(BB+XX*BBA)/(1.+XX)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
20 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 10
ELSE
C
C - THIS SECTION TREATS THE NLTE CASE
C
30 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 40 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(1.0-EM(I)/ODVI)*(BB+XX*BBA)/(1.+XX)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
40 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 30
C
ENDIF
ELSE
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER TWICE:
C
C - FOR TANGENT PATHS AND FOR THE CASE OF THE REFLECTED ATMOSPHERE
C
IF (IHIRAC.NE.4) THEN
C
C - THIS SECTION TREATS THE LTE CASE
C
50 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 60 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EMX = (1.-TR(I))/(1.+XX)
EM(I) = EMX*(BB+XX*BBA)
EMB(I) = EMX*(BB+XX*BBB)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBB = BBB+BBDLB
60 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 50
C
ELSE
C
C - THIS SECTION TREATS THE CASE OF NLTE
C
70 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 80 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EMX = (1.-TR(I))*(1.0-EM(I)/ODVI)/(1.+XX)
EM(I) = EMX*(BB+XX*BBA)
EMB(I) = EMX*(BB+XX*BBB)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBB = BBB+BBDLB
80 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 70
C
ENDIF
ENDIF
C
RETURN
C
900 FORMAT ('0EMISSION AND TRANSMISSION (MOLECULAR) ')
905 FORMAT ('0EMISSION AND TRANSMISSION (AEROSOLS EFFECTS INCLUDED)')
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE EMINIT (NPTS,MFILE,JPATHL,TBND) 2,26
C
IMPLICIT REAL*8 (V)
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON NEWEM(2410),NEWTR(2410)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYER,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMIHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYDUM,YI1,YID(10),LSTWDF
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILA,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
CHARACTER*40 CEXT,CYID
C
DIMENSION EMLAYB(2410)
DIMENSION XFILHD(2),OPNLHD(2)
DIMENSION EMLAYR(2),TRLAYR(2)
C
EQUIVALENCE (XFILHD(1),XID(1)) , (OPNLHD(1),V1PO)
EQUIVALENCE (NEWEM(1),EMLAYR(1)) , (NEWTR(1),TRLAYR(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
REAL NEWEM,NEWTR
C
C
C *********************************************************************
C **** THIS SUBROUTINE COMPUTES THE EMISSION FOR THE FIRST LAYER ****
C *********************************************************************
C
C TBND IS THE BOUNDARY BLACK BODY TEMPERATUE
C
C IPATHL =-1 IS FOR THE LOOKING DOWN CASE WITH REFLECTED ATMOSPHERE
C IPATHL = 0 IS FOR THE HORIZONTAL PATH CASE (HOMOGENEOUS LAYER)
C IPATHL = 1 IS FOR THE LOOKING DOWN CASE (TO DENSER LAYERS)
C IPATHL = 2 IS FOR THE SYMMETRIC TANGENT PATH CASE
C IPATHL = 3 IS FOR THE LOOKING UP CASE (TO LESS DENSE LAYERS
C
CALL CPUTIM (TIME)
C
C ** NOTE ON IPATHL =2
C THE TANGENT MERGE ROUTINES ARE DIVIDED INTO ANTERIOR (1ST)
C AND POSTERIOR (2ND) LAYER CROSSINGS. THIS RECOGNITION IS
C TRIGGERED BY THE VALUE OF "IANT"
C
C IF IANT.EQ. 1 THEN POSTERIOR MERGE
C IF IANT.EQ. 0 THEN NORMAL MERGE
C IF IANT.EQ.-1 THEN ANTERIOR MERGE
C
WRITE (IPR,900) TIME
NPANLS = 0
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
C FOR AEROSOL RUNS, MOVE EXTID INTO YID
C
IF (IAERSL.GT.0) THEN
WRITE (CEXT,'(10A4)') EXTID
WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CYID(19:40) = CEXT(19:40)
READ (CYID,'(5A8)') (YID(I),I=3,7)
ENDIF
C
C IF BOUNDARY PROPERTIES ARE SUPPLIED, AND DOWNWARD LOOKING
C CASE; SET IPATHL TO REFLECTED ATMOSPHERE CASE
C
IF (IBPROP.EQ.1.AND.IPATHL.EQ.1) IPATHL = -1
IEMIT = 1
FACT = 1.
TIMEM = 0.0
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
DO 10 MOL = 1, NMOL
WK(MOL) = WK(MOL)*FACT
10 CONTINUE
WBROAD = WBROAD*FACT
LAYR1 = LAYER
WRITE (IPR,905) LAYR1,LAYER,KFILE,MFILE
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
XKT = TAVE/RADCN2
XKTBND = TBND/RADCN2
IF (IPATHL.EQ.-1) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ENDIF
IF (IPATHL.EQ.0) THEN
XKTA = 0.
XKTB = 0.
ENDIF
IF (IPATHL.EQ.1) THEN
XKTA = TZU/RADCN2
XKTB = 0.
ENDIF
IF (IPATHL.EQ.2) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ENDIF
IF (IPATHL.EQ.3) THEN
XKTA = TZL/RADCN2
XKTB = 0.
ENDIF
WRITE (IPR,910) IPATHL,IANT
C
20 CONTINUE
C
CALL CPUTIM (TIMEM1)
CALL EMIN (V1PO,V2PO,DVPO,NLIMO,KFILE,EMLAYR,EMLAYB,TRLAYR,KEOF,
* NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 80
VI = V1PO-DVPO
VIDVBD = VI
VIDVEM = VI
VIDVRF = VI
BBLAST = -1.
EMLAST = -1.
RFLAST = -1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) THEN
DO 30 J = 1, NLIMO
TRJ = TRLAYR(J)
NEWEM(J) = EMLAYR(J)+EMLAYB(J)*TRJ
TRLAYR(J) = TRLAYR(J)*TRJ
30 CONTINUE
ELSEIF ((IPATHL.EQ.1).AND.(TBND.GT.0.)) THEN
C
NLIM1 = 0
NLIM2 = 0
EMDUM = 0.
BBDUM = 0.
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMDUM)
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBDUM)
IEMBB = 0
IF (VIDVBD.GT.VIDVEM) IEMBB = 1
C
40 NLIM1 = NLIM2+1
C
VI = V1PO+FLOAT(NLIM1-1)*DVPO
IF (IEMBB.EQ.0) THEN
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDV,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
ELSE
EMISIV = EMISFN(VI,DVPO,VIDV,EMDEL,EMLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
ENDIF
C
IF (VIDV.GE.9.E+4) THEN
NLIM2 = NLIMO+1
ELSE
NLIM2 = (VIDV-V1PO)/DVPO+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIMO)
C
DO 50 J = NLIM1, NLIM2
V=VI+FLOAT(J-1)*DVPO
NEWEM(J) = EMLAYR(J)+TRLAYR(J)*BB*EMISIV
C
C Increment interpolation value
C
EMISIV = EMISIV+EMDEL
BB = BB+BBDEL
50 CONTINUE
C
IF (NLIM2.LT.NLIMO) GO TO 40
C
ELSEIF ((IPATHL.EQ.-1).AND.(TBND.GT.0.)) THEN
C
NLIM1 = 0
NLIM2 = 0
EMDUM = 0.
RFDUM = 0.
BBDUM = 0.
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMDUM)
REFLCT = REFLFN(VI,DVPO,VIDVRF,RFDEL,RFDUM)
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBDUM)
IEMBB = 0
IF (VIDVEM.LT.VIDVRF.AND.VIDVEM.LT.VIDVBD) IEMBB = 1
IF (VIDVRF.LT.VIDVEM.AND.VIDVRF.LT.VIDVBD) IEMBB = 2
C
60 NLIM1 = NLIM2+1
C
VI = V1PO+FLOAT(NLIM1-1)*DVPO
IF (IEMBB.EQ.0) THEN
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDV,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
VIDVRF = -VIDV
REFLCT = REFLFN(VI,DVPO,VIDVRF,RFDEL,RFLAST)
ELSEIF (IEMBB.EQ.1) THEN
EMISIV = EMISFN(VI,DVPO,VIDV,EMDEL,EMLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
VIDVRF = -VIDV
REFLCT = REFLFN(VI,DVPO,VIDVRF,RFDEL,RFLAST)
ELSE
REFLCT = REFLFN(VI,DVPO,VIDV,RFDEL,RFLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
ENDIF
C
IF (VIDV.GE.9.E+4) THEN
NLIM2 = NLIMO+1
ELSE
NLIM2 = (VIDV-V1PO)/DVPO+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIMO)
C
DO 70 J = NLIM1, NLIM2
V=VI+FLOAT(J-1)*DVPO
NEWEM(J) = EMLAYR(J)+EMLAYB(J)*REFLCT*TRLAYR(J)+
* TRLAYR(J)*BB*EMISIV
C
C Increment interpolation value
C
EMISIV = EMISIV+EMDEL
REFLCT = REFLCT+RFDEL
BB = BB+BBDEL
70 CONTINUE
C
IF (NLIM2.LT.NLIMO) GO TO 60
C
ENDIF
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,NEWEM,NEWTR,MFILE,NPTS,NPANLS)
GO TO 20
80 CALL CPUTIM (TIME1)
TIME = TIME1-TIME
WRITE (IPR,915) TIME,TIMEM
C
RETURN
C
900 FORMAT (' TIME AT THE START OF --EMINIT-- ',F10.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,
* '0 INPUT FILE =',I5,' OUTPUT FILE =',I5)
910 FORMAT ('0 IPATHL AND IANT',2I5)
915 FORMAT (' TIME REQUIRED FOR --EMINIT-- ',F10.3,
* ' --EMIN-- ',F10.3)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE RADMRG (NPTS,LFILE,MFILE,JPATHL,TBND) 1,18
C
IMPLICIT REAL*8 (V)
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 8 APRIL 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON RADN(2410),TRAN(2410),RADO(0:5000)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYDUM,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /XPANEL/ V1P,V2P,DVP,NLIM,RMIN,RMAX,NPNLXP,NSHIFT,NPTSS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /XME/ TRAO(0:5000)
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
DIMENSION RADLYB(2410)
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
DIMENSION A1(10),A2(10),A3(10),A4(10)
DIMENSION RADLYR(2),TRALYR(2),RADOI(2),TRAOI(2)
DIMENSION WKSAV(35)
C
CHARACTER*40 CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO)
EQUIVALENCE (RADO(1),RADOI(1)) , (TRAO(1),TRAOI(1)),
* (RADN(1),RADLYR(1)) , (TRAN(1),TRALYR(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
DATA NDIM / 2410 /,ND2 / 5000 /
C
C
C
C ************************************************************
C ****** THIS SUBROUTINE DOES LAYER MERGE FOR RADIANCE ******
C ************************************************************
C
C IPATHL =-1 IS FOR THE LOOKING DOWN CASE FOR REFLECTED ATMOSPHERE
C IPATHL = 1 IS FOR THE LOOKING DOWN CASE (TO DENSER LAYERS)
C IPATHL = 2 IS FOR THE SYMMETRIC TANGENT PATH CASE
C IPATHL = 3 IS FOR THE LOOKING UP CASE (TO LESS DENSE LAYERS)
C
C
C ** NOTE ON IPATHL = 2
C THE TANGENT MERGE ROUTINES ARE DIVIDED INTO ANTERIOR (1ST)
C AND POSTERIOR (2ND) LAYER CROSSINGS THIS RECOGNITION IS
C TRIGGERED BY THE VALUE OF "IANT"
C
C IF IANT.EQ. 1 THEN POSTERIOR MERGE
C IF IANT.EQ. 0 THEN NORMAL MERGE
C IF IANT.EQ.-1 THEN ANTERIOR MERGE
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
TIMEM = 0.0
TIMRD = 0.0
TIMOT = 0.0
C
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
LAY1SV = LAYR1
DVL = DV
PL = PAVE
TL = TAVE
WTOTL = 0.
C
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
C
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
C
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
C IF BOUNDARY PROPERTIES ARE SUPPLIED, AND DOWNWARD LOOKING
C CASE; SET IPATHL TO REFLECTED ATMOSPHERE CASE
C
IF (IBPROP.EQ.1.AND.IPATHL.EQ.1) IPATHL = -1
TAVK = TAVE
DVK = DV
FACT = 1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
C
IF (DVL.EQ.DVK) THEN
ITYPE = 0
ELSEIF (DVL.GT.DVK) THEN
ITYPE = DVK/(DVL-DVK)+0.5
ELSE
C
C DVL.LT.DVK
C
ITYPE = -INT(DVL/(DVK-DVL)+0.5)
ENDIF
IF (ITYPE.LT.0) STOP ' RADMRG; ITYPE LT 0 '
C
WTOTK = 0.
LAYR1 = LAY1SV
WRITE (IPR,905) LAYR1,LAYER,KFILE,LFILE,MFILE
IEMIT = 1
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+WK(MOL)*FACT
WK(MOL) = WK(MOL)*FACT+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+WBROAD*FACT
WBROAD = WBROAD*FACT+WN2SAV
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
SECANT = 0.
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
XKT = TAVK/RADCN2
C
WRITE (IPR,910) IPATHL,IANT
C
IF (IPATHL.EQ.-1) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ELSEIF (IPATHL.EQ.1) THEN
XKTA = TZU/RADCN2
XKTB = 0.
ELSEIF (IPATHL.EQ.2) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ELSEIF (IPATHL.EQ.3) THEN
XKTA = TZL/RADCN2
XKTB = 0.
ELSE
STOP ' RADMRG; IPATHL '
ENDIF
C
ATYPE = ITYPE
LL = ITYPE+1
AP = 1.0/(ATYPE+1.0)
C
C A1, A2, A3 AND A4 ARE THE CONSTANTS
C FOR THE LAGRANGE 4 POINT INTERPOLATION
C
DO 30 JPG = 1, ITYPE
APG = JPG
IPL = JPG+1
P = 1.0-(AP*APG)
A1(IPL) = -P*(P-1.0)*(P-2.0)/6.0
A2(IPL) = (P**2-1.0)*(P-2.0)*0.5
A3(IPL) = -P*(P+1.0)*(P-2.0)*0.5
A4(IPL) = P*(P**2-1.0)/6.0
30 CONTINUE
C
C *** BEGINNING OF LOOP THAT DOES MERGE ***
C
NPE = 0
RADO(0) = 0.0
TRAO(0) = 0.0
V1PO = 0.0
V2PO = 0.0
DVPO = 0.0
C
40 CONTINUE
C
CALL CPUTIM (TIMEM1)
CALL EMIN (V1P,V2P,DVP,NLIM,KFILE,RADLYR,RADLYB,TRALYR,KEOF,
* NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 80
II = 1
C
IF (V2PO.LE.V2P+DVPO) THEN
50 CALL CPUTIM (TIMEM1)
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 60
CALL BUFIN (LFILE,LEOF,RADOI(NPE+1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAOI(NPE+1),NLIMO)
CALL CPUTIM (TIMEM2)
TIMRD = TIMRD+TIMEM2-TIMEM1
NPE = NLIMO+NPE
IF (V2PO.LE.V2P+DVPO) GO TO 50
ENDIF
C
C ZERO POINT OF FIRST PANEL
C
60 IF (RADO(0).EQ.0.0.AND.TRAO(0).EQ.0.0) THEN
RADO(0) = RADO(1)
TRAO(0) = TRAO(1)
ENDIF
C
C END POINT OF LAST PANEL
C
IF (V2PO+DVPO.GE.V2) THEN
RADO(NPE+1) = RADO(NPE)
RADO(NPE+2) = RADO(NPE)
TRAO(NPE+1) = TRAO(NPE)
TRAO(NPE+2) = TRAO(NPE)
ENDIF
C
NPL = 1
C
C NPL IS LOCATION OF FIRST ELEMENT ON ARRAYS RADO AND TRAO
C
CALL RADNN (RADN,TRAN,RADO,TRAO,RADLYB,NLIM,NDIM,ND2,V1P,DVP,
* IPATHL,A1,A2,A3,A4,LL,NPL)
C
CALL CPUTIM (TIMEM1)
C
IF (TBND.GT.0.) CALL EMBND (V1P,V2P,DVP,NLIM,RADN,TRAN,TBND)
C
CALL EMOUT (V1P,V2P,DVP,NLIM,RADN,TRAN,MFILE,NPTS,NPANLS)
CALL CPUTIM (TIMEM2)
TIMOT = TIMOT+TIMEM2-TIMEM1
C
C NPL IS NOW LOCATION OF FIRST ELEMENT TO BE USED FOR NEXT PASS
C
IPL = -1
DO 70 NL = NPL, NPE
IPL = IPL+1
RADO(IPL) = RADO(NL)
TRAO(IPL) = TRAO(NL)
70 CONTINUE
C
NPE = IPL
C
GO TO 40
80 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,915) TIME1,TIM,TIMEM,TIMRD,TIMOT
C
RETURN
C
C
900 FORMAT ('0 THE TIME AT THE START OF RADMRG IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5)
910 FORMAT ('0 IPATHL AND IANT',2I5)
915 FORMAT ('0 THE TIME AT THE END OF RADMRG IS ',F12.3/F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE - EMIN - ',
* F12.3,' - READ - ',F12.3,' - EMOUT - ',F12.3)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE RADNN (RADLYR,TRALYR,RADO,TRAO,RADLYB,NLIM,NDIM,ND2, 2,2
* V1P,DVP,IPATHL,A1,A2,A3,A4,LL,NPL)
C
IMPLICIT REAL*8 (V)
C
C THIS SUBROUTINE CALCULATES THE NEW RADIANCE AND TRANSMISSION
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 8 APRIL 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM: R.D. WORSHAM
C S.A. CLOUGH
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
DIMENSION RADLYR(NDIM),TRALYR(NDIM),RADO(0:ND2),TRAO(0:ND2),
* RADLYB(NDIM),A1(*),A2(*),A3(*),A4(*)
C
LLM1 = LL-1
LLM1 = MAX(LLM1,1)
C
C LOOPING OVER POINTS WITH SAME WEIGHTS
C
DO 110 NL = 1, LL
IPL = (NPL+NL-1)-LLM1
IF (NL.GT.1) IPL = IPL-1
C
IF (NL.EQ.1) THEN
C
C EXACT FREQUENCY - NO INTERPOLATION
C
IF (IPATHL.EQ.1) THEN
C
DO 10 I = NL, NLIM, LL
IPL = IPL+LLM1
RADLYR(I) = TRALYR(I)*RADO(IPL)+RADLYR(I)
TRALYR(I) = TRALYR(I)*TRAO(IPL)
10 CONTINUE
C
ELSEIF (IPATHL.EQ.2) THEN
C
DO 20 I = NL, NLIM, LL
IPL = IPL+LLM1
TRTEMP = TRALYR(I)*TRAO(IPL)
RADLYR(I) = TRALYR(I)*RADO(IPL)+RADLYR(I)+
* RADLYB(I)*TRTEMP
TRALYR(I) = TRALYR(I)*TRTEMP
20 CONTINUE
C
ELSEIF (IPATHL.EQ.3) THEN
C
DO 30 I = NL, NLIM, LL
IPL = IPL+LLM1
RADLYR(I) = RADO(IPL)+RADLYR(I)*TRAO(IPL)
TRALYR(I) = TRALYR(I)*TRAO(IPL)
30 CONTINUE
C
ELSEIF (IPATHL.EQ.-1) THEN
C
VI = V1P-DVP
DVI = DVP*FLOAT(LL)
VIDVRF = VI
RFLAST = -1.
NLIM1 = 0
NLIM2 = NL-1
C
40 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
REFLCT = REFLFN(VI,DVI,VIDVRF,RFDEL,RFLAST)
C
IF (VIDVRF.GE.9.E+4) THEN
NLIM2 = NLIM+1
ELSE
NLIM2 = (VIDVRF+DVI-DVP-V1P)/DVP+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIM)
C
C Test to make sure LL divides evenly into (NLIM2-NLIM1+1)
C
NRMNDR = MOD(NLIM2-NLIM1,LL)
IF ((NRMNDR.NE.0).AND.(NLIM2+(LL-NRMNDR).LT.2400)) THEN
NLIM2 = NLIM2+(LL-NRMNDR)
ENDIF
C
DO 50 I = NLIM1, NLIM2, LL
IPL = IPL+LLM1
TRTEMP = TRALYR(I)*TRAO(IPL)
RADLYR(I) = TRALYR(I)*RADO(IPL)+RADLYR(I)+
* RADLYB(I)*TRAO(IPL)*TRTEMP*REFLCT
TRALYR(I) = TRTEMP
C
C Increment interpolation values
C
REFLCT = REFLCT+RFDEL
ILAST = I
50 CONTINUE
C
IF (NLIM2.LT.NLIM) THEN
NLIM2 = ILAST + LL-1
GO TO 40
ENDIF
C
ENDIF
C
C NOT EXACT FREQUENCY - INTERPOLATE RESULT
C
ELSE
C
A1N = A1(NL)
A2N = A2(NL)
A3N = A3(NL)
A4N = A4(NL)
C
IF (IPATHL.EQ.1) THEN
DO 60 I = NL, NLIM, LL
IPL = IPL+LLM1
C
C INTERPOLATE THE OLD RADIANCE
C
RADLYR(I) = TRALYR(I)*(A1N*RADO(IPL-1)+A2N*RADO(IPL)+
* A3N*RADO(IPL+1)+A4N*RADO(IPL+2))+RADLYR(I)
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRALYR(I) = TRALYR(I)*(A1N*TRAO(IPL-1)+A2N*TRAO(IPL)+
* A3N*TRAO(IPL+1)+A4N*TRAO(IPL+2))
60 CONTINUE
C
ELSEIF (IPATHL.EQ.2) THEN
C
DO 70 I = NL, NLIM, LL
IPL = IPL+LLM1
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRTEMP = TRALYR(I)*(A1N*TRAO(IPL-1)+A2N*TRAO(IPL)+
* A3N*TRAO(IPL+1)+A4N*TRAO(IPL+2))
C
C INTERPOLATE THE OLD RADIANCE
C
RADLYR(I) = TRALYR(I)*(A1N*RADO(IPL-1)+A2N*RADO(IPL)+
* A3N*RADO(IPL+1)+A4N*RADO(IPL+2))+
* RADLYR(I)+RADLYB(I)*TRTEMP
TRALYR(I) = TRALYR(I)*TRTEMP
70 CONTINUE
C
ELSEIF (IPATHL.EQ.3) THEN
C
DO 80 I = NL, NLIM, LL
IPL = IPL+LLM1
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRAOI = A1N*TRAO(IPL-1)+A2N*TRAO(IPL)+
* A3N*TRAO(IPL+1)+A4N*TRAO(IPL+2)
C
C INTERPOLATE THE OLD RADIANCE
C
RADLYR(I) = A1N*RADO(IPL-1)+A2N*RADO(IPL)+
* A3N*RADO(IPL+1)+A4N*RADO(IPL+2)+
* RADLYR(I)*TRAOI
TRALYR(I) = TRALYR(I)*TRAOI
80 CONTINUE
C
ELSEIF (IPATHL.EQ.-1) THEN
C
VI = V1P-DVP
DVI = DVP*FLOAT(LL)
VIDVRF = VI
RFLAST = -1
NLIM1 = 0
NLIM2 = NL-1
C
90 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
REFLCT = REFLFN(VI,DVI,VIDVRF,RFDEL,RFLAST)
C
IF (VIDVRF.GE.9.E+4) THEN
NLIM2 = NLIM+1
ELSE
NLIM2 = (VIDVRF+DVI-DVP-V1P)/DVP+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIM)
C
C Test to make sure LL divides evenly into (NLIM2-NLIM1+1)
C
NRMNDR = MOD(NLIM2-NLIM1,LL)
IF ((NRMNDR.NE.0).AND.(NLIM2+(LL-NRMNDR).LT.2400)) THEN
NLIM2 = NLIM2+(LL-NRMNDR)
ENDIF
DO 100 I = NLIM1, NLIM2, LL
IPL = IPL+LLM1
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRAOI = A1N*TRAO(IPL-1)+A2N*TRAO(IPL)+
* A3N*TRAO(IPL+1)+A4N*TRAO(IPL+2)
C
TRTEMP = TRALYR(I)*TRAOI
C
C INTERPOLATE THE OLD RADIANCE
C
RADLYR(I) = TRALYR(I)*(A1N*RADO(IPL-1)+A2N*RADO(IPL)+
* A3N*RADO(IPL+1)+A4N*RADO(IPL+2))+
* RADLYR(I)+RADLYB(I)*TRAOI*TRTEMP*REFLCT
TRALYR(I) = TRTEMP
C
C Increment interpolation values
C
REFLCT = REFLCT+RFDEL
ILAST = I
100 CONTINUE
C
IF (NLIM2.LT.NLIM) THEN
NLIM2 = ILAST + LL -1
GO TO 90
ENDIF
C
ENDIF
C
ENDIF
C
110 CONTINUE
C
NPL = IPL
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE RADINT (NPTS,LFILE,MFILE,JPATHL,TBND) 1,30
C
IMPLICIT REAL*8 (V)
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 5 APRIL 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON RADN(2410),TRAN(2410),RADLYR(-1:4818)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYDUM,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /XPANEL/ V1P,V2P,DVP,NLIM,RMIN,RMAX,NPNLXP,NSHIFT,NPTSS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /XMI/ TRALYR(-1:4818)
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
DIMENSION A1(0:100),A2(0:100),A3(0:100),A4(0:100)
DIMENSION RADO(2),TRAO(2)
DIMENSION WKSAV(35)
C
DIMENSION RADLYB(-1:4818)
C
CHARACTER*40 CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO)
EQUIVALENCE (RADN(1),RADO(1)) , (TRAN(1),TRAO(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
C ************************************************************
C ****** THIS SUBROUTINE DOES LAYER MERGE FOR EMISSION ******
C ****** USING FOUR POINT GENERAL INTERPOLATION ******
C ************************************************************
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
TIMEM = 0.0
TIMRD = 0.0
TIMTB = 0.0
TIMOT = 0.0
C
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
DVL = DV
LAY1SV = LAYR1
PL = PAVE
TL = TAVE
WTOTL = 0.
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
XKT = TAVE/RADCN2
XKTA = TZU/RADCN2
XKTB = 0.
DVK = DV
LAYR1 = LAY1SV
FACT = 1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
ATYPE = 9.999E09
IF (DVK.EQ.DVL) ATYPE = 0.
IF (DVL.GT.DVK) ATYPE = DVK/(DVL-DVK)+0.5
IF (DVL.LT.DVK) ATYPE = -DVL/(DVK-DVL)-0.5
C
C IF (ATYPE .GT. 0) STOP ' RADINT; ATYPE GT 0 '
C
WTOTK = 0.
WRITE (IPR,905) LAYR1,LAYER,KFILE,LFILE,MFILE,ATYPE
IEMIT = 1
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+WK(MOL)*FACT
WK(MOL) = WK(MOL)*FACT+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+WBROAD*FACT
WBROAD = WBROAD*FACT+WN2SAV
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
SECANT = 0.
DV = DVL
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
C
IF (ATYPE.EQ.0.) THEN
C
C 1/1 RATIO ONLY
C
30 CONTINUE
CALL CPUTIM (TIMEM1)
CALL EMIN (V1P,V2P,DVP,NLIM,KFILE,RADLYR(1),RADLYB(1),
* TRALYR(1),KEOF,NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 110
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
CALL BUFIN (LFILE,LEOF,RADO(1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAO(1),NLIMO)
CALL CPUTIM (TIMEM3)
TIMRD = TIMRD+TIMEM3-TIMEM2
DO 40 I = 1, NLIM
RADN(I) = RADO(I)+RADLYR(I)*TRAO(I)
TRAN(I) = TRALYR(I)*TRAO(I)
40 CONTINUE
CALL CPUTIM (TIMEM1)
IF (TBND.GT.0.)
* CALL EMBND (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,TBND)
CALL CPUTIM (TIMEM2)
TIMTB = TIMTB+TIMEM2-TIMEM1
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,MFILE,NPTS,NPANLS)
CALL CPUTIM (TIMEM3)
TIMOT = TIMOT+TIMEM3-TIMEM2
GO TO 30
C
ENDIF
C
C ALL RATIOS EXCEPT 1/1
C
DO 50 JP = 0, 100
APG = JP
P = 0.01*APG
C
C THE FOLLOW ARE THE CONSTANTS FOR THE LAGRANGE 4 POINT
C INTERPOLATION
C
A1(JP) = -P*(P-1.0)*(P-2.0)/6.0
A2(JP) = (P**2-1.0)*(P-2.0)*0.5
A3(JP) = -P*(P+1.0)*(P-2.0)*0.5
A4(JP) = P*(P**2-1.0)/6.0
50 CONTINUE
C
C *** BEGINNING OF LOOP THAT DOES MERGE ***
C
NPE = 0
RADLYR(0) = 0.0
TRALYR(0) = 0.0
V1P = 0.0
V2P = 0.0
DVP = 0.0
V1PO = 0.0
V2PO = 0.0
DVPO = 0.0
KEOF = 1
C
60 CONTINUE
C
CALL CPUTIM (TIMEM1)
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 110
CALL BUFIN (LFILE,LEOF,RADO(1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAO(1),NLIMO)
CALL CPUTIM (TIMEM2)
TIMRD = TIMRD+TIMEM2-TIMEM1
II = 1
C
IF (V2P.LE.V2PO+DVP .AND. KEOF.GT.0) THEN
70 CALL CPUTIM (TIMEM2)
CALL EMIN (V1P,V2P,DVP,NLIM,KFILE,RADLYR(NPE+1),RADLYB(NPE+1),
* TRALYR(NPE+1),KEOF,NPANLS)
CALL CPUTIM (TIMEM3)
TIMEM = TIMEM+TIMEM3-TIMEM2
IF (KEOF.LE.0) GO TO 80
V1P = V1P-FLOAT(NPE)*DVP
NPE = NLIM+NPE
IF (V2P.LE.V2PO+DVP) GO TO 70
ENDIF
C
C ZERO POINT OF FIRST PANEL
C
80 IF (RADLYR(0).EQ.0.0.AND.TRALYR(0).EQ.0.0) THEN
TRALYR(-1) = TRALYR(1)
RADLYR(-1) = RADLYR(1)
RADLYB(-1) = RADLYB(1)
TRALYR(0) = TRALYR(1)
RADLYR(0) = RADLYR(1)
RADLYB(0) = RADLYB(1)
ENDIF
C
C END POINT OF LAST PANEL
C
IF (V2P+DVP.GE.V2) THEN
TRALYR(NPE+1) = TRALYR(NPE)
RADLYR(NPE+1) = RADLYR(NPE)
RADLYB(NPE+1) = RADLYB(NPE)
TRALYR(NPE+2) = TRALYR(NPE)
RADLYR(NPE+2) = RADLYR(NPE)
RADLYB(NPE+2) = RADLYB(NPE)
ENDIF
C
C NPL IS LOCATION OF FIRST ELEMENT ON ARRAYS RADO AND TRAO
C
NPL = 1
C
RATDV = DVL/DVK
C
C FJJ IS OFFSET BY 2. FOR ROUNDING PURPOSES
C
FJ1DIF = (V1PO-V1P)/DVP+1.+2.
C
C ***** BEGINNING OF LOOP THAT DOES MERGE *****
C
DO 90 II = 1, NLIMO
C
FJJ = FJ1DIF+RATDV*FLOAT(II-1)
JJ = IFIX(FJJ)-2
C
JP = (FJJ-FLOAT(JJ))*100.-199.5
C
C INTERPOLATE THE OLD EMISSION
C
RADN(II) = RADO(II)+(A1(JP)*RADLYR(JJ-1)+A2(JP)*RADLYR(JJ)+
* A3(JP)*RADLYR(JJ+1)+A4(JP)*RADLYR(JJ+2))*TRAO(II)
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRAN(II) = (A1(JP)*TRALYR(JJ-1)+A2(JP)*TRALYR(JJ)+
* A3(JP)*TRALYR(JJ+1)+A4(JP)*TRALYR(JJ+2))*TRAO(II)
C
90 CONTINUE
C
NPL = JJ-1
C
CALL CPUTIM (TIMEM1)
IF (TBND.GT.0.) CALL EMBND (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,TBND)
CALL CPUTIM (TIMEM2)
TIMTB = TIMTB+TIMEM2-TIMEM1
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,MFILE,NPTS,NPANLS)
CALL CPUTIM (TIMEM3)
TIMOT = TIMOT+TIMEM3-TIMEM2
C
C NPL IS NOW LOCATION OF FIRST ELEMENT TO BE USED FOR NEXT PASS
C
IPL = -2
DO 100 NL = NPL, NPE
IPL = IPL+1
TRALYR(IPL) = TRALYR(NL)
RADLYR(IPL) = RADLYR(NL)
RADLYB(IPL) = RADLYB(NL)
100 CONTINUE
C
V1P = V1P+FLOAT(NPL+1)*DVP
NPE = IPL
C
GO TO 60
110 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,910) TIME1,TIM,TIMEM,TIMRD,TIMTB,TIMOT
C
RETURN
C
900 FORMAT ('0 THE TIME AT THE START OF RADINT IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5,' WITH XTYPE=',
* G15.5)
910 FORMAT ('0 THE TIME AT THE END OF RADINT IS ',F12.3/F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE - EMIN - ',F12.3,
* ' - READ - ',F12.3,' - EMBND - ',F12.3,' - EMOUT - ',
* F12.3)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE EMBND (V1PO,V2PO,DVPO,NLIMO,NEWEM,NEWTR,TBND) 7,6
C
IMPLICIT REAL*8 (V)
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 9 APRIL 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
DIMENSION NEWEM(*),NEWTR(*)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
REAL NEWEM,NEWTR
C
XKTBND = TBND/RADCN2
VI = V1PO-DVPO
VIDVBD = VI
VIDVEM = VI
BBLAST = -1.
EMLAST = -1.
NLIM1 = 0
NLIM2 = 0
EMDUM = 0.
BBDUM = 0.
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMDUM)
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBDUM)
IEMBB = 0
IF (VIDVBD.GT.VIDVEM) IEMBB = 1
C
10 NLIM1 = NLIM2+1
C
VI = V1PO+FLOAT(NLIM1-1)*DVPO
IF (IEMBB.EQ.0) THEN
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDV,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
ELSE
EMISIV = EMISFN(VI,DVPO,VIDV,EMDEL,EMLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
ENDIF
C
IF (VIDV.GE.9.E+4) THEN
NLIM2 = NLIMO+1
ELSE
NLIM2 = (VIDV-V1PO)/DVPO+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIMO)
C
DO 20 J = NLIM1, NLIM2
NEWEM(J) = NEWEM(J)+NEWTR(J)*EMISIV*BB
C
C Increment interpolation values
C
BB = BB+BBDEL
EMISIV = EMISIV+EMDEL
20 CONTINUE
C
IF (NLIM2.LT.NLIMO) GO TO 10
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE EMOUT (V1P,V2P,DVP,NLIM,NEWEM,NEWTR,MFILE,NPTS,NPANLS) 12,3
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE EMOUT OUTPUTS MERGED EMISSION AND TRANSMITTANCE RESULT
C TO MFILE
C
COMMON /BUFPNL/ V1PBF,V2PBF,DVPBF,NLIMBF
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
DIMENSION PNLHDR(2)
DIMENSION NEWEM(*),NEWTR(*)
C
EQUIVALENCE (PNLHDR(1),V1PBF)
C
REAL NEWEM,NEWTR
C
NPANLS = NPANLS+1
V1PBF = V1P
V2PBF = V2P
DVPBF = DVP
NLIMBF = NLIM
C
CALL BUFOUT (MFILE,PNLHDR(1),NPHDRF)
CALL BUFOUT (MFILE,NEWEM(1),NLIMBF)
CALL BUFOUT (MFILE,NEWTR(1),NLIMBF)
C
IF (NPTS.GT.0) THEN
IF (NPANLS.EQ.1) WRITE (IPR,900)
WRITE (IPR,905)
NNPTS = NPTS
IF (NPTS.GT.(NLIMBF/2)+1) NNPTS = (NLIMBF/2)+1
JEND = NLIMBF-NNPTS+1
DO 10 J = 1, NNPTS
VJ = V1PBF+FLOAT(J-1)*DVPBF
K = J+JEND-1
VK = V1PBF+FLOAT(K-1)*DVPBF
WRITE (IPR,910) J,VJ,NEWEM(J),NEWTR(J),
* K,VK,NEWEM(K),NEWTR(K)
10 CONTINUE
ENDIF
C
RETURN
C
900 FORMAT ('0 ','LOCATION WAVENUMBER',2X,'RADIANCE',7X,
* 'TRANSMITTANCE',22X,'LOCATION WAVENUMBER',2X,
* 'RADIANCE',7X,'TRANSMITTANCE')
905 FORMAT (' ')
910 FORMAT (I8,2X,F12.6,1P2E15.7,20X,I8,2X,0PF12.6,1P2E15.7)
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE OPNODF(NLAYER,LAYER,PTHODL,HFMODL,IEMIT) 4
C
C This subroutine opens file for calculating the radiance using
C precalculated optical depths
C (IEMIT = 1,IMRG=A/12,B/22,C/32,40,41)
C
LOGICAL OP
CHARACTER*57 FILE1
CHARACTER*55 PTHODL
CHARACTER*11 CFORM
CHARACTER*10 HFMODL
C
C -------------------------
C Common block for analytic derivative
C -------------------------
COMMON /ADRFIL/ KODFIL,KODTOT,KTEMP,KFILAD
C -------------------------
C
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
C 123456789-123456789-123456789-123456789-123456789-1234567
DATA FILE1 /
* ' '/
DATA CFORM / 'UNFORMATTED' /
C
C For calculation of analytic derivative, open previously
C calculated optical depth files as unit KODFIL. Otherwise,
C open as KFILE.
C
IF (IEMIT.EQ.3) THEN
KOPEN = KODFIL
ELSE
KOPEN = KFILE
ENDIF
C
WRITE(IPR,910) LAYER,NLAYER
INQUIRE (UNIT=KOPEN,OPENED=OP)
IF (OP) CLOSE (KOPEN)
WRITE(FILE1,HFMODL) PTHODL,LAYER
OPEN(UNIT=KOPEN,FILE=FILE1,FORM=CFORM,STATUS='OLD')
C
C Write procedure
C
WRITE(IPR,900) FILE1
C
RETURN
C
900 FORMAT (' Opened layer optical depth file: ',A57)
910 FORMAT ('LAYER ',I5,' OF ',I5,':')
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE OPNRAD(NLAYER,LAYER) 2
C
C This subroutine opens file for calculating the layer radiances
C (IEMIT = 1; IMRG=45,46)
C
LOGICAL OP
CHARACTER*57 FILE1
CHARACTER*11 CFORM
CHARACTER*55 PTHRAD
CHARACTER*10 HFMRAD
C
C Common block for layer radiances
C -------------------------
COMMON /RADLAY/ PTHRAD,HFMRAD
C -------------------------
C
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
C 123456789-123456789-123456789-123456789-123456789-1234567
DATA FILE1 /
* ' '/
DATA CFORM / 'UNFORMATTED' /
C
INQUIRE (UNIT=NFILE,OPENED=OP)
IF (OP) CLOSE (NFILE)
WRITE(FILE1,HFMRAD) PTHRAD,LAYER
OPEN(UNIT=NFILE,FILE=FILE1,FORM=CFORM,STATUS='UNKNOWN')
C
C Write procedure
C
WRITE(IPR,900) FILE1
C
RETURN
C
900 FORMAT (' Opened layer radiance file: ',A57,/)
910 FORMAT ('LAYER ',I5,' OF ',I5,':')
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE OPNDRV(NLAYER,LAYER,LAYTOT) 2
C
C This subroutine opens file for calculating the layer derivatives
C (IEMIT = 3)
C
LOGICAL OP
CHARACTER*57 FILE1,FILE2,FILE3
CHARACTER*11 CFORM
CHARACTER*55 CDUM1,PTHODI,PTHODT,PTHRDR
CHARACTER*10 HFMODI,HFMODT,HFMRDR
C
C Common block for analytic derivatives
C -------------------------
COMMON /ADRPNM/ CDUM1,PTHODI,PTHODT,PTHRDR
COMMON /ADRFRM/ HFMODI,HFMODT,HFMRDR
COMMON /ADRFIL/ KODFIL,KODTOT,KTEMP,KFILAD
C -------------------------
C
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
C
C 123456789-123456789-123456789-123456789-123456789-1234567
DATA FILE1 /
* ' '/,
* FILE2 /
* ' '/,
* FILE3 /
* ' '/
DATA CFORM / 'UNFORMATTED' /
C
WRITE(IPR,910) LAYER,NLAYER
INQUIRE (UNIT=KODFIL,OPENED=OP)
IF (OP) CLOSE (KODFIL)
WRITE(FILE1,HFMODI) PTHODI,LAYER
OPEN(UNIT=KODFIL,FILE=FILE1,FORM=CFORM,STATUS='OLD')
C
IF (LAYER.NE.NLAYER) THEN
INQUIRE (UNIT=KODTOT,OPENED=OP)
IF (OP) CLOSE (KODTOT)
WRITE(FILE2,HFMODT) PTHODT,LAYTOT
OPEN(UNIT=KODTOT,FILE=FILE2,FORM=CFORM,STATUS='OLD')
ELSE
FILE2 = 'NOT USED'
ENDIF
C
INQUIRE (UNIT=KTEMP,OPENED=OP)
IF (OP) CLOSE (KTEMP)
OPEN(UNIT=KTEMP,FILE='mono_ad.TMP',FORM=CFORM,STATUS='unknown')
C
INQUIRE (UNIT=KFILAD,OPENED=OP)
IF (OP) CLOSE (KFILAD)
WRITE(FILE3,HFMRDR) PTHRDR,LAYER
OPEN(UNIT=KFILAD,FILE=FILE3,FORM=CFORM,STATUS='UNKNOWN')
C
C Write procedure
C
WRITE(IPR,900) FILE1,FILE2,FILE3
C
RETURN
C
900 FORMAT (' Opened layer optical depth file: ',A57,/,
* ' Opened accumulated optical depth file: ',A57,/,
* ' Opened layer analytic derivative file: ',A57)
910 FORMAT ('LAYER ',I5,' OF ',I5,':')
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE EMDM (V1P,V2P,DVP,NLIM,KFILE,EM,EMB,TR,KEOF,NPANLS) 2,70
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE EMDM INPUTS OPTICAL DEPTH VALUES FROM KFILE AND
C CALCULATES SOURCE FUNCTION FOR THE LAYER.
C THIS VERSION WORKS FOR AEROSOLS AND NLTE.
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KDUMY,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /LAMCHN/ ONEPL,ONEMI,EXPMIN,ARGMIN
COMMON /BUFPNL/ V1PBF,V2PBF,DVPBF,NLIMBF
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
COMMON /EMDMSV/ BBSAV(2410),BBASAV(2410),XXSAV(2410)
C
DIMENSION PNLHDR(2),EM(*),EMB(*),TR(*)
C
EQUIVALENCE (FSCDID(1),IHIRAC) , (FSCDID(2),ILBLF4)
EQUIVALENCE (PNLHDR(1),V1PBF)
EQUIVALENCE (FSCDID(4),IAERSL)
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) RETURN
CALL BUFIN (KFILE,KEOF,TR(1),NLIMBF)
C
C TR CONTAINS THE OPTICAL DEPTHS AT THIS STAGE
C
IF (IHIRAC.EQ.4) CALL BUFIN (KFILE,KEOF,EM(1),NLIMBF)
C
C EM CONTAINS THE OPTICAL DEPTH CORRECTIONS FOR NLTE AT THIS STAGE
C
IF (NPANLS.LT.1.AND.IAERSL.EQ.0) WRITE (IPR,900)
IF (NPANLS.LT.1.AND.IAERSL.NE.0) WRITE (IPR,905)
C
EXT = 0.
ADEL = 0.
RADFN0 = 0.
RDEL = 0.
BB = 0.
BBDEL = 0.
BBA = 0.
BBDLA = 0.
BBB = 0.
BBDLB = 0.
C
V1P = V1PBF
V2P = V2PBF
DVP = DVPBF
NLIM = NLIMBF
VI = V1P-DVP
VIDV = VI
VIBB = VI
VAER = VI
VDUM = VI
BBLAST = -1.
BBLXTA = -2.
BBLXTB = -3.
RDLAST = -1.
BBDUM = -4.
RDDUM = -1.
NLIM1 = 0
NLIM2 = 0
C
AA = 0.278
C
IF (IAERSL.NE.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBDUM)
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDDUM)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
IAFBB = 0
IF (VITST.LT.VAER.AND.VITST.LT.VIBB) IAFBB = 1
IF (VAER.LT.VITST.AND.VAER.LT.VIBB) IAFBB = 2
ELSE
IAFBB = -1
ENDIF
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER ONLY ONCE
C
C - WITH XKTA=0 THIS ALGORITHM REVERTS TO THE ORIGINAL
C
IF (XKTB.LE.0.) THEN
C
C - THIS SECTION TREATS THE LTE CASE
C
IF (IHIRAC.NE.4) THEN
C
10 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 20 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(BB+XX*BBA)/(1.+XX)
C
C Store BB, BBA, and XX for derivative source term
C
XXSAV(I) = XX
BBSAV(I) = BB
BBASAV(I) = BBA
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
20 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 10
ELSE
C
C - THIS SECTION TREATS THE NLTE CASE
C
30 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 40 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(1.0-EM(I)/ODVI)*(BB+XX*BBA)/(1.+XX)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
40 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 30
C
ENDIF
ELSE
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER TWICE:
C
C - FOR TANGENT PATHS AND FOR THE CASE OF THE REFLECTED ATMOSPHERE
C
IF (IHIRAC.NE.4) THEN
C
C - THIS SECTION TREATS THE LTE CASE
C
50 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 60 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EMX = (1.-TR(I))/(1.+XX)
EM(I) = EMX*(BB+XX*BBA)
EMB(I) = EMX*(BB+XX*BBB)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBB = BBB+BBDLB
60 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 50
C
ELSE
C
C - THIS SECTION TREATS THE CASE OF NLTE
C
70 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 80 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EMX = (1.-TR(I))*(1.0-EM(I)/ODVI)/(1.+XX)
EM(I) = EMX*(BB+XX*BBA)
EMB(I) = EMX*(BB+XX*BBB)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBB = BBB+BBDLB
80 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 70
C
ENDIF
ENDIF
C
RETURN
C
900 FORMAT ('0EMISSION AND TRANSMISSION (MOLECULAR) ')
905 FORMAT ('0EMISSION AND TRANSMISSION (AEROSOLS EFFECTS INCLUDED)')
C
END
C
C ---------------------------------------------------------------
SUBROUTINE EMDT (V1P,V2P,DVP,NLIM,KFILE,EM,EMB,TR,KEOF,NPANLS) 1,74
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE EMDT inputs optical depth values from kfile and
C calculates source function for the layer.
C
C This version is used for analytic temperature derivatives.
C Non-LTE and limb not yet implemented.
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C IMPLEMENTATION: P.D. Brown
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KDUMY,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /LAMCHN/ ONEPL,ONEMI,EXPMIN,ARGMIN
COMMON /BUFPNL/ V1PBF,V2PBF,DVPBF,NLIMBF
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
COMMON /EMDMSV/ BBSAV(2410),BBASAV(2410),XXSAV(2410)
COMMON /EMDTSV/ BBDSAV(2410)
C
DIMENSION PNLHDR(2),EM(*),EMB(*),TR(*)
C
EQUIVALENCE (FSCDID(1),IHIRAC) , (FSCDID(2),ILBLF4)
EQUIVALENCE (PNLHDR(1),V1PBF)
EQUIVALENCE (FSCDID(4),IAERSL)
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) RETURN
CALL BUFIN (KFILE,KEOF,TR(1),NLIMBF)
C
C TR contains the optical depths at this stage
C
IF (IHIRAC.EQ.4) CALL BUFIN (KFILE,KEOF,EM(1),NLIMBF)
C
C EM contains the optical depth corrections for nlte at this stage
C
IF (NPANLS.LT.1.AND.IAERSL.EQ.0) WRITE (IPR,900)
IF (NPANLS.LT.1.AND.IAERSL.NE.0) WRITE (IPR,905)
C
EXT = 0.
ADEL = 0.
RADFN0 = 0.
RDEL = 0.
BB = 0.
BBDEL = 0.
BBA = 0.
BBDLA = 0.
BBB = 0.
BBDLB = 0.
BBD = 0.
BBADDL = 0.
C
V1P = V1PBF
V2P = V2PBF
DVP = DVPBF
NLIM = NLIMBF
VI = V1P-DVP
VIDV = VI
VIBB = VI
VAER = VI
VDUM = VI
BBLAST = -1.
BBLXTA = -2.
BBLXTB = -3.
RDLAST = -1.
BBDUM = -4.
RDDUM = -1.
NLIM1 = 0
NLIM2 = 0
C
AA = 0.278
C
IF (IAERSL.NE.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBDUM)
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDDUM)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
IAFBB = 0
IF (VITST.LT.VAER.AND.VITST.LT.VIBB) IAFBB = 1
IF (VAER.LT.VITST.AND.VAER.LT.VIBB) IAFBB = 2
ELSE
IAFBB = -1
ENDIF
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER ONLY ONCE
C
C - WITH XKTA=0 THIS ALGORITHM REVERTS TO THE ORIGINAL
C
IF (XKTB.LE.0.) THEN
C
C - THIS SECTION TREATS THE LTE CASE
C
IF (IHIRAC.NE.4) THEN
C
10 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
BBD = BBAD(BB,VI,DVP,V2P,XKT,VIDD,BBADDL,BBADOL)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
BBD = BBAD(BB,VI,DVP,V2P,XKT,VIDD,BBADDL,BBADOL)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
BBD = BBAD(BB,VI,DVP,V2P,XKT,VIDD,BBADDL,BBADOL)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
BBD = BBAD(BB,VI,DVP,V2P,XKT,VIDD,BBADDL,BBADOL)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 20 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
C
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(BB+XX*BBA)/(1.+XX)
C
C Store BB, BBA, BBD, and XX for derivative source term
C
XXSAV(I) = XX
BBSAV(I) = BB
BBASAV(I) = BBA
BBDSAV(I) = BBD
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBD = BBD+BBADDL
20 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 10
ELSE
C
C - THIS SECTION TREATS THE NLTE CASE
C
30 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 40 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
C
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(1.0-EM(I)/ODVI)*(BB+XX*BBA)/(1.+XX)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
40 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 30
C
ENDIF
ELSE
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER TWICE:
C
C - FOR TANGENT PATHS AND FOR THE CASE OF THE REFLECTED ATMOSPHERE
C
IF (IHIRAC.NE.4) THEN
C
C - THIS SECTION TREATS THE LTE CASE
C
50 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 60 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EMX = (1.-TR(I))/(1.+XX)
EM(I) = EMX*(BB+XX*BBA)
EMB(I) = EMX*(BB+XX*BBB)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBB = BBB+BBDLB
60 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 50
C
ELSE
C
C - THIS SECTION TREATS THE CASE OF NLTE
C
70 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
IF (XKTB.GT.0.) THEN
VIBB = -VIDV
BBB = BBFN(VI,DVP,V2P,XKTB,VIBB,BBDLB,BBLXTB)
ELSE
BBB = BB
BBDLB = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 80 I = NLIM1, NLIM2
ODVI = TR(I)+EXT*RADFN0
C
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EMX = (1.-TR(I))*(1.0-EM(I)/ODVI)/(1.+XX)
EM(I) = EMX*(BB+XX*BBA)
EMB(I) = EMX*(BB+XX*BBB)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
BBB = BBB+BBDLB
80 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 70
C
ENDIF
ENDIF
C
RETURN
C
900 FORMAT ('0EMISSION AND TRANSMISSION (MOLECULAR) ')
905 FORMAT ('0EMISSION AND TRANSMISSION (AEROSOLS EFFECTS INCLUDED)')
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE EMADL1 (NPTS,MFILE,JPATHL,TBND) 1,33
C
IMPLICIT REAL*8 (V)
C
C Calculates radiance and radiance derivative for first layer
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C IMPLEMENTATION: P.D. Brown
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON NEWEM(2410),NEWTR(2410)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYER,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMIHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYDUM,YI1,YID(10),LSTWDF
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILA,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
COMMON /ADRLYR/ KSUBL(0:5000),OPDT(0:5000)
COMMON /ADRFIL/ KODFIL,KODTOT,KTEMP,KFILAD
COMMON /IADFLG/ IANDER,NSPCRT
C
CHARACTER*40 CEXT,CYID
C
DIMENSION EMLAYB(2410)
DIMENSION XFILHD(2),OPNLHD(2),XFHDUM(2)
DIMENSION EMLAYR(2),TRLAYR(2)
DIMENSION XKMOLC(2),ODACUM(2)
DIMENSION RPRIME(2410),OLDEM(0:5000)
C
EQUIVALENCE (XFILHD(1),XID(1)) , (OPNLHD(1),V1PO)
EQUIVALENCE (NEWEM(1),EMLAYR(1)) , (NEWTR(1),TRLAYR(1)),
* (KSUBL(1),XKMOLC(1)) , (OPDT(1),ODACUM(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
REAL NEWEM,NEWTR
C
DATA NDIM / 2410 /,ND2 / 5000 /
C
C TBND IS THE BOUNDARY BLACK BODY TEMPERATUE
C
C IPATHL =-1 IS FOR THE LOOKING DOWN CASE WITH REFLECTED ATMOSPHERE
C IPATHL = 0 IS FOR THE HORIZONTAL PATH CASE (HOMOGENEOUS LAYER)
C IPATHL = 1 IS FOR THE LOOKING DOWN CASE (TO DENSER LAYERS)
C IPATHL = 2 IS FOR THE SYMMETRIC TANGENT PATH CASE
C IPATHL = 3 IS FOR THE LOOKING UP CASE (TO LESS DENSE LAYERS
C
CALL CPUTIM (TIME)
C
C ** NOTE ON IPATHL =2
C THE TANGENT MERGE ROUTINES ARE DIVIDED INTO ANTERIOR (1ST)
C AND POSTERIOR (2ND) LAYER CROSSINGS. THIS RECOGNITION IS
C TRIGGERED BY THE VALUE OF "IANT"
C
C IF IANT.EQ. 1 THEN POSTERIOR MERGE
C IF IANT.EQ. 0 THEN NORMAL MERGE
C IF IANT.EQ.-1 THEN ANTERIOR MERGE
C
WRITE (IPR,900) TIME
NPANLS = 0
C
C Read in file headers for layer absorptance coefficients and
C layer optical depths and total optical depths (if there is more th
C one layer between the present layer and the observer)
C
CALL BUFIN (KFILE,KEOF,XFHDUM(1),NFHDRF)
IF (LAYER.LT.NLAYER) CALL BUFIN (KODTOT,KEOF,XFILHD(1),NFHDRF)
CALL BUFIN (KODFIL,KEOF,XFILHD(1),NFHDRF)
C
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
C FOR AEROSOL RUNS, MOVE EXTID INTO YID
C
IF (IAERSL.GT.0) THEN
WRITE (CEXT,'(10A4)') EXTID
WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CYID(19:40) = CEXT(19:40)
READ (CYID,'(5A8)') (YID(I),I=3,7)
ENDIF
C
C IF BOUNDARY PROPERTIES ARE SUPPLIED, AND DOWNWARD LOOKING
C CASE; SET IPATHL TO REFLECTED ATMOSPHERE CASE
C
IF (IBPROP.EQ.1.AND.IPATHL.EQ.1) IPATHL = -1
IEMIT = 1
FACT = 1.
TIMEM = 0.0
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
DO 10 MOL = 1, NMOL
WK(MOL) = WK(MOL)*FACT
10 CONTINUE
WBROAD = WBROAD*FACT
LAYR1 = LAYER
WRITE (IPR,905) LAYR1,LAYER,KFILE,MFILE
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
CALL BUFOUT (KTEMP,XFILHD(1),NFHDRF)
DVXM = DV
XKT = TAVE/RADCN2
XKTBND = TBND/RADCN2
IF (IPATHL.EQ.-1) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ENDIF
IF (IPATHL.EQ.0) THEN
XKTA = 0.
XKTB = 0.
ENDIF
IF (IPATHL.EQ.1) THEN
XKTA = TZU/RADCN2
XKTB = 0.
ENDIF
IF (IPATHL.EQ.2) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ENDIF
IF (IPATHL.EQ.3) THEN
XKTA = TZL/RADCN2
XKTB = 0.
ENDIF
WRITE (IPR,910) IPATHL,IANT
C
20 CONTINUE
C
C
C Input emission and transmission, and calculate layer
C source function
C
C Call EMDT for temperature retrieval
C Call EMDM for molecular retrieval
C
C
CALL CPUTIM (TIMEM1)
IF (NSPCRT.EQ.29) THEN
CALL EMDT (V1PO,V2PO,DVPO,NLIMO,KODFIL,EMLAYR,EMLAYB,
* TRLAYR,KEOF,NPANLS)
ELSEIF (NSPCRT.GT.0) THEN
CALL EMDM (V1PO,V2PO,DVPO,NLIMO,KODFIL,EMLAYR,EMLAYB,
* TRLAYR,KEOF,NPANLS)
ENDIF
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 80
VI = V1PO-DVPO
VIDVBD = VI
VIDVEM = VI
VIDVRF = VI
BBLAST = -1.
EMLAST = -1.
RFLAST = -1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) THEN
DO 30 J = 1, NLIMO
TRJ = TRLAYR(J)
NEWEM(J) = EMLAYR(J)+EMLAYB(J)*TRJ
TRLAYR(J) = TRLAYR(J)*TRJ
30 CONTINUE
ELSEIF (((IPATHL.EQ.1).AND.(TBND.GT.0.)).OR.(IPATHL.EQ.3)) THEN
C
NLIM1 = 0
NLIM2 = 0
EMDUM = 0.
BBDUM = 0.
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMDUM)
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBDUM)
IEMBB = 0
IF (VIDVBD.GT.VIDVEM) IEMBB = 1
C
40 NLIM1 = NLIM2+1
C
VI = V1PO+FLOAT(NLIM1-1)*DVPO
IF (IEMBB.EQ.0) THEN
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDV,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
ELSE
EMISIV = EMISFN(VI,DVPO,VIDV,EMDEL,EMLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
ENDIF
C
IF (VIDV.GE.9.E+4) THEN
NLIM2 = NLIMO+1
ELSE
NLIM2 = (VIDV-V1PO)/DVPO+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIMO)
C
DO 50 J = NLIM1, NLIM2
OLDEM(J) = BB*EMISIV
NEWEM(J) = EMLAYR(J)+TRLAYR(J)*OLDEM(J)
C
C Increment interpolation values
C
EMISIV = EMISIV+EMDEL
BB = BB+BBDEL
50 CONTINUE
C
IF (NLIM2.LT.NLIMO) GO TO 40
C
ELSEIF ((IPATHL.EQ.-1).AND.(TBND.GT.0.)) THEN
C
NLIM1 = 0
NLIM2 = 0
EMDUM = 0.
RFDUM = 0.
BBDUM = 0.
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMDUM)
REFLCT = REFLFN(VI,DVPO,VIDVRF,RFDEL,RFDUM)
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBDUM)
IEMBB = 0
IF (VIDVEM.LT.VIDVRF.AND.VIDVEM.LT.VIDVBD) IEMBB = 1
IF (VIDVRF.LT.VIDVEM.AND.VIDVRF.LT.VIDVBD) IEMBB = 2
C
60 NLIM1 = NLIM2+1
C
VI = V1PO+FLOAT(NLIM1-1)*DVPO
IF (IEMBB.EQ.0) THEN
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDV,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
VIDVRF = -VIDV
REFLCT = REFLFN(VI,DVPO,VIDVRF,RFDEL,RFLAST)
ELSEIF (IEMBB.EQ.1) THEN
EMISIV = EMISFN(VI,DVPO,VIDV,EMDEL,EMLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
VIDVRF = -VIDV
REFLCT = REFLFN(VI,DVPO,VIDVRF,RFDEL,RFLAST)
ELSE
REFLCT = REFLFN(VI,DVPO,VIDV,RFDEL,RFLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
ENDIF
C
IF (VIDV.GE.9.E+4) THEN
NLIM2 = NLIMO+1
ELSE
NLIM2 = (VIDV-V1PO)/DVPO+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIMO)
C
DO 70 J = NLIM1, NLIM2
NEWEM(J) = EMLAYR(J)+EMLAYB(J)*REFLCT*TRLAYR(J)+
* TRLAYR(J)*BB*EMISIV
C
C Increment interpolation values
C
BB = BB+BBDEL
EMISIV = EMISIV+EMDEL
REFLCT = REFLCT+RFDEL
70 CONTINUE
C
IF (NLIM2.LT.NLIMO) GO TO 60
C
ENDIF
C
C Output radiance to MFILE
C
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,NEWEM,NEWTR,MFILE,NPTS,NPANLS)
C
C Combine terms of analytic layer radiance derivative
C
C -------------------------------
C Analytic Derivative calculation
C -------------------------------
C
C
C Call TDERIV for temperature retrieval
C Call ADERIV for molecular retrieval
C
C If molecular retrieval:
C Input layer absorptance coefficients and accumulated
C transmittances and calculate layer derivatives
C
C If temperature retrieval:
C Input Planck function derivative and layer
C transmittances and calculate layer derivatives
C
IF (NSPCRT.EQ.29) THEN
CALL TDERIV (KFILE,KODTOT,RPRIME,OLDEM,TRLAYR,KSUBL,XKMOLC,
* OPDT,ODACUM,NLIMO,NDIM,ND2,IPATHL,LAYER,NLAYER,
* v1po,dvpo)
ELSEIF (NSPCRT.GT.0) THEN
CALL ADERIV (KFILE,KODTOT,RPRIME,OLDEM,TRLAYR,KSUBL,XKMOLC,
* OPDT,ODACUM,NLIMO,NDIM,ND2,IPATHL,LAYER,NLAYER,
* v1po,dvpo)
ENDIF
C
C Output monochromatic radiance derivative to KTEMP
C
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,RPRIME,NEWTR,KTEMP,NPTS,NPANLS)
C
GO TO 20
80 CALL CPUTIM (TIME1)
TIME = TIME1-TIME
WRITE (IPR,915) TIME,TIMEM
C
RETURN
C
900 FORMAT (' TIME AT THE START OF --EMADL1-- ',F10.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,
* '0 INPUT FILE =',I5,' OUTPUT FILE =',I5)
910 FORMAT ('0 IPATHL AND IANT',2I5)
915 FORMAT (' TIME REQUIRED FOR --EMADL1-- ',F10.3,
* ' --EMDM-- ',F10.3)
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE EMADMG (NPTS,LFILE,MFILE,JPATHL,TBND) 1,24
C
IMPLICIT REAL*8 (V)
C
C Merges layer radiances when calculating layer
C radiance derivatives.
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON RADN(2410),TRAN(2410),RADO(0:5000)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYER,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYDUM,YI1,YID(10),LSTWDF
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /XPANEL/ V1P,V2P,DVP,NLIM,RMIN,RMAX,NPNLXP,NSHIFT,NPTSS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /XME/ TRAO(0:5000)
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
COMMON /ADRLYR/ KSUBL(0:5000),OPDT(0:5000)
COMMON /ADRFIL/ KODFIL,KODTOT,KTEMP,KFILAD
COMMON /IADFLG/ IANDER,NSPCRT
C
DIMENSION RADLYB(2410)
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
DIMENSION A1(10),A2(10),A3(10),A4(10)
DIMENSION RADLYR(2),TRALYR(2),RADOI(2),TRAOI(2)
DIMENSION WKSAV(35)
DIMENSION XKMOLC(2),ODACUM(2)
DIMENSION RPRIME(2410)
C
CHARACTER*40 CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO)
EQUIVALENCE (RADO(1),RADOI(1)) , (TRAO(1),TRAOI(1)),
* (KSUBL(1),XKMOLC(1)) , (OPDT(1),ODACUM(1)),
* (RADN(1),RADLYR(1)) , (TRAN(1),TRALYR(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
DATA A1 /10*0.0/, A2 /10*0.0/, A3 /10*0.0/, A4 /10*0.0/
DATA NDIM / 2410 /,ND2 / 5000 /
C
C
C IPATHL =-1 IS FOR THE LOOKING DOWN CASE FOR REFLECTED ATMOSPHERE
C IPATHL = 1 IS FOR THE LOOKING DOWN CASE (TO DENSER LAYERS)
C IPATHL = 2 IS FOR THE SYMMETRIC TANGENT PATH CASE
C IPATHL = 3 IS FOR THE LOOKING UP CASE (TO LESS DENSE LAYERS)
C
C
C ** NOTE ON IPATHL = 2
C THE TANGENT MERGE ROUTINES ARE DIVIDED INTO ANTERIOR (1ST)
C AND POSTERIOR (2ND) LAYER CROSSINGS THIS RECOGNITION IS
C TRIGGERED BY THE VALUE OF "IANT"
C
C IF IANT.EQ. 1 THEN POSTERIOR MERGE
C IF IANT.EQ. 0 THEN NORMAL MERGE
C IF IANT.EQ.-1 THEN ANTERIOR MERGE
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
TIMEM = 0.0
TIMRD = 0.0
TIMOT = 0.0
C
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
LAY1SV = LAYR1
DVL = DV
PL = PAVE
TL = TAVE
WTOTL = 0.
C
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
C
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
C
C Read in file headers for layer absorptance coefficients, layer
C optical depths, and total optical depths (if there is more than
C one layer between the present layer and the observer)
C
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (LAYER.LT.NLAYER) CALL BUFIN (KODTOT,KEOF,XFILHD(1),NFHDRF)
CALL BUFIN (KODFIL,KEOF,XFILHD(1),NFHDRF)
C
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
C
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
C IF BOUNDARY PROPERTIES ARE SUPPLIED, AND DOWNWARD LOOKING
C CASE; SET IPATHL TO REFLECTED ATMOSPHERE CASE
C
IF (IBPROP.EQ.1.AND.IPATHL.EQ.1) IPATHL = -1
TAVK = TAVE
DVK = DV
FACT = 1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
C
IF (DVL.EQ.DVK) THEN
ITYPE = 0
ELSEIF (DVL.GT.DVK) THEN
ITYPE = DVK/(DVL-DVK)+0.5
ELSE
C
C DVL.LT.DVK
C
ITYPE = -INT(DVL/(DVK-DVL)+0.5)
ENDIF
IF (ITYPE.LT.0) STOP ' EMADMG; ITYPE LT 0 '
C
WTOTK = 0.
LAYR1 = LAY1SV
WRITE (IPR,905) LAYR1,LAYER,KODFIL,LFILE,MFILE
IEMIT = 1
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+WK(MOL)*FACT
WK(MOL) = WK(MOL)*FACT+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+WBROAD*FACT
WBROAD = WBROAD*FACT+WN2SAV
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
SECANT = 0.
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
CALL BUFOUT (KTEMP,XFILHD(1),NFHDRF)
DVXM = DV
XKT = TAVK/RADCN2
C
WRITE (IPR,910) IPATHL,IANT
C
IF (IPATHL.EQ.-1) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ELSEIF (IPATHL.EQ.1) THEN
XKTA = TZU/RADCN2
XKTB = 0.
ELSEIF (IPATHL.EQ.2) THEN
XKTA = TZU/RADCN2
XKTB = TZL/RADCN2
ELSEIF (IPATHL.EQ.3) THEN
XKTA = TZL/RADCN2
XKTB = 0.
ELSE
STOP ' EMADMG; IPATHL '
ENDIF
C
ATYPE = ITYPE
LL = ITYPE+1
AP = 1.0/(ATYPE+1.0)
C
C *** BEGINNING OF LOOP THAT DOES MERGE ***
C
NPE = 0
RADO(0) = 0.0
TRAO(0) = 0.0
V1PO = 0.0
V2PO = 0.0
DVPO = 0.0
C
40 CONTINUE
C
C
C Input emission and transmission, and calculate layer
C source function
C
CALL CPUTIM (TIMEM1)
CALL EMDM (V1P,V2P,DVP,NLIM,KODFIL,RADLYR,RADLYB,TRALYR,KEOF,
* NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 80
C
II = 1
C
IF (V2PO.LE.V2P+DVPO) THEN
50 CALL CPUTIM (TIMEM1)
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 60
CALL BUFIN (LFILE,LEOF,RADOI(NPE+1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAOI(NPE+1),NLIMO)
CALL CPUTIM (TIMEM2)
TIMRD = TIMRD+TIMEM2-TIMEM1
NPE = NLIMO+NPE
IF (V2PO.LE.V2P+DVPO) GO TO 50
ENDIF
C
C ZERO POINT OF FIRST PANEL
C
60 IF (RADO(0).EQ.0.0.AND.TRAO(0).EQ.0.0) THEN
RADO(0) = RADO(1)
TRAO(0) = TRAO(1)
ENDIF
C
C END POINT OF LAST PANEL
C
IF (V2PO+DVPO.GE.V2) THEN
RADO(NPE+1) = RADO(NPE)
RADO(NPE+2) = RADO(NPE)
TRAO(NPE+1) = TRAO(NPE)
TRAO(NPE+2) = TRAO(NPE)
ENDIF
C
C -------------------------------
C Analytic Derivative calculation
C -------------------------------
C
C Call TDERIV for temperature retrieval
C Call ADERIV for molecular retrieval
C
C If molecular retrieval:
C Input layer absorptance coefficients and accumulated
C transmittances and calculate layer derivatives
C
C If temperature retrieval:
C Input Planck function derivative and layer
C transmittances and calculate layer derivatives
C
IF (NSPCRT.EQ.29) THEN
CALL TDERIV (KFILE,KODTOT,RPRIME,RADO,TRALYR,KSUBL,XKMOLC,
* OPDT,ODACUM,NLIM,NDIM,ND2,IPATHL,LAYER,NLAYER,
* v1p,dvp)
ELSEIF (NSPCRT.GT.0) THEN
CALL ADERIV (KFILE,KODTOT,RPRIME,RADO,TRALYR,KSUBL,XKMOLC,
* OPDT,ODACUM,NLIM,NDIM,ND2,IPATHL,LAYER,NLAYER,
* v1p,dvp)
ENDIF
C
C Output monochromatic radiance derivatives to KTEMP
C
CALL EMOUT (V1P,V2P,DVP,NLIM,RPRIME,TRAN,KTEMP,NPTS,NPANLS)
C
C --------------------------
C Layer Radiance Calculation
C --------------------------
C
NPL = 1
C
C NPL IS LOCATION OF FIRST ELEMENT ON ARRAYS RADO AND TRAO
C Combine terms of layer radiative transfer
C
CALL RADNN (RADN,TRAN,RADO,TRAO,RADLYB,NLIM,NDIM,ND2,V1P,DVP,
* IPATHL,A1,A2,A3,A4,LL,NPL)
C
CALL CPUTIM (TIMEM1)
C
IF (TBND.GT.0.) CALL EMBND (V1P,V2P,DVP,NLIM,RADN,TRAN,TBND)
C
C Output radiance to MFILE
C
CALL EMOUT (V1P,V2P,DVP,NLIM,RADN,TRAN,MFILE,NPTS,NPANLS)
C
CALL CPUTIM (TIMEM2)
TIMOT = TIMOT+TIMEM2-TIMEM1
C
C NPL IS NOW LOCATION OF FIRST ELEMENT TO BE USED FOR NEXT PASS
C
IPL = -1
DO 70 NL = NPL, NPE
IPL = IPL+1
RADO(IPL) = RADO(NL)
TRAO(IPL) = TRAO(NL)
70 CONTINUE
C
NPE = IPL
C
GO TO 40
C
C End of loop over panels
C
80 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,915) TIME1,TIM,TIMEM,TIMRD,TIMOT
C
RETURN
C
C
900 FORMAT ('0 THE TIME AT THE START OF EMADMG IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5)
910 FORMAT ('0 IPATHL AND IANT',2I5)
915 FORMAT ('0 THE TIME AT THE END OF EMADMG IS ',F12.3/F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE - EMDM - ',
* F12.3,' - READ - ',F12.3,' - EMOUT - ',F12.3)
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE ADERIV (KFILE,KODTOT,RPRIME,RADO,TRALYR,KSUBL,XKMOLC, 2,4
* OPDT,ODACUM,NLIM,NDIM,ND2,IPATHL,LAYER,NLAYER,
* v1po,dvpo)
C
C This subroutine inputs abosrption coefficient values from
C KFILE and total transmittance from KODTOT (if there is more
C than one layer between the present layer and the observer),
C and then calculates the radiance derivatives
C
IMPLICIT REAL*8 (V)
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
REAL KSUBL(0:ND2)
C
DIMENSION RADO(0:ND2),OPDT(0:ND2)
DIMENSION RPRIME(NDIM)
DIMENSION TRALYR(*)
DIMENSION PNLHDR(2),XKMOLC(2),ODACUM(2)
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /BUFPNL/ V1P,V2P,DVP,NLIMP
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KDUMY,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /EMDMSV/ BBSAV(2410),BBASAV(2410),XXSAV(2410)
C
EQUIVALENCE (PNLHDR(1),V1P)
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) THEN
WRITE(*,*) 'End of KFILE ',KFILE
GOTO 10
ENDIF
C
C Set number of words to be input as 2400 or the number of
C points between the absolute end point V2 and the panel
C starting point V1P, whichever is smaller. If you already
C have the point (from a previous panel), then skip over
C remaining panels.
C
XLIMP = (V2-V1P)/DVP + 1.00001
IF (XLIMP.GT.2400) THEN
NLIMP = 2400
ELSEIF (XLIMP.LT.1.) THEN
GOTO 10
ELSE
NLIMP = XLIMP
ENDIF
C
C Read in absorptance coefficients
C
CALL BUFIN (KFILE,KEOF,XKMOLC(1),NLIMP)
C
C Read in total optical depths if LAYER < NLAYER
C
IF (LAYER.LT.NLAYER) THEN
CALL BUFIN (KODTOT,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) THEN
WRITE(IPR,900) KODTOT,KFILE
STOP 'IN SUBROUTINE ADERIV: SEE OUTPUT FILE'
ENDIF
C
C Set number of words to be input as 2400 or the number of
C points between the absolute end point V2 and the panel
C starting point V1P, whichever is smaller. If you already
C have the point (from a previous panel), then skip over
C remaining panels.
C
XLIMP = (V2-V1P)/DVP + 1.00001
IF (XLIMP.GT.2400) THEN
NLIMP = 2400
ELSEIF (XLIMP.LT.1.) THEN
GOTO 10
ELSE
NLIMP = XLIMP
ENDIF
C
C Read in total optical depths
C
CALL BUFIN (KODTOT,KEOF,ODACUM(1),NLIMP)
ENDIF
C
C Calculate layer derivatives,
C
C dR/du = T{n-1}*[d{nonsource]/du + d{source}/du], where
C
C d{nonsource}/du = -Km*exp{-tau}*Rn,
C
C and
C
C d{source}/du = (Km/(1+a*tau)) *
C [exp{-tau} *
C (Bbar + B*a*tau + a*(B-Bbar)/(1+a*tau)) +
C a*(B-Bbar)/(1+a*tau)],
C
C such that
C Km is layer absorbtance coefficients
C tau is layer optical depths
C Rn is the radiance of the layer
C previously done
C
C Numerical substitutions:
C
C Y = 1/(1+a*tau)
C YZ = a*(B-Bbar)/(1+a*tau)
C AA = a
C XXSAV = a*tau
C BBSAV = Bbar
C BBASAV = B
C
C
C When calculating the derivative of the layer nearest the observer,
C omit the total accumulated transmittance, TRACCM
C
10 AA = 0.278
c write(*,*) 'ipathl = ',ipathl
IF (IPATHL.EQ.1.OR.IPATHL.EQ.3) THEN
c write(*,*) 'layer = ',layer,' nlayer = ',nlayer
IF (LAYER.NE.NLAYER) THEN
DO 20 I = 1, NLIM
vtest = v1po + dvpo*(i-1)
Y = 1./(1.+XXSAV(I))
YZ = Y*AA*(BBASAV(I)-BBSAV(I))
TRACCM = EXP(-OPDT(I))
SOURCE = Y*KSUBL(I)*
* (TRALYR(I)*(BBSAV(I)+BBASAV(I)*XXSAV(I)-YZ)+YZ)
SRCNON = KSUBL(I)*TRALYR(I)*RADO(I)
RPRIME(I) = TRACCM*(SOURCE - SRCNON)
c if ((i.le.10).and.(layer.eq.1)) then
c write(*,*) '#1 i, gnu = ',i,vtest
c write(*,*) ' y = ',y
c write(*,*) ' yz = ',yz
c write(*,*) ' traccm = ',traccm
c write(*,*) ' source = ',source
c write(*,*) ' ksubl = ',ksubl(i)
c write(*,*) ' tralyr = ',tralyr(i)
c write(*,*) ' xxsav = ',xxsav(i)
c write(*,*) ' srcnon = ',srcnon
c write(*,*) ' ksubl = ',ksubl(i)
c write(*,*) ' tralyr = ',tralyr(i)
c write(*,*) ' rado = ',rado(i)
c write(*,*) ' rprime = ',rprime(i)
c endif
20 CONTINUE
ELSE
DO 30 I = 1, NLIM
vtest = v1po + dvpo*(i-1)
Y = 1./(1.+XXSAV(I))
YZ = Y*AA*(BBASAV(I)-BBSAV(I))
SOURCE = Y*KSUBL(I)*
* (TRALYR(I)*(BBSAV(I)+BBASAV(I)*XXSAV(I)-YZ)+YZ)
SRCNON = KSUBL(I)*TRALYR(I)*RADO(I)
RPRIME(I) = SOURCE - SRCNON
c if ((i.le.10).and.(layer.eq.1)) then
c write(*,*) '#2 i, gnu = ',i,vtest
c write(*,*) ' y = ',y
c write(*,*) ' yz = ',yz
c write(*,*) ' source = ',source
c write(*,*) ' ksubl = ',ksubl(i)
c write(*,*) ' tralyr = ',tralyr(i)
c write(*,*) ' xxsav = ',xxsav(i)
c write(*,*) ' srcnon = ',srcnon
c write(*,*) ' ksubl = ',ksubl(i)
c write(*,*) ' tralyr = ',tralyr(i)
c write(*,*) ' rado = ',rado(i)
c write(*,*) ' rprime = ',rprime(i)
c endif
30 CONTINUE
ENDIF
c write(*,*) '**************************'
ELSEIF (IPATHL.EQ.2) THEN
STOP 'ADERIV NOT SET FOR IPATHL = 2'
c ELSEIF (IPATHL.EQ.3) THEN
c STOP 'ADERIV NOT SET FOR IPATHL = 3'
ELSEIF (IPATHL.EQ.-1) THEN
STOP 'ADERIV NOT SET FOR IPATHL = -1'
ENDIF
C
RETURN
C
900 FORMAT ('KODTOT, ',I2.2,', reached end prior to end of KFILE, ',
* I2.2)
C
END
C
C ---------------------------------------------------------------
C
SUBROUTINE TDERIV (KFILE,KODTOT,RPRIME,RADO,TRALYR,KSUBL,XKMOLC, 2,2
* OPDT,ODACUM,NLIM,NDIM,ND2,IPATHL,LAYER,NLAYER,
* v1po,dvpo)
C
C This subroutine combines the Planck function derivative
C (calculated in SUBROUTINE EMDT) and the layer transmittance
C and then calculates the radiance derivatives with respect to
C temperature
C
IMPLICIT REAL*8 (V)
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
REAL KSUBL(0:ND2)
C
DIMENSION RADO(0:ND2),OPDT(0:ND2)
DIMENSION RPRIME(NDIM)
DIMENSION TRALYR(*)
DIMENSION PNLHDR(2),XKMOLC(2),ODACUM(2)
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /BUFPNL/ V1P,V2P,DVP,NLIMP
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KDUMY,KPANEL,LINFIL,NFILE,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /EMDMSV/ BBSAV(2410),BBASAV(2410),XXSAV(2410)
COMMON /EMDTSV/ BBDSAV(2410)
C
EQUIVALENCE (PNLHDR(1),V1P)
C
c CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
c IF (KEOF.LE.0) THEN
c WRITE(*,*) 'End of KFILE ',KFILE
c GOTO 10
c ENDIF
C
C Set number of words to be input as 2400 or the number of
C points between the absolute end point V2 and the panel
C starting point V1P, whichever is smaller. If you already
C have the point (from a previous panel), then skip over
C remaining panels.
C
XLIMP = (V2-V1P)/DVP + 1.00001
IF (XLIMP.GT.2400) THEN
NLIMP = 2400
ELSEIF (XLIMP.LT.1.) THEN
GOTO 10
ELSE
NLIMP = XLIMP
ENDIF
C
C Read in absorptance coefficients
C
c CALL BUFIN (KFILE,KEOF,XKMOLC(1),NLIMP)
C
C Read in total optical depths if LAYER < NLAYER
C
IF (LAYER.LT.NLAYER) THEN
CALL BUFIN (KODTOT,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) THEN
WRITE(IPR,900) KODTOT,KFILE
STOP 'IN SUBROUTINE ADERIV: SEE OUTPUT FILE'
ENDIF
C
C Set number of words to be input as 2400 or the number of
C points between the absolute end point V2 and the panel
C starting point V1P, whichever is smaller. If you already
C have the point (from a previous panel), then skip over
C remaining panels.
C
XLIMP = (V2-V1P)/DVP + 1.00001
IF (XLIMP.GT.2400) THEN
NLIMP = 2400
ELSEIF (XLIMP.LT.1.) THEN
GOTO 10
ELSE
NLIMP = XLIMP
ENDIF
C
C Read in total optical depths
C
CALL BUFIN (KODTOT,KEOF,ODACUM(1),NLIMP)
ENDIF
C
C Calculate layer derivatives,
C
C dR/dt = Tt*[(1-T)*dB/dt], where
C
C dB/dt = (h*c*gnu)/(k*t*t) * 1/(1-exp{-beta*gnu})*B(gnu,t)
C
C
C When calculating the derivative of the layer nearest the observer,
C omit the total accumulated transmittance, TRACCM
C
10 AA = 0.278
IF (IPATHL.EQ.1.OR.IPATHL.EQ.3) THEN
IF (LAYER.NE.NLAYER) THEN
DO 20 I = 1, NLIM
vtest = v1po + dvpo*(i-1)
c Y = 1./(1.+XXSAV(I))
c YZ = Y*AA*(BBASAV(I)-BBSAV(I))
TRACCM = EXP(-OPDT(I))
SOURCE = BBDSAV(I)*(1.-TRALYR(I))
c SRCNON =
RPRIME(I) = TRACCM*SOURCE
20 CONTINUE
ELSE
DO 30 I = 1, NLIM
vtest = v1po + dvpo*(i-1)
c Y = 1./(1.+XXSAV(I))
c YZ = Y*AA*(BBASAV(I)-BBSAV(I))
SOURCE = BBDSAV(I)*(1.-TRALYR(I))
c SRCNON =
RPRIME(I) = SOURCE
30 CONTINUE
ENDIF
ELSEIF (IPATHL.EQ.2) THEN
STOP 'TDERIV NOT SET FOR IPATHL = 2'
c ELSEIF (IPATHL.EQ.3) THEN
c STOP 'TDERIV NOT SET FOR IPATHL = 3'
ELSEIF (IPATHL.EQ.-1) THEN
STOP 'TDERIV NOT SET FOR IPATHL = -1'
ENDIF
C
RETURN
C
900 FORMAT ('KODTOT, ',I2.2,', reached end prior to end of KFILE, ',
* I2.2)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE FLXIN (V1P,V2P,DVP,NLIM,KFILE,EM,TR,KEOF,NPANLS) 4,19
C
IMPLICIT REAL*8 (V)
C
C SUBROUTINE FLXIN INPUTS OPTICAL DEPTH VALUES FROM KFILE AND
C CALCULATES FLUX FOR THE LAYER. THIS VERSION WORKS FOR AEROSOLS.
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C M.J. IACONO
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KDUMY,KPANEL,LINFIL,NFILA,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /LAMCHN/ ONEPL,ONEMI,EXPMIN,ARGMIN
COMMON /BUFPNL/ V1PBF,V2PBF,DVPBF,NLIMBF
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
DIMENSION PNLHDR(2),EM(*),TR(*)
C
EQUIVALENCE (FSCDID(1),IHIRAC) , (FSCDID(2),ILBLF4)
EQUIVALENCE (PNLHDR(1),V1PBF)
EQUIVALENCE (FSCDID(4),IAERSL)
C
CALL BUFIN (KFILE,KEOF,PNLHDR(1),NPHDRF)
IF (KEOF.LE.0) RETURN
CALL BUFIN (KFILE,KEOF,TR(1),NLIMBF)
C
C TR CONTAINS THE OPTICAL DEPTHS AT THIS STAGE
C
IF (IHIRAC.EQ.4) STOP ' IHIRAC=4 FLXIN '
C
C EM CONTAINS THE OPTICAL DEPTH CORRECTIONS FOR NLTE AT THIS STAGE
C
IF (NPANLS.LT.1.AND.IAERSL.EQ.0) WRITE (IPR,900)
IF (NPANLS.LT.1.AND.IAERSL.NE.0) WRITE (IPR,905)
C
EXT = 0.
ADEL = 0.
RADFN0 = 0.
RDEL = 0.
BB = 0.
BBDEL = 0.
BBA = 0.
BBDLA = 0.
C
V1P = V1PBF
V2P = V2PBF
DVP = DVPBF
NLIM = NLIMBF
VI = V1P-DVP
VIDV = VI
VIBB = VI
VAER = VI
VDUM = VI
BBLAST = -1.
BBLXTA = -2.
RDLAST = -1.
BBDUM = -4.
RDDUM = -1.
NLIM1 = 0
NLIM2 = 0
C
AA = 0.278
C
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBDUM)
IF (IAERSL.NE.0) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDDUM)
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
IAFBB = 0
IF (VITST.LT.VAER.AND.VITST.LT.VIBB) IAFBB = 1
IF (VAER.LT.VITST.AND.VAER.LT.VIBB) IAFBB = 2
ELSE
IAFBB = -1
ENDIF
C
C - THIS SECTION TREATS THE CASE WHERE THE LAYER CONTRIBUTES
C TO THE RADIATIVE TRANSFER ONLY ONCE
C
C - WITH XKTA=0 THIS ALGORITHM REVERTS TO THE ORIGINAL
C
IF (XKTB.GT.0.) STOP ' XKTB GT 0. FLXIN '
C
10 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFBB.EQ.-1) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
ELSEIF (IAFBB.EQ.0) THEN
BB = BBFN(VI,DVP,V2P,XKT,VIDV,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.1) THEN
RADFN0 = RADFNI(VI,DVP,XKT,VIDV,RDEL,RDLAST)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSEIF (IAFBB.EQ.2) THEN
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VIBB = -VIDV
BB = BBFN(VI,DVP,V2P,XKT,VIBB,BBDEL,BBLAST)
IF (XKTA.GT.0.) THEN
VIBB = -VIDV
BBA = BBFN(VI,DVP,V2P,XKTA,VIBB,BBDLA,BBLXTA)
ELSE
BBA = BB
BBDLA = BBDEL
ENDIF
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 20 I = NLIM1, NLIM2
ODVI = SECNT*TR(I)+EXT*RADFN0
C
XX = AA*ODVI
C
TR(I) = EXP(-ODVI)
EM(I) = (1.-TR(I))*(BB+XX*BBA)/(1.+XX)
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
BB = BB+BBDEL
BBA = BBA+BBDLA
20 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 10
C
RETURN
C
900 FORMAT ('0EMISSION AND TRANSMISSION (MOLECULAR) ')
905 FORMAT ('0EMISSION AND TRANSMISSION (AEROSOLS EFFECTS INCLUDED)')
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE FLINIT (NPTS,MFILE,JPATHL,TBND) 2,14
C
IMPLICIT REAL*8 (V)
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C M.J. IACONO
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON NEWEM(2410),NEWTR(2410)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYER,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMIHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYDUM,YI1,YID(10),LSTWDF
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILA,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
CHARACTER*40 CEXT,CYID
C
REAL NEWEM,NEWTR
C
DIMENSION XFILHD(2),OPNLHD(2)
DIMENSION EMLAYR(2),TRLAYR(2)
C
EQUIVALENCE (XFILHD(1),XID(1)) , (OPNLHD(1),V1PO)
EQUIVALENCE (NEWEM(1),EMLAYR(1)) , (NEWTR(1),TRLAYR(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
C
C *******************************************************************
C *** THIS SUBROUTINE COMPUTES THE EMISSION FOR THE FIRST LAYER ***
C *******************************************************************
C
C TBND IS THE BOUNDARY BLACK BODY TEMPERATUE
C
C IPATHL = 1 IS FOR THE DOWNWELLING FLUX CASE
C IPATHL = 3 IS FOR THE UPWELLING FLUX CASE
C
CALL CPUTIM (TIME)
C
WRITE (IPR,900) TIME
NPANLS = 0
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
C FOR AEROSOL RUNS, MOVE EXTID INTO YID
C
IF (IAERSL.GT.0) THEN
WRITE (CEXT,'(10A4)') EXTID
WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CYID(19:40) = CEXT(19:40)
READ (CYID,'(5A8)') (YID(I),I=3,7)
ENDIF
C
C READ IN DIRECTION COSINE
C
READ (IRD,905) DIRCOS
SECNT = 1.0/DIRCOS
SECANT = SECNT
SECNT0 = SECNT
WRITE (IPR,910) DIRCOS
C
C IF BOUNDARY PROPERTIES ARE SUPPLIED, AND DOWNWARD LOOKING
C CASE; SET IPATHL TO REFLECTED ATMOSPHERE CASE
C
IF (IBPROP.EQ.1.AND.IPATHL.EQ.1) IPATHL = -1
IEMIT = 1
FACT = 1.
TIMEM = 0.0
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
DO 10 MOL = 1, NMOL
WK(MOL) = WK(MOL)*FACT
10 CONTINUE
WBROAD = WBROAD*FACT
LAYR1 = LAYER
WRITE (IPR,915) LAYR1,LAYER,KFILE,MFILE
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
XKT = TAVE/RADCN2
XKTBND = TBND/RADCN2
IF (IPATHL.EQ.-1) STOP ' IPATH=-1 '
IF (IPATHL.EQ.0) STOP ' IPATH=0 '
IF (IPATHL.EQ.1) THEN
XKTA = TZL/RADCN2
XKTB = 0.
ENDIF
IF (IPATHL.EQ.2) STOP ' IPATH=2 '
IF (IPATHL.EQ.3) THEN
XKTA = TZU/RADCN2
XKTB = 0.
ENDIF
WRITE (IPR,920) IPATHL,IANT
C
20 CONTINUE
C
CALL CPUTIM (TIMEM1)
CALL FLXIN (V1PO,V2PO,DVPO,NLIMO,KFILE,EMLAYR,TRLAYR,KEOF,NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 50
VI = V1PO-DVPO
VIDVBD = VI
VIDVEM = VI
VIDVRF = VI
BBLAST = -1.
EMLAST = -1.
IF ((IPATHL.EQ.3).AND.(TBND.GT.0.)) THEN
C
NLIM1 = 0
NLIM2 = 0
EMDUM = 0.
BBDUM = 0.
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMDUM)
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBDUM)
IEMBB = 0
IF (VIDVBD.GT.VIDVEM) IEMBB = 1
C
30 NLIM1 = NLIM2+1
C
VI = V1PO+FLOAT(NLIM1-1)*DVPO
IF (IEMBB.EQ.0) THEN
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDV,BBDEL,BBLAST)
VIDVEM = -VIDV
EMISIV = EMISFN(VI,DVPO,VIDVEM,EMDEL,EMLAST)
ELSE
EMISIV = EMISFN(VI,DVPO,VIDV,EMDEL,EMLAST)
VIDVBD = -VIDV
BB = BBFN(VI,DVPO,V2PO,XKTBND,VIDVBD,BBDEL,BBLAST)
ENDIF
C
IF (VIDV.GE.9.E+4) THEN
NLIM2 = NLIMO+1
ELSE
NLIM2 = (VIDV-V1PO)/DVPO+1.001
ENDIF
NLIM2 = MIN(NLIM2,NLIMO)
C
DO 40 J = NLIM1, NLIM2
NEWEM(J) = EMLAYR(J)+TRLAYR(J)*BB*EMISIV
C
C Increment interpolation values
C
BB = BB+BBDEL
EMISIV = EMISIV+EMDEL
40 CONTINUE
C
IF (NLIM2.LT.NLIMO) GO TO 30
C
ENDIF
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,NEWEM,NEWTR,MFILE,NPTS,NPANLS)
GO TO 20
50 CALL CPUTIM (TIME1)
TIME = TIME1-TIME
WRITE (IPR,925) TIME,TIMEM
C
RETURN
C
900 FORMAT (' TIME AT THE START OF --FLINIT-- ',F10.3)
905 FORMAT (F10.8)
910 FORMAT ('0 ***** DIR. COSINE ***** ',/,7X,F10.8)
915 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,
* '0 INPUT FILE =',I5,' OUTPUT FILE =',I5)
920 FORMAT ('0 IPATHL AND IANT',2I5)
925 FORMAT (' TIME REQUIRED FOR --FLINIT-- ',F10.3,' --FLXIN-- ',
* F10.3)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE FLUXUP (NPTS,LFILE,MFILE,JPATHL,TBND) 1,18
C
IMPLICIT REAL*8 (V)
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C M.J. IACONO
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON RADN(2410),TRAN(2410),RADO(0:5000)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYDUM,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /BNDPRP/ TMPBND,BNDEMI(3),BNDRFL(3),IBPROP
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /XPANEL/ V1P,V2P,DVP,NLIM,RMIN,RMAX,NPNLXP,NSHIFT,NPTSS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILA,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /XME/ TRAO(0:5000)
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
DIMENSION A1(10),A2(10),A3(10),A4(10)
DIMENSION RADLYR(2),TRALYR(2),RADOI(2),TRAOI(2)
DIMENSION WKSAV(35)
C
CHARACTER*40 CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO)
EQUIVALENCE (RADO(1),RADOI(1)) , (TRAO(1),TRAOI(1))
EQUIVALENCE (RADN(1),RADLYR(1)) , (TRAN(1),TRALYR(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
DATA NDIM / 2410 /,ND2 / 5000 /
C
C
C
C ************************************************************
C ****** THIS SUBROUTINE DOES LAYER MERGE FOR RADIANCE ******
C ************************************************************
C
C IPATHL = 3 IS FOR THE UPWELLING FLUX CASE
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
TIMEM = 0.0
TIMRD = 0.0
TIMOT = 0.0
C
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
SECNT = SECANT
LAY1SV = LAYR1
DVL = DV
PL = PAVE
TL = TAVE
WTOTL = 0.
C
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
C
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
C
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
TAVK = TAVE
DVK = DV
FACT = 1.
C
IF (DVL.EQ.DVK) THEN
ITYPE = 0
ELSEIF (DVL.GT.DVK) THEN
ITYPE = DVK/(DVL-DVK)+0.5
ELSE
C
C DVL.LT.DVK
C
ITYPE = -INT(DVL/(DVK-DVL)+0.5)
ENDIF
IF (ITYPE.LT.0) STOP ' FLUXUP; ITYPE LT 0 '
C
WTOTK = 0.
LAYR1 = LAY1SV
WRITE (IPR,905) LAYR1,LAYER,KFILE,LFILE,MFILE
IEMIT = 1
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+WK(MOL)*FACT
WK(MOL) = WK(MOL)*FACT+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+WBROAD*FACT
WBROAD = WBROAD*FACT+WN2SAV
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
SECANT = SECNT
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
XKT = TAVK/RADCN2
C
WRITE (IPR,910) IPATHL,IANT
C
IF (IPATHL.EQ.3) THEN
XKTA = TZU/RADCN2
XKTB = 0.
ELSE
STOP ' FLUXUP; IPATHL '
ENDIF
C
ATYPE = ITYPE
LL = ITYPE+1
AP = 1.0/(ATYPE+1.0)
C
C A1, A2, A3 AND A4 ARE THE CONSTANTS
C FOR THE LAGRANGE 4 POINT INTERPOLATION
C
DO 30 JPG = 1, ITYPE
APG = JPG
IPL = JPG+1
P = 1.0-(AP*APG)
A1(IPL) = -P*(P-1.0)*(P-2.0)/6.0
A2(IPL) = (P**2-1.0)*(P-2.0)*0.5
A3(IPL) = -P*(P+1.0)*(P-2.0)*0.5
A4(IPL) = P*(P**2-1.0)/6.0
30 CONTINUE
C
C *** BEGINNING OF LOOP THAT DOES MERGE ***
C
NPE = 0
RADO(0) = 0.0
TRAO(0) = 0.0
V1PO = 0.0
V2PO = 0.0
DVPO = 0.0
C
40 CONTINUE
C
CALL CPUTIM (TIMEM1)
CALL FLXIN (V1P,V2P,DVP,NLIM,KFILE,RADLYR,TRALYR,KEOF,NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 80
II = 1
C
IF (V2PO.LE.V2P+DVPO) THEN
50 CALL CPUTIM (TIMEM1)
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 60
CALL BUFIN (LFILE,LEOF,RADOI(NPE+1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAOI(NPE+1),NLIMO)
CALL CPUTIM (TIMEM2)
TIMRD = TIMRD+TIMEM2-TIMEM1
NPE = NLIMO+NPE
IF (V2PO.LE.V2P+DVPO) GO TO 50
ENDIF
C
C ZERO POINT OF FIRST PANEL
C
60 IF (RADO(0).EQ.0.0.AND.TRAO(0).EQ.0.0) THEN
RADO(0) = RADO(1)
TRAO(0) = TRAO(1)
ENDIF
C
C END POINT OF LAST PANEL
C
IF (V2PO+DVPO.GE.V2) THEN
RADO(NPE+1) = RADO(NPE)
RADO(NPE+2) = RADO(NPE)
TRAO(NPE+1) = TRAO(NPE)
TRAO(NPE+2) = TRAO(NPE)
ENDIF
C
NPL = 1
C
C NPL IS LOCATION OF FIRST ELEMENT ON ARRAYS RADO AND TRAO
C
CALL FLUXNN (RADN,TRAN,RADO,TRAO,NLIM,NDIM,ND2,V1P,DVP,IPATHL,
* A1,A2,A3,A4,LL,NPL)
C
CALL CPUTIM (TIMEM1)
C
IF (TBND.GT.0.) CALL EMBND (V1P,V2P,DVP,NLIM,RADN,TRAN,TBND)
C
CALL EMOUT (V1P,V2P,DVP,NLIM,RADN,TRAN,MFILE,NPTS,NPANLS)
CALL CPUTIM (TIMEM2)
TIMOT = TIMOT+TIMEM2-TIMEM1
C
C NPL IS NOW LOCATION OF FIRST ELEMENT TO BE USED FOR NEXT PASS
C
IPL = -1
DO 70 NL = NPL, NPE
IPL = IPL+1
RADO(IPL) = RADO(NL)
TRAO(IPL) = TRAO(NL)
70 CONTINUE
C
NPE = IPL
C
GO TO 40
80 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,915) TIME1,TIM,TIMEM,TIMRD,TIMOT
C
RETURN
C
900 FORMAT ('0 THE TIME AT THE START OF FLUXUP IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5)
910 FORMAT ('0 IPATHL AND IANT',2I5)
915 FORMAT ('0 THE TIME AT THE END OF FLUXUP IS ',F12.3,/,F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE - FLXIN - ',F12.3,
* ' - READ - ',F12.3,' - EMOUT - ',F12.3)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE FLUXNN (RADLYR,TRALYR,RADO,TRAO,NLIM,NDIM,ND2,V1P,DVP, 1
* IPATHL,A1,A2,A3,A4,LL,NPL)
C
IMPLICIT REAL*8 (V)
C
C THIS SUBROUTINE CALCULATES THE NEW RADIANCE AND TRANSMISSION
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM: R.D. WORSHAM
C S.A. CLOUGH
C M.J. IACONO
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
DIMENSION RADLYR(NDIM),TRALYR(NDIM),RADO(0:ND2),TRAO(0:ND2),
* A1(*),A2(*),A3(*),A4(*)
C
LLM1 = LL-1
LLM1 = MAX(LLM1,1)
C
C LOOPING OVER POINTS WITH SAME WEIGHTS
C
DO 30 NL = 1, LL
IPL = (NPL+NL-1)-LLM1
IF (NL.GT.1) IPL = IPL-1
C
IF (NL.EQ.1) THEN
C
C EXACT FREQUENCY - NO INTERPOLATION
C
DO 10 I = NL, NLIM, LL
IPL = IPL+LLM1
RADLYR(I) = TRALYR(I)*RADO(IPL)+RADLYR(I)
TRALYR(I) = TRALYR(I)*TRAO(IPL)
10 CONTINUE
C
C NOT EXACT FREQUENCY - INTERPOLATE RESULT
C
ELSE
C
A1N = A1(NL)
A2N = A2(NL)
A3N = A3(NL)
A4N = A4(NL)
C
DO 20 I = NL, NLIM, LL
IPL = IPL+LLM1
C
C INTERPOLATE THE OLD RADIANCE
C
RADLYR(I) = TRALYR(I)*(A1N*RADO(IPL-1)+A2N*RADO(IPL)+
* A3N*RADO(IPL+1)+A4N*RADO(IPL+2))+RADLYR(I)
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRALYR(I) = TRALYR(I)*(A1N*TRAO(IPL-1)+A2N*TRAO(IPL)+
* A3N*TRAO(IPL+1)+A4N*TRAO(IPL+2))
20 CONTINUE
C
C
ENDIF
C
30 CONTINUE
C
NPL = IPL
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE FLUXDN (NPTS,LFILE,MFILE,JPATHL,TBND) 1,30
C
IMPLICIT REAL*8 (V)
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
C LAST MODIFICATION: 14 AUGUST 1991
C
C IMPLEMENTATION: R.D. WORSHAM
C
C ALGORITHM REVISIONS: S.A. CLOUGH
C M.J. IACONO
C R.D. WORSHAM
C J.L. MONCET
C
C
C ATMOSPHERIC AND ENVIRONMENTAL RESEARCH INC.
C 840 MEMORIAL DRIVE, CAMBRIDGE, MA 02139
C
C----------------------------------------------------------------------
C
C WORK SUPPORTED BY: THE ARM PROGRAM
C OFFICE OF ENERGY RESEARCH
C DEPARTMENT OF ENERGY
C
C
C SOURCE OF ORIGINAL ROUTINE: AFGL LINE-BY-LINE MODEL
C
C FASCOD3
C
C
CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
C
COMMON RADN(2410),TRAN(2410),RADLYR(-1:4818)
COMMON /MANE/ P0,TEMP0,NLAYER,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYDUM,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /CONSTS/ PI,PLANCK,BOLTZ,CLIGHT,AVOG,RADCN1,RADCN2
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /EMHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYER ,YI1,YID(10),LSTWDF
COMMON /OPANL/ V1PO,V2PO,DVPO,NLIMO
COMMON /XPANEL/ V1P,V2P,DVP,NLIM,RMIN,RMAX,NPNLXP,NSHIFT,NPTSS
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILA,IAFIL,IEXFIL,
* NLTEFL,LNFIL4,LNGTH4
COMMON /XMI/ TRALYR(-1:4818)
COMMON /RMRG/ XKT,XKTA,XKTB,SECNT
C
DIMENSION XFILHD(2),PNLHDR(2),OPNLHD(2)
DIMENSION A1(0:100),A2(0:100),A3(0:100),A4(0:100)
DIMENSION RADO(2),TRAO(2)
DIMENSION WKSAV(35)
C
CHARACTER*40 CYID
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHDR(1),V1P),
* (OPNLHD(1),V1PO)
EQUIVALENCE (RADN(1),RADO(1)) , (TRAN(1),TRAO(1)),
* (FSCDID(4),IAERSL) , (FSCDID(5),IEMIT),
* (FSCDID(7),IPLOT) , (FSCDID(8),IPATHL),
* (FSCDID(16),LAYR1)
C
C ************************************************************
C ****** THIS SUBROUTINE DOES LAYER MERGE FOR EMISSION ******
C ****** USING FOUR POINT GENERAL INTERPOLATION ******
C ************************************************************
C
CALL CPUTIM (TIME)
WRITE (IPR,900) TIME
NPANLS = 0
TIMEM = 0.0
TIMRD = 0.0
TIMTB = 0.0
TIMOT = 0.0
C
CALL BUFIN (LFILE,LEOF,XFILHD(1),NFHDRF)
SECNT = SECANT
DVL = DV
LAY1SV = LAYR1
PL = PAVE
TL = TAVE
WTOTL = 0.
DO 10 MOL = 1, NMOL
WTOTL = WTOTL+WK(MOL)
WKSAV(MOL) = WK(MOL)
10 CONTINUE
WTOTL = WTOTL+WBROAD
WN2SAV = WBROAD
C
C FOR AEROSOL RUNS, MOVE YID (LFILE) INTO YID (MFILE)
C
IF (IAERSL.GT.0) WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CALL BUFIN (KFILE,KEOF,XFILHD(1),NFHDRF)
IF (IAERSL.GT.0) READ (CYID,'(5A8)') (YID(I),I=3,7)
IF (JPATHL.GE.1) IPATHL = JPATHL
PLAY = PAVE
TLAY = TAVE
C
IF (IPATHL.NE.1) STOP ' FLUXDN: IPATHL .NE. 1 '
C
XKT = TAVE/RADCN2
XKTA = TZL/RADCN2
XKTB = 0.
DVK = DV
LAYR1 = LAY1SV
FACT = 1.
IF (IPATHL.EQ.2.AND.IANT.EQ.0) FACT = 2.
ATYPE = 9.999E09
IF (DVK.EQ.DVL) ATYPE = 0.
IF (DVL.GT.DVK) ATYPE = DVK/(DVL-DVK)+0.5
IF (DVL.LT.DVK) ATYPE = -DVL/(DVK-DVL)-0.5
C
C IF (ATYPE .GT. 0) STOP ' FLUXDN; ATYPE GT 0 '
C
WTOTK = 0.
WRITE (IPR,905) LAYR1,LAYER,KFILE,LFILE,MFILE,ATYPE
IEMIT = 1
DO 20 MOL = 1, NMOL
WTOTK = WTOTK+WK(MOL)*FACT
WK(MOL) = WK(MOL)*FACT+WKSAV(MOL)
20 CONTINUE
WTOTK = WTOTK+WBROAD*FACT
WBROAD = WBROAD*FACT+WN2SAV
PAVE = (PL*WTOTL+PAVE*WTOTK)/(WTOTL+WTOTK)
TAVE = (TL*WTOTL+TAVE*WTOTK)/(WTOTL+WTOTK)
SECANT = SECNT
DV = DVL
C
C WK IS NOW THE ACCUMULATED SUM OF THE COLUMN DENSITIES
C
CALL BUFOUT (MFILE,XFILHD(1),NFHDRF)
DVXM = DV
C
IF (ATYPE.EQ.0.) THEN
C
C 1/1 RATIO ONLY
C
C
30 CONTINUE
CALL CPUTIM (TIMEM1)
CALL FLXIN (V1P,V2P,DVP,NLIM,KFILE,RADLYR(1),TRALYR(1),KEOF,
* NPANLS)
CALL CPUTIM (TIMEM2)
TIMEM = TIMEM+TIMEM2-TIMEM1
IF (KEOF.LE.0) GO TO 110
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
CALL BUFIN (LFILE,LEOF,RADO(1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAO(1),NLIMO)
CALL CPUTIM (TIMEM3)
TIMRD = TIMRD+TIMEM3-TIMEM2
DO 40 I = 1, NLIM
RADN(I) = RADO(I)*TRALYR(I)+RADLYR(I)
TRAN(I) = TRALYR(I)*TRAO(I)
40 CONTINUE
CALL CPUTIM (TIMEM1)
IF (TBND.GT.0.)
* CALL EMBND (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,TBND)
C
CALL CPUTIM (TIMEM2)
TIMTB = TIMTB+TIMEM2-TIMEM1
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,MFILE,NPTS,NPANLS)
CALL CPUTIM (TIMEM3)
TIMOT = TIMOT+TIMEM3-TIMEM2
GO TO 30
C
ENDIF
C
C ALL RATIOS EXCEPT 1/1
C
DO 50 JP = 0, 100
APG = JP
P = 0.01*APG
C
C THE FOLLOW ARE THE CONSTANTS FOR THE LAGRANGE 4 POINT
C INTERPOLATION
C
A1(JP) = -P*(P-1.0)*(P-2.0)/6.0
A2(JP) = (P**2-1.0)*(P-2.0)*0.5
A3(JP) = -P*(P+1.0)*(P-2.0)*0.5
A4(JP) = P*(P**2-1.0)/6.0
50 CONTINUE
C
C *** BEGINNING OF LOOP THAT DOES MERGE ***
C
NPE = 0
RADLYR(0) = 0.0
TRALYR(0) = 0.0
V1P = 0.0
V2P = 0.0
DVP = 0.0
V1PO = 0.0
V2PO = 0.0
DVPO = 0.0
KEOF = 1
C
60 CONTINUE
C
CALL CPUTIM (TIMEM1)
CALL BUFIN (LFILE,LEOF,OPNLHD(1),NPHDRF)
IF (LEOF.LE.0) GO TO 110
CALL BUFIN (LFILE,LEOF,RADO(1),NLIMO)
CALL BUFIN (LFILE,LEOF,TRAO(1),NLIMO)
CALL CPUTIM (TIMEM2)
TIMRD = TIMRD+TIMEM2-TIMEM1
II = 1
C
IF (V2P.LE.V2PO+DVP .AND. KEOF.GT.0) THEN
70 CALL CPUTIM (TIMEM2)
CALL FLXIN (V1P,V2P,DVP,NLIM,KFILE,RADLYR(NPE+1),
* TRALYR(NPE+1),KEOF,NPANLS)
CALL CPUTIM (TIMEM3)
TIMEM = TIMEM+TIMEM3-TIMEM2
IF (KEOF.LE.0) GO TO 80
V1P = V1P-FLOAT(NPE)*DVP
NPE = NLIM+NPE
IF (V2P.LE.V2PO+DVP) GO TO 70
ENDIF
C
C ZERO POINT OF FIRST PANEL
C
80 IF (RADLYR(0).EQ.0.0.AND.TRALYR(0).EQ.0.0) THEN
TRALYR(-1) = TRALYR(1)
RADLYR(-1) = RADLYR(1)
TRALYR(0) = TRALYR(1)
RADLYR(0) = RADLYR(1)
ENDIF
C
C END POINT OF LAST PANEL
C
IF (V2P+DVP.GE.V2) THEN
TRALYR(NPE+1) = TRALYR(NPE)
RADLYR(NPE+1) = RADLYR(NPE)
TRALYR(NPE+2) = TRALYR(NPE)
RADLYR(NPE+2) = RADLYR(NPE)
ENDIF
C
C NPL IS LOCATION OF FIRST ELEMENT ON ARRAYS RADO AND TRAO
C
NPL = 1
C
RATDV = DVL/DVK
C
C FJJ IS OFFSET BY 2. FOR ROUNDING PURPOSES
C
FJ1DIF = (V1PO-V1P)/DVP+1.+2.
C
C ***** BEGINNING OF LOOP THAT DOES MERGE *****
C
DO 90 II = 1, NLIMO
C
FJJ = FJ1DIF+RATDV*FLOAT(II-1)
JJ = IFIX(FJJ)-2
C
JP = (FJJ-FLOAT(JJ))*100.-199.5
C
C INTERPOLATE THE OLD TRANSMISSION
C
TRNLYR = A1(JP)*TRALYR(JJ-1)+A2(JP)*TRALYR(JJ)+
* A3(JP)*TRALYR(JJ+1)+A4(JP)*TRALYR(JJ+2)
C
C INTERPOLATE THE OLD EMISSION
C
RADN(II) = RADO(II)*TRNLYR+(A1(JP)*RADLYR(JJ-1)+
* A2(JP)*RADLYR(JJ)+A3(JP)*RADLYR(JJ+1)+
* A4(JP)*RADLYR(JJ+2))
C
TRAN(II) = TRNLYR*TRAO(II)
C
90 CONTINUE
C
NPL = JJ-1
C
CALL CPUTIM (TIMEM1)
IF (TBND.GT.0.) CALL EMBND (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,TBND)
C
CALL CPUTIM (TIMEM2)
TIMTB = TIMTB+TIMEM2-TIMEM1
CALL EMOUT (V1PO,V2PO,DVPO,NLIMO,RADN,TRAN,MFILE,NPTS,NPANLS)
CALL CPUTIM (TIMEM3)
TIMOT = TIMOT+TIMEM3-TIMEM2
C
C NPL IS NOW LOCATION OF FIRST ELEMENT TO BE USED FOR NEXT PASS
C
IPL = -2
DO 100 NL = NPL, NPE
IPL = IPL+1
TRALYR(IPL) = TRALYR(NL)
RADLYR(IPL) = RADLYR(NL)
100 CONTINUE
C
V1P = V1P+FLOAT(NPL+1)*DVP
NPE = IPL
C
GO TO 60
110 CONTINUE
C
CALL CPUTIM (TIME1)
TIM = TIME1-TIME
WRITE (IPR,910) TIME1,TIM,TIMEM,TIMRD,TIMTB,TIMOT
C
RETURN
C
900 FORMAT ('0 THE TIME AT THE START OF FLUXDN IS ',F12.3)
905 FORMAT ('0 INITIAL LAYER',I5,' FINAL LAYER',I5,/,'0 FILE ',I5,
* ' MERGED WITH FILE ',I5,' ONTO FILE',I5,' WITH XTYPE=',
* G15.5)
910 FORMAT ('0 THE TIME AT THE END OF FLUXDN IS ',F12.3/F12.3,
* ' SECS WERE REQUIRED FOR THIS MERGE - FLXIN - ',F12.3,
* ' - READ - ',F12.3,' - EMBND - ',F12.3,' - EMOUT - ',
* F12.3)
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE GETEXT (IEXFIL,LYRNOW,IEMITT) 3,6
C
IMPLICIT REAL*8 (V)
C
C THIS SUBROUTINE HAS BEEN MODIFIED TO ALSO READ IN THE ASYMMETRY
C FACTORS AND TO SEARCH FOR THE PROPER LAYER IF DESIRED.
C
C JAN 1986 (A.E.R. INC.)
C
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFIL,IEXFLD,
* NLTEFL,LNFIL4,LNGTH4
C
C ROUTINE TO BUFFER IN AEROSOL ABSORPTION AND SCATTERRING
C FROM TAPE 20 INTO COMMON ABSORB SCATTR, AND ASYMT
C
COMMON /PNLHDR/ V1P,V2P,DVP,NLIM,LDUM
COMMON /MANE/ P0,TEMP0,NLAYRS,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYRFX,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDA/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
COMMON /ABSORA/ V1ABS,V2ABS,DVABS,NPTABS,ABSRB(2025)
COMMON /SCATTA/ V1SC,V2SC,DVSC,NPTSC,SCTTR(2025)
COMMON /ASYMMA/ V1AS,V2AS,DVAS,NPTAS,ASYMT(2025)
DIMENSION APNLHD(2),AFILHD(2)
C
CHARACTER*40 CEXT
C
EQUIVALENCE (APNLHD(1),V1P) , (AFILHD(1),XID(1))
C
IF (IEMITT.EQ.0) THEN
CALL BUFIN (IEXFIL,IEOF,AFILHD(1),NFHDRF)
C
C MOVE YID INTO EXTID
C
WRITE (CEXT,'(5A8)') (YID(I),I=3,7)
READ (CEXT,'(10A4)') EXTID
ENDIF
C
IF (LYRNOW.EQ.1.AND.IEMITT.EQ.1) THEN
REWIND IEXFIL
CALL BUFIN (IEXFIL,IEOF,AFILHD(1),NFHDRF)
C
C MOVE YID INTO EXTID
C
WRITE (CEXT,'(5A8)') (YID(I),I=3,7)
READ (CEXT,'(10A4)') EXTID
ENDIF
C
LAYER = 0
C
10 DO 20 I = 1, 2025
ABSRB(I) = 0.
SCTTR(I) = 0.
ASYMT(I) = 0.
20 CONTINUE
C
CALL BUFIN (IEXFIL,IEOF,APNLHD(1),NPHDRF)
C
LAYER = LAYER+1
IF (IEOF.LE.0) STOP 'GETEXT; IEXFIL EMPTY'
C
CALL BUFIN (IEXFIL,IEOF,ABSRB(1),NLIM)
CALL BUFIN (IEXFIL,IEOF,SCTTR(1),NLIM)
CALL BUFIN (IEXFIL,IEOF,ASYMT(1),NLIM)
C
V1ABS = V1P
V1SC = V1P
V1AS = V1P
V2ABS = V2P
V2SC = V2P
V2AS = V2P
DVABS = DVP
DVSC = DVP
DVAS = DVP
NPTABS = NLIM
NPTSC = NLIM
NPTAS = NLIM
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE ADARSL (NNPTS,IEXFIL,MFILE,IAFIL,IEMIT) 1,12
C
IMPLICIT REAL*8 (V)
C
C ROUTINE TO ADD ABSORPTION AND SCATTERING TO THE TRANSMITTANCE
C VALUES AT EACH POINT. THE AEROSOL VALUES ARE STORED IN
C COMMON ABSORB AND COMMON SCATTR.
C
PARAMETER (MXFSC=200, MXLAY=MXFSC+3,MXZMD=3400,
* MXPDIM=MXLAY+MXZMD,IM2=MXPDIM-2,MXMOL=35,MXTRAC=22)
C
COMMON R1(2410)
COMMON /ABSPNL/ V1P,V2P,DVP,NLIM,NSHFT,NPTS
COMMON /MANE/ P0,TEMP0,NLAYRS,DVXM,H2OSLF,WTOT,ALBAR,ADBAR,AVBAR,
* AVFIX,LAYRFX,SECNT0,SAMPLE,DVSET,ALFAL0,AVMASS,
* DPTMIN,DPTFAC,ALTAV,AVTRAT,TDIFF1,TDIFF2,ALTD1,
* ALTD2,ANGLE,IANT,LTGNT,LH1,LH2,IPFLAG,PLAY,TLAY,
* EXTID(10)
COMMON /MSACCT/ IOD,IDIR,ITOP,ISURF,MSPTS,MSPANL(MXLAY),
* MSPNL1(MXLAY),
* MSLAY1,ISFILE,JSFILE,KSFILE,LSFILE,MSFILE,IEFILE,
* JEFILE,KEFILE
COMMON /IFIL/ IRD,IPR,IPU,NOPR,NFHDRF,NPHDRF,NFHDRL,NPHDRL,
* NLNGTH,KFILE,KPANEL,LINFIL,NFILE,IAFLD,IEXFLD,
* NLTEFL,LNFIL4,LNGTH4
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WBROAD,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
C
EQUIVALENCE (XFILHD(1),XID(1)) , (PNLHD(1),V1P)
EQUIVALENCE (FSCDID(4),IAERSL)
C
CHARACTER*40 CEXT,CYID
C
DIMENSION PNLHD(4),XFILHD(2)
C
XKT0 = 0.6951*296.
CALL GETEXT (IEXFIL,LAYRS,IEMIT)
CALL BUFIN (MFILE,IEOF,XFILHD(1),NFHDRF)
C
C FOR AEROSOL RUNS, MOVE EXTID INTO YID
C
WRITE (CEXT,'(10A4)') EXTID
WRITE (CYID,'(5A8)') (YID(I),I=3,7)
CYID(19:40) = CEXT(19:40)
READ (CYID,'(5A8)') (YID(I),I=3,7)
C
CALL BUFOUT (IAFIL,XFILHD(1),NFHDRF)
NPANLS = 0
IF (NOPR.EQ.0) WRITE (IPR,900) XID,(YID(N),N=1,2)
10 CALL BUFIN (MFILE,IEOF,PNLHD(1),NPHDRF)
IF (IEOF.LE.0) GO TO 40
CALL BUFIN (MFILE,IEOF,R1(1),NLIM)
IF (NPANLS.EQ.0) VIDV = V1P-DVP
VAER = VIDV
VITST = VIDV
NPTS = NNPTS
RDLAST = -1.
NLIM1 = 0
NLIM2 = 0
RDDUM = 0.
AF = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
RDF = RADFNI(VI,DVP,XKT0,VITST,RDEL,RDDUM)
IAFRD = 0
IF (VITST.GT.VAER) IAFRD = 1
C
20 NLIM1 = NLIM2+1
C
VI = V1P+FLOAT(NLIM1-1)*DVP
IF (IAFRD.EQ.0) THEN
RADFN0 = RADFNI(VI,DVP,XKT0,VIDV,RDEL,RDLAST)
VAER = -VIDV
EXT = AERF(VI,DVP,VAER,ADEL,TAUSCT,TDEL,GFACT,GDEL)
ELSE
EXT = AERF(VI,DVP,VIDV,ADEL,TAUSCT,TDEL,GFACT,GDEL)
VITST = -VIDV
RADFN0 = RADFNI(VI,DVP,XKT0,VITST,RDEL,RDLAST)
ENDIF
C
NLIM2 = (VIDV-V1P)/DVP+1.001
NLIM2 = MIN(NLIM2,NLIM)
C
DO 30 K = NLIM1, NLIM2
R1(K) = R1(K)+EXT*RADFN0
C
C Increment interpolation values
C
EXT = EXT+ADEL
RADFN0 = RADFN0+RDEL
30 CONTINUE
C
IF (NLIM2.LT.NLIM) GO TO 20
C
CALL ABSOUT (V1P,V2P,DVP,NLIM,1,IAFIL,NPTS,R1,NPANLS)
GO TO 10
40 CONTINUE
IAERSL = 2
RETURN
C
900 FORMAT ('1',5X,' AEROSOLS'//1X,10A8,2X,2(1X,A8,1X)//,5X,
* 'FILE 20 AEROSOL EXTINCTIONS ADDED TO FILE 12 SENT TO ',
* 'FILE 14')
C
END
C
C ----------------------------------------------------------------
C
FUNCTION AERF (VI,DVI,VINXT,ADEL,TAUSCT,TDEL,GFACT,GDEL) 46
C
IMPLICIT REAL*8 (V)
C
C THIS FUNCTION CORRELATES THE AEROSOL FREQ. WITH THE LBLRTM
C FREQ. AND ADDS THE ABSORPTION TO THE
C SCATTERING TO PRODUCE THE EXTINCTION
C
C THIS FUNCTION HAS BEEN MODIFIED TO RETURN THE SCATTERING
C SEPARATELY, AND TO ALSO RETURN THE ASYMMETRY FACTOR.
C
C JAN 1986 (A.E.R. INC.)
C
COMMON /ABSORA/ V1ABS,V2ABS,DVABS,NPTABS,ABSRB(2025)
COMMON /SCATTA/ V1SC,V2SC,DVSC,NPTSC,SCTTR(2025)
COMMON /ASYMMA/ V1AS,V2AS,DVAS,NPTAS,ASYMT(2025)
C
character*8 XID, HMOLID, YID
real*8 SECANT, XALTZ
C
COMMON /FILHDR/ XID(10),SECANT,PAVE,TAVE,HMOLID(60),XALTZ(4),
* WK(60),PZL,PZU,TZL,TZU,WN2 ,DV ,V1 ,V2 ,TBOUND,
* EMISIV,FSCDID(17),NMOL,LAYRS ,YI1,YID(10),LSTWDF
C
EQUIVALENCE (XFILHD(1),XID(1))
DIMENSION XFILHD(2)
C
VDIF = VI-V1ABS
IAER = VDIF/DVABS+1.00001
VAER = V1ABS+DVABS*FLOAT(IAER-1)
DERIVS = (SCTTR(IAER+1)-SCTTR(IAER))/DVABS
DERIVA = (ASYMT(IAER+1)-ASYMT(IAER))/DVABS
DERIV = DERIVS+(ABSRB(IAER+1)-ABSRB(IAER))/DVABS
C
C TAUSCT IS THE SCATTERING TERM
C
TAUSCT = SCTTR(IAER)+DERIVS*(VI-VAER)
C
C GFACT IS THE ASYMMETRY FACTOR
C
GFACT = ASYMT(IAER)+DERIVA*(VI-VAER)
AERF = SCTTR(IAER)+ABSRB(IAER)+DERIV*(VI-VAER)
C
C ADEL, TDEL & GDEL ARE THE CHANGE PER DVI
C
ADEL = DERIV*DVI
TDEL = DERIVS*DVI
GDEL = DERIVA*DVI
C
C Set next point for interpolation to be the point
C corresponding to the next element of the SCTTR,
C ASYMT, & ABSRB arrays
C
VINXT = VAER+DVABS
C
RETURN
C
END
C
C ----------------------------------------------------------------
C
SUBROUTINE LINTCO(V1,Z1,V2,Z2,VINT,ZINT,ZDEL) 2
C
C Linearly interpolates emission and reflection values which
C are directly read in from ASCII files
C
IMPLICIT REAL*8 (V)
C
C ZDEL is the slope of the line
C
ZDEL = (Z2-Z1)/(V2-V1)
C
C ZCEPT is the intercept for V = 0.0
C
ZCEPT = Z1 - ZDEL*V1
C
C Calculate ZINT value at VINT
C
ZINT = ZDEL*VINT + ZCEPT
C
RETURN
C
END
C ----------------------------------------------------------------