subroutine ptz(nlev,ntype,nmodel,pd,pd8) C C ** THIS PROGRAM CALCULATES TEMPERATURES ,H2O MIXING RATIOS ** C ** AND O3 MIXING RATIOS BY USING AN ANALYTICAL ** C ** FUNCTION WHICH APPROXIMATES ** C ** THE US STANDARD (1976). THIS IS ** C ** CALCULATED IN FUNCTION 'ANTEMP', WHICH IS CALLED BY THE ** C ** MAIN PROGRAM. THE FORM OF THE ANALYTICAL FUNCTION WAS ** C ** SUGGESTED TO ME IN 1971 BY RICHARD S. LINDZEN. ** C C*****THIS VERSION IS ONLY USABLE FOR 1976 US STD ATM AND OBTAINS C QUANTITIES FOR CO2 INTERPOLATION AND INSERTION INTO OPERA- C TIONAL RADIATION CODES C C definitions: c pd,pd8: pressures (mb) for data levels. pd is for the case c where p(sfc)=1013.25 mb; pd8 applies when p(sfc)=810.6 mb. c in either case, index 1 is at p=0 mb; index (nlev+2) is at the sfc. c prsint: same as pd, but with indices reversed,index 1 at c the surface, and index (nlev+1) at the top (nonzero) data level. c press: pressures used for quadratures (4 quad. pts.) c indices run as in prsint. c prout: same as pd. c tquad: temperature at quad. pts. dimension index increases c from the surface. c tmpint: temperature at quad. pts, saved over both height c and quadrature index. values come from tquad. c temp: layer-averaged temperature (bet. 2 values of c prsint or pd) obtained by simpson's rule. lowest index is c the surface layer. c tmpout: as in temp, but with indices reversed (index 1= c temp at top data level),others are layer averages c tmpmid: temp. at mid-point of layer (defined by adjacent c pd's. NOT averaged. index 1 is at top and is top data level. c tmpflx: temp. at layer boundaries, with index 1 at top c and given by top data level. c alt: height values (km) generated by ANTEMP. lowest c index = surface. c a: altitudes saves over both height and quaarature index. c wmxint: h2o mass mixing ratio in g/g for quad. pts. c indices as in tmpint. c wmix: layer-averaged h2o mixing ratios,with indices c as in temp. c wmxout: wmix with indices reversed, and index 1=top data c level. c omxint: o3 mass mixing ratio in g/g for quad. pts. c indices as in tmpint. c o3mix: layer-averaged o3 mixing ratios,with indices c as in temp. c omxout: omix with indices reversed, and index 1=top data c level. implicit real(selected_real_kind(15)) (a-h,o-z) CHARACTER*20 PROFIL dimension pd(200),pd8(200) dimension prsint(200),press(200),prout(200) dimension tquad(200),tmpint(200,4),temp(200) dimension tmpout(200),tmpmid(200),tmpflx(200) dimension alt(200),a(200,4) dimension wmxint(200,4),wmix(200),wmxout(200) dimension omxint(200,4),o3mix(200),omxout(200) dimension gtemp(200) C C DATA PROFIL/ 6 'US STANDARD 1976'/ DATA PSMAX/1013.250/ C nl=nlev nlp=nl+1 nlp2=nl+2 c DELZAP=0.5 R=8.31432 G0=9.80665 ZMASS=28.9644 AA=6356.766 ALT(1)=0.0 zero = 0.0 Tquad(1)=ANTEMP(6,zero) C*******DETERMINE THE PRESSURES (PRESS) PSTAR=PSMAX C IF (NTYPE.EQ.1) CALL SKYP(nlev,pd,pd8,gtemp) IF (NTYPE.EQ.2) CALL SKY80P(nlev,pd,pd8,gtemp) IF (NTYPE.EQ.0) CALL SIGP(nlev,ntype,nmodel,pd,pd8,gtemp) IF (NTYPE.EQ.3) CALL SIGP(nlev,ntype,nmodel,pd,pd8,gtemp) PD(NLP2)=PSTAR DO 40 N=1,NLP PRSINT(N)=PD(NLP2+1-N) 40 CONTINUE C *** CALCULATE TEMPS FOR SEVERAL PRESSURES TO DO QUADRATURE DO 504 NQ=1,4 DO 505 N=2,NLP 505 PRESS(N)=PRSINT(N)+0.25*(NQ-1)*(PRSINT(N-1)-PRSINT(N)) PRESS(1)=PRSINT(1) C********************* c print*, ' PRESS ', (press(n),n=1,nlp) DO 100 N=1,NLEV C C ** ESTABLISH COMPUTATATIONAL LEVELS BETWEEN USER LEVELS AT ** C ** INTERVALS OF APPROXIMATELY 'DELZAP' KM. ** C DLOGP=7.0*LOG(PRESS(N)/PRESS(N+1)) NINT=DLOGP/DELZAP NINT=NINT+1 ZNINT=NINT G=G0 DZ=R*DLOGP/(7.0*ZMASS*G*ZNINT) HT=ALT(N) C C ** CALCULATE HEIGHT AT NEXT USER LEVEL BY MEANS OF ** C ** RUNGE-KUTTA INTEGRATION. ** C DO 200 M=1,NINT RK1=ANTEMP(6,HT)*DZ RK2=ANTEMP(6,HT+0.5*RK1)*DZ RK3=ANTEMP(6,HT+0.5*RK2)*DZ RK4=ANTEMP(6,HT+RK3)*DZ HT=HT+0.16666667*(RK1+RK2+RK2+RK3+RK3+RK4) 200 CONTINUE ALT(N+1)=HT Tquad(N+1)=ANTEMP(6,HT) c print*, n, alt(n), alt(n+1), tquad(n+1) 100 CONTINUE c print*, ' TQUAD ', (tquad(n),n=1,nlp) DO 506 N=1,NLP TMPINT(N,NQ)=Tquad(N) A(N,NQ)=ALT(N) 506 CONTINUE DO 523 N=1,NLP CALL MIXRAT(6,PRESS(N),TMPINT(N,NQ),WMXINT(N,NQ)) WMXINT(N,NQ)=WMXINT(N,NQ)*1.0E-6 CALL OZONE(6,PRESS(N),TMPINT(N,NQ),OMXINT(N,NQ)) 523 CONTINUE 504 CONTINUE C ***APPLY SIMPSON'S RULE DO 507 N=2,NLP TEMP(N)=1./12.*(TMPINT(N-1,1)+TMPINT(N,1)+4.*TMPINT(N,2)+ 1 2.*TMPINT(N,3)+4.*TMPINT(N,4)) WMIX(N)=1./12.*(WMXINT(N-1,1)+WMXINT(N,1)+4.*WMXINT(N,2)+ 1 2.*WMXINT(N,3)+4.*WMXINT(N,4)) O3MIX(N)=1./12.*(OMXINT(N-1,1)+OMXINT(N,1)+4.*OMXINT(N,2)+ 1 2.*OMXINT(N,3)+4.*OMXINT(N,4)) 507 CONTINUE C***OUTPUT FOR LINE-BY-LINE CALCS TMPOUT(1)=TMPINT(NLP,1) WMXOUT(1)=WMXINT(NLP,1) OMXOUT(1)=OMXINT(NLP,1) TMPMID(1)=TMPINT(NLP,1) DO 520 I=2,NLP2 TMPOUT(I)=TEMP(NLP2+1-I) WMXOUT(I)=WMIX(NLP2+1-I) OMXOUT(I)=O3MIX(NLP2+1-I) TMPMID(I)=TMPINT(NLP2+1-I,3) 520 CONTINUE DO 521 I=1,NLP2 PROUT(I)=PD(I) 521 CONTINUE DO 5221 I=2,NLP2 TMPFLX(I)=TMPINT(NLP2+1-I,1) 5221 CONTINUE TMPFLX(1)=TMPFLX(2) C C ** CALCULATE WATER MIXING RATIO USING LUTHER PROGRAM C DO 300 N=1,NLP CALL MIXRAT(6,PRSINT(N),TMPINT(N,1),WMIX(N)) WMIX(N)=WMIX(N)*1.0E-06 300 CONTINUE C C CALCULATE OZONE MIXING RATIO USING SCHWARZKOPF PROGRAM C DO 400 N=1,NLP CALL OZONE(6,PRSINT(N),TMPINT(N,1),O3MIX(N)) 400 CONTINUE C c WRITE (10,101) PROFIL 101 FORMAT (1X,A20) c WRITE (10,201) 201 FORMAT(5X,' HEIGHT TEMPERATURE PRESSURE R(H2O) 1 R(O3)') c WRITE (10,202) A(1,1),TMPINT(1,1),PRSINT(1),WMIX(1),O3MIX(1) DO 210 N=2,NLP c WRITE (10,203) TMPOUT(NLP2+1-N),WMXOUT(NLP2+1-N),OMXOUT(NLP2+1-N) c WRITE (10,202) A(N,1),TMPINT(N,1),PRSINT(N),WMIX(N),O3MIX(N) 210 CONTINUE 203 FORMAT (1X,14X,F14.6,14X,2E14.6) 202 FORMAT(1X,2F14.6,E14.6,2E14.6) C DETERMINE WHAT FORMATS TO USE*** NREP=NLP/5 NREM=NLP-5*NREP IF (NREM.EQ.0) THEN NREP=NREP-1 NREM=5 ENDIF NREP2=NL/5 NREM2=NL-5*NREP2 IF (NREM2.EQ.0) THEN NREP2=NREP2-1 NREM2=5 ENDIF NREP3=NLP/4 NREM3=NLP-4*NREP3 IF (NREM3.EQ.0) THEN NREP3=NREP3-1 NREM3=4 ENDIF C***OUTPUT TEMPERATURES DO 800 IOUT=1,2 IF (IOUT.EQ.1) THEN WRITE (16,701) c WRITE (10,701) ELSE WRITE (16,702) c WRITE (10,702) ENDIF 701 FORMAT (6X,'DATA DTEMP /') 702 FORMAT (6X,'DATA STEMP /') NF=0 IF (NREP.NE.0) THEN DO 801 NR=1,NREP NS=5*(NR-1)+1 NF=NS+4 WRITE (16,656) (TMPINT(NLP2-N,1),N=NS,NF) if (iout.eq.2) WRITE (10,1656) (TMPINT(NLP2-N,1),N=NS,NF) 801 CONTINUE ENDIF NF=NF+1 IF (NREM.EQ.1) THEN WRITE (16,616) (TMPINT(NLP2-N,1),N=NF,NLP) if (iout.eq.2) WRITE (10,1616) (TMPINT(NLP2-N,1),N=NF,NLP) ENDIF IF (NREM.EQ.2) THEN WRITE (16,626) (TMPINT(NLP2-N,1),N=NF,NLP) if (iout.eq.2) WRITE (10,1626) (TMPINT(NLP2-N,1),N=NF,NLP) ENDIF IF (NREM.EQ.3) THEN WRITE (16,636) (TMPINT(NLP2-N,1),N=NF,NLP) if (iout.eq.2) WRITE (10,1636) (TMPINT(NLP2-N,1),N=NF,NLP) ENDIF IF (NREM.EQ.4) THEN WRITE (16,646) (TMPINT(NLP2-N,1),N=NF,NLP) if (iout.eq.2) WRITE (10,1646) (TMPINT(NLP2-N,1),N=NF,NLP) ENDIF IF (NREM.EQ.5) THEN WRITE (16,606) (TMPINT(NLP2-N,1),N=NF,NLP) if (iout.eq.2) WRITE (10,1606) (TMPINT(NLP2-N,1),N=NF,NLP) ENDIF 616 FORMAT (5X,'*',1X,F12.6,'/') 626 FORMAT (5X,'*',1X,F12.6,',',F12.6,'/') 636 FORMAT (5X,'*',1X,F12.6,',',F12.6,',',F12.6,'/') 646 FORMAT (5X,'*',1X,F12.6,',',F12.6,',',F12.6,',', 1 F12.6,'/') 656 FORMAT (5X,'*',1X,F12.6,',',F12.6,',',F12.6,',', 1 F12.6,',',F12.6,',') 606 FORMAT (5X,'*',1X,F12.6,',',F12.6,',',F12.6,',', 1 F12.6,',',F12.6,'/') 800 CONTINUE 1616 FORMAT (F13.6) 1626 FORMAT (2F13.6) 1636 FORMAT (3F13.6) 1646 FORMAT (4F13.6) 1656 FORMAT (5F13.6) 1606 FORMAT (5F13.6) C***OUTPUT GTEMP WRITE (16,706) c WRITE (10,706) 706 FORMAT (6X,'DATA GTEMP /') NF=0 IF (NREP3.NE.0) THEN DO 805 NR=1,NREP3 NS=4*(NR-1)+1 NF=NS+3 WRITE (16,648) (GTEMP(N),N=NS,NF) WRITE (10,1648) (GTEMP(N),N=NS,NF) 805 CONTINUE ENDIF NF=NF+1 IF (NREM3.EQ.1) THEN WRITE (16,618) (GTEMP(N),N=NF,NLP) WRITE (10,1618) (GTEMP(N),N=NF,NLP) ENDIF IF (NREM3.EQ.2) THEN WRITE (16,628) (GTEMP(N),N=NF,NLP) WRITE (10,1628) (GTEMP(N),N=NF,NLP) ENDIF IF (NREM3.EQ.3) THEN WRITE (16,638) (GTEMP(N),N=NF,NLP) WRITE (10,1638) (GTEMP(N),N=NF,NLP) ENDIF IF (NREM3.EQ.4) THEN WRITE (16,608) (GTEMP(N),N=NF,NLP) WRITE (10,1608) (GTEMP(N),N=NF,NLP) ENDIF C***OUTPUT WMIX IN GM/GM NF=0 WRITE (16,703) 703 FORMAT (6X,'DATA RR /') IF (NREP2.NE.0) THEN DO 802 NR=1,NREP2 NS=5*(NR-1)+1 NF=NS+4 WRITE (16,657) (WMIX(NLP2-N),N=NS,NF) 802 CONTINUE ENDIF NF=NF+1 IF (NREM2.EQ.1) THEN WRITE (16,617) (WMIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.2) THEN WRITE (16,627) (WMIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.3) THEN WRITE (16,637) (WMIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.4) THEN WRITE (16,647) (WMIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.5) THEN WRITE (16,607) (WMIX(NLP2-N),N=NF,NL) ENDIF 617 FORMAT (5X,'*',1X,E12.6,'/') 627 FORMAT (5X,'*',1X,E12.6,',',E12.6,'/') 637 FORMAT (5X,'*',1X,E12.6,',',E12.6,',',E12.6,'/') 647 FORMAT (5X,'*',1X,E12.6,',',E12.6,',',E12.6,',', 1 E12.6,'/') 657 FORMAT (5X,'*',1X,E12.6,',',E12.6,',',E12.6,',', 1 E12.6,',',E12.6,',') 607 FORMAT (5X,'*',1X,E12.6,',',E12.6,',',E12.6,',', 1 E12.6,',',E12.6,'/') C***OUTPUT O3MIX IN GM/GM NF=0 WRITE (16,704) 704 FORMAT (6X,'DATA QQO3 /') IF (NREP2.NE.0) THEN DO 803 NR=1,NREP2 NS=5*(NR-1)+1 NF=NS+4 WRITE (16,657) (O3MIX(NLP2-N),N=NS,NF) 803 CONTINUE ENDIF NF=NF+1 IF (NREM2.EQ.1) THEN WRITE (16,617) (O3MIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.2) THEN WRITE (16,627) (O3MIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.3) THEN WRITE (16,637) (O3MIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.4) THEN WRITE (16,647) (O3MIX(NLP2-N),N=NF,NL) ENDIF IF (NREM2.EQ.5) THEN WRITE (16,607) (O3MIX(NLP2-N),N=NF,NL) ENDIF C***OUTPUT PROUT IN CGS DO 629 N=2,NLP2 PROUT(N)=PROUT(N)*1.0E3 629 CONTINUE WRITE (16,705) 705 FORMAT (6X,'DATA PPRESS /') NF=0 IF (NREP3.NE.0) THEN DO 804 NR=1,NREP3 NS=4*(NR-1)+2 NF=NS+3 WRITE (16,648) (PROUT(N),N=NS,NF) 804 CONTINUE ENDIF NF=NF+1 IF (NREM3.EQ.1) THEN WRITE (16,618) (PROUT(N),N=NF,NLP2) ENDIF IF (NREM3.EQ.2) THEN WRITE (16,628) (PROUT(N),N=NF,NLP2) ENDIF IF (NREM3.EQ.3) THEN WRITE (16,638) (PROUT(N),N=NF,NLP2) ENDIF IF (NREM3.EQ.4) THEN WRITE (16,608) (PROUT(N),N=NF,NLP2) ENDIF 618 FORMAT (5X,'*',1X,E15.9,'/') 628 FORMAT (5X,'*',1X,E15.9,',',E15.9,'/') 638 FORMAT (5X,'*',1X,E15.9,',',E15.9,',',E15.9,'/') 648 FORMAT (5X,'*',1X,E15.9,',',E15.9,',',E15.9,',',E15.9,',') 608 FORMAT (5X,'*',1X,E15.9,',',E15.9,',',E15.9,',',E15.9,'/') 1618 FORMAT (E15.9) 1628 FORMAT (2E16.9) 1638 FORMAT (3E16.9) 1648 FORMAT (4E16.9) 1608 FORMAT (4E16.9) return END