/[MITgcm]/MITgcm_contrib/jscott/igsm/src/radia.F
ViewVC logotype

Annotation of /MITgcm_contrib/jscott/igsm/src/radia.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.2 - (hide annotations) (download)
Tue Aug 22 20:25:52 2006 UTC (18 years, 11 months ago) by jscott
Branch: MAIN
Changes since 1.1: +1 -1 lines
changed AGRID.COM -> AGRID.h

1 jscott 1.1 #include "ctrparam.h"
2     SUBROUTINE RADIA 5001.
3     C**** 5002.
4     C**** THIS SUBROUTINES ADDS THE RADIATION HEATING TO THE TEMPERATURES 5003.
5     C**** 5004.
6    
7     #include "BD2G04.COM"
8     #include "chem_para"
9     #include "chem_com"
10    
11     COMMON U,V,T,P,Q 5006.
12     COMMON/ADDALB/BVSURFA,XVSURFA,BNSURFA,XNSURFA
13     COMMON/WORK1/CONV(IM0,JM0,LM0),PK(IM0,JM0,LM0),PREC(IM0,JM0),
14     & TPREC(IM0,JM0), 5007.
15     * COSZ1(IO0,JM0),COSZ2(IO0,JM0),COSZA(IO0,JM0), 5008.
16     * TRINCG(IO0,JM0),BTMPW(IO0,JM0),SNFS(IO0,JM0,4),TNFS(IO0,JM0,4), 5009.
17     * TRHRS(IO0,JM0,3),SRHRS(IO0,JM0,3),ALB(IO0,JM0,9) 5010.
18     COMMON/WORK2/CLDSS(IM0,JM0,LM0),CLDMC(IM0,JM0,LM0), 5011.
19     * TOTCLD(36) 5012.
20     DIMENSION TRNFP0(JM0),TRNFP1(JM0),ALBJ(JM0,9)
21     real ODATA2(JM0,2),GDATA2(JM0,14),BDATA2(JM0,2),FDATA2(JM0,2),
22     * RQT2(JM0,3)
23     common/SURRAD/TRSURF(JM0,4),SRSURF(JM0,4)
24     common/FORAERSOL/FORSULF,FORBC
25     logical FORSULF,FORBC
26     C COMMON/WORK4/ IS BEING USED BY THE RADIATION ROUTINES 5013.
27     C 5014.
28     C RADCOM: CONTROL/INPUT PARAMETERS 5015.
29     C 5016.
30     COMMON/RADCOM/VADATA(11,4,3),DGLAT(46),DGLON(72),TMINSR,FULGAS(18)5017.
31     A ,FRACSL,RATQSL,FOGTSL,PTLISO,TLGRAD,TKCICE,FGOLDU(18)5018.
32     B ,FLONO3,FRAYLE,FCLDTR,FCLDSR,FALGAE,FMARCL,FEMTRA(6) 5019.
33     C ,WETTRA,WETSRA,DMOICE,DMLICE,LICETK,NTRCE,FZASRA(6) 5020.
34     D ,ID5(5),ITR(4),IMG(2),ILG(2),LAPGAS,KWVCON,NORMS0,NV 5021.
35     E ,KEEPRH,KEEPAL,ISOSCT,IHGSCT,KFRACC,KGASSR,KAERSR 5022.
36     F ,MARCLD,LAYTOP,LMR,LMRP,JMLAT,IMLON,KFORCE,LASTVC 5023.
37     C 5024.
38     C BASIC RADCOM INPUT DATA 5025.
39     C 5026.
40     G ,PLE(40),HLB(40),TLB(40),TLT(40),TL(40),U0GAS(40,9) 5027.
41     H ,ULGAS(40,9),TRACER(40,4),RTAU(40),QL(40),RHL(40) 5028.
42     I ,POCEAN,PEARTH,POICE,PLICE,AGESN,SNOWE,SNOWOI,SNOWLI 5029.
43     J ,TGO,TGE,TGOI,TGLI,TS,WS,WEARTH,ZOICE,FSPARE(200) 5030.
44     K ,S0,COSZ,PVT(11),BXA(153),SRBXAL(15,2),FRC(5),LUXGAS 5031.
45     L ,JYEARR,JDAYR,JLAT,ILON,MEANAL,KALVIS,ISPARE(25),SGPS5032.
46     C 5033.
47     C BASIC RADCOM OUTPUT DATA 5034.
48     C 5035.
49     M ,TRDFLB(40),TRUFLB(40),TRNFLB(40),TRFCRL(40),TRSLCR 5036.
50     N ,SRDFLB(40),SRUFLB(40),SRNFLB(40),SRFHRL(40),SRSLHR 5037.
51     O ,SRIVIS,SROVIS,PLAVIS,SRINIR,SRONIR,PLANIR,SRXATM(4) 5038.
52     P ,SRDVIS,SRUVIS,ALBVIS,SRDNIR,SRUNIR,ALBNIR,FSRNFG(4) 5039.
53     Q ,SRTVIS,SRRVIS,SRAVIS,SRTNIR,SRRNIR,SRANIR,FTRUFG(4) 5040.
54     R ,TRDFGW,TRUFGW,TRUFTW,BTEMPW,TRDFSL,TRUFSL,DTRUFG(4) 5041.
55     S ,TRSLTS,TRSLTG,TRSLWV,TRSLBS,TTRUFG,LBOTCL,LTOPCL 5042.
56     DIMENSION COE(39) 5043.
57     LOGICAL POLE,DC25,HPRNT,WRCLD,CLDFEED 5044.
58     #if ( defined OCEAN_3D )
59 jscott 1.2 #include "AGRID.h"
60 jscott 1.1 #endif
61    
62     #if ( defined CLM )
63     #include "CLM.COM"
64     #endif
65     c
66    
67     common /ATCO2/atm_co2(jm0)
68     common/conprn/HPRNT
69     common/COMCLD/READGHG,PCLOUD,WRCLD,NWRCLD,NWRCL,INYEAR,JNDAY
70     &,CFAEROSOL,ALFA,CFBC,cfvolaer
71     common/ghstfor/GHSFALB,GHSF,ALBCF,FVOLADD,STRARFOR,S0FOR,CO2FOR,
72     & CO2F,ghostfv(LM0+1),ghostf(LM0+1,JM0),ghflux(LM0+1,JM0)
73     logical GHSFALB,GHSF,GSOEQ,STRARFOR,S0FOR,CO2FOR
74     common/eqgso/GSOEQ
75     common/cldfdb/coefcl(3),CLDFEED
76     common/aexpc/AEXP,ISTRT1
77     common/ SNOWALB/FRSNALB,III
78     common/ S0XR/S0RATE,CFS0X
79     common/ BACKGRGHG/GHGBGR(5)
80     COMMON/CO2TRND/ALFFOR,CO2TR,YEARGT,CO2IN,INYRAD
81     dimension CLDSSF(JM0,LM0),CLDMCF(JM0,LM0)
82     &,BSO4LAND(JM0),BSO4OCEAN(JM0),BSO4TOTAL(JM0)
83     &,bcodmn(JM0,LM0,12)
84     dimension DSWSRF(jm0),DLWSRF(jm0),DSWVIS(jm0),DSWNIR(jm0)
85     integer PCLOUD
86     common/TSUR/TSURFC(JM0,0:13),TSURFT(JM0),TSURFD(JM0),DTSURF(JM0)
87     *,cfcld(JM0,3)
88     CHARACTER*4 JMNTHF,JMLAST
89     DATA JMLAST /'LAST'/
90     DATA TF/273.16/,TCIR/258.16/,STBO/.567257E-7/,IFIRST/1/,JDLAST/-9/5045.
91     DATA IRFIRST /1/
92     C **** CLEAR SKY
93     dimension SRHRCL(JM0),TRHRCL(JM0),ALBCL(JM0),SNP1CL(JM0),
94     *SNP0CL(JM0),TRINCL(JM0),TRP0CL(JM0),TRP1CL(JM0)
95     dimension RTAU0(LM0),TOTCLD0(LM0)
96     common/clrsk/CLEAR(JM0),NCLR(JM0),AJCLR(JM0,12),BJCLR(JM0,12),
97     * CJCLR(JM0,12)
98     integer CLEAR
99     C AJCLR
100     C 1 SW INC AT P0 RD (AJ(1))
101     C 2 SW ABS BELOW P0 RD (AJ(2))
102     C 3 SW ABS BELOW P1 RD (AJ(3))
103     C 4 SW ABS AT Z0 RD (AJ(6))
104     C 5 SW INC AT Z0 RD (AJ(5))
105     C 6 LW INC AT Z0 RD (AJ(67))
106     C 7 NET LW AT Z0 SF (AJ(9))
107     C 8 NET LW AT P0 RD (AJ(7))
108     C 9 NET LW AT P1 RD (AJ(8))
109     C 10 NET RAD AT P0 DG (AJ(10))
110     C 11 NET RAD AT P1 DG (AJ(11))
111     C 12 NET RAD AT Z0 DG (AJ(12))
112     C **** CLEAR SKY
113     C**** 5046.
114     C**** FDATA 2 LAND COVERAGE (1) 5047.
115     C**** 3 RATIO OF LAND ICE COVERAGE TO LAND COVERAGE (1) 5048.
116     C**** 5049.
117     C**** ODATA 1 OCEAN TEMPERATURE (C) 5050.
118     C**** 2 RATIO OF OCEAN ICE COVERAGE TO WATER COVERAGE (1) 5051.
119     C**** 5052.
120     C**** GDATA 1 OCEAN ICE SNOW AMOUNT (KG/M**2) 5053.
121     C**** 2 EARTH SNOW AMOUNT (KG/M**2) 5054.
122     C**** 3 OCEAN ICE TEMPERATURE OF FIRST LAYER (C) 5055.
123     C**** 4 EARTH TEMPERATURE OF FIRST LAYER (C) 5056.
124     C**** 5 EARTH WATER OF FIRST LAYER (KG/M**2) 5057.
125     C**** 6 EARTH ICE OF FIRST LAYER (KG/M**2) 5058.
126     C**** 11 AGE OF SNOW (DAYS) 5059.
127     C**** 12 LAND ICE SNOW AMOUNT (KG/M**2) 5060.
128     C**** 13 LAND ICE TEMPERATURE OF FIRST LAYER (C) 5061.
129     C**** 5062.
130     C**** BLDATA 1 COMPOSITE SURFACE WIND MAGNITUDE (M/S) 5063.
131     C**** 2 COMPOSITE SURFACE AIR TEMPERATURE (K) 5064.
132     C**** 5 FREE 5065.
133     C**** 5066.
134     C**** VDATA 1-8 EARTH RATIOS FOR THE 8 VEGETATION TYPES (1) 5067.
135     C**** 9 WATER FIELD CAPACITY OF FIRST LAYER (KG/M**2) 5068.
136     C**** 5069.
137     dimension GHGBGR0(5)
138    
139     IF(MODRD.EQ.0) IDACC(2)=IDACC(2)+1 5070.
140     IF (IFIRST.NE.1) GO TO 50 5071.
141     if(GHSF)then
142     do j=1,jm
143     ghflux(1,j)=ghostf(1,j)
144     do l=2,lm0+1
145     ghflux(l,j)= ghflux(l-1,j)+ghostf(l,j)
146     enddo
147     enddo
148     do l=1,lm0+1
149     print *,l,ghflux(l,1)
150     enddo
151     endif
152     JDAYR=JNDAY
153     JYEARR=INYEAR
154     c GHGBGR(1)=GHGBGR0(1)
155     c GHGBGR(2)=GHGBGR0(2)
156     c GHGBGR(3)=GHGBGR0(3)
157     c GHGBGR(4)=GHGBGR0(4)*1000.
158     c GHGBGR(5)=GHGBGR0(5)*1000.
159     c print *,'Background GHGs'
160     print '(5E12.4)',GHGBGR
161     nreadcld=0
162     nrbyyr=24*365/5
163     nrcldmax=20*nrbyyr
164     c print *,' CLOUDS for ',nrcldmax/nrbyyr,' years'
165     print *,'RADIA'
166     KTREND=-CO2
167     JDAY00=-1
168     if(CLDFEED)then
169     c print *,' radia.f coefcl=',coefcl
170     print *,' for low and middle clouds',coefcl(1)
171     print *,' for top clouds',coefcl(2)
172     print *,' for MC clouds',coefcl(3)
173     endif
174     print *,' READGHG=',READGHG
175     print *,' CFAEROSOL=',CFAEROSOL
176     #if ( defined PREDICTED_BC)
177     print *,'With black carbon forcing'
178     print *,' CFBC=',CFBC
179     nrdbc=0
180     read(769),bcodmn
181     ! do nm=1,12
182     ! do l=1,lm
183     ! do j=1,jm
184     ! bcodmn(j,l,nm)=bcod(1,j,l)
185     ! enddo
186     ! enddo
187     ! enddo
188     #endif
189     do j=1,jm
190     BSO4TOTAL(j)=1.0
191     BSO4LAND(j)=1.0
192     BSO4OCEAN(j)=1.0
193     enddo
194     print *,' separate caclulations for land and ocean'
195     DC25=.TRUE.
196     ! DC25=.FALSE.
197     if(DC25)then
198     print *,' with DC'
199     else
200     print *,' without DC'
201     print *,' subroutine COSZR'
202     end if
203     if(abs(PCLOUD-3.).gt.1.5.and..NOT.WRCLD)IFIRST=0 5072.
204     LMP1=LM+1 5072.1
205     DTCNDS=NCNDS*DT 5073.
206     C**** SET THE CONTROL PARAMETERS FOR THE RADIATION 5074.
207     JMLAT=JM 5074.1
208     if(JM.ne.24) then
209     DO J=1,JMLAT
210     DGLAT(J)=acos(COSP(J))*360./TWOPI
211     if(J.le.JMLAT/2)DGLAT(J)=-DGLAT(J)
212     END DO
213     endif
214     IMLON=IO 5074.2
215     LMR=LM+3 5075.
216     COEX=.01*GRAV*KAPA/RGAS 5076.
217     PSFMPT=PSF-PTOP 5077.
218     DO 30 L=1,LM 5078.
219     COE(L)=DTCNDS*COEX/DSIG(L) 5079.
220     30 PLE(L)=SIGE(L)*(PSF-PTOP)+PTOP 5080.
221     PLE(LMP1)=PTOP 5081.
222     PLE(LM+2)=.5*PTOP 5082.
223     PLE(LMR)=.2*PTOP 5083.
224     PLE(LMR+1)=1.E-5 5084.
225     DO 40 LR=LMP1,LMR 5085.
226     COE(LR)=DT*NRAD*COEX/(PLE(LR)-PLE(LR+1)) 5086.
227     QL(LR)=.3E-5 5087.
228     40 RTAU(LR)=0. 5088.
229     DPMICE=10. 5089.
230     C S0X=1. 5089.1
231     KSTREND=S0RATE*100.
232     print *,'S0RATE=',S0RATE,' KSTREND=',KSTREND
233     print *,'From first call of radia CFS0X=',CFS0X
234     #if ( defined FIXED_FOR )
235     print *,'Run with fixed forcing'
236     #endif
237     C
238     if(GSOEQ)then
239     if(KSTREND.ne.0)then
240     print *,'Wrong setting of GSOEQ and KSTREND'
241     print *,' GSOEQ=',GSOEQ,' KSTREND=',KSTREND
242     c stop
243     endif
244     if(YEARGT.ne.1860)then
245     print *,'Wrong setting of GSOEQ and YEARGT'
246     print *,' GSOEQ=',GSOEQ,' YEARGT=',YEARGT
247     stop
248     endif
249     C For an equilibrium run for 1860 conditions for GSO runs
250     print *,
251     & 'An equilibrium run for 1860 conditions for GSO runs'
252     S0X=1365.3596/1367.
253     print *,'Sconstant=',S0X*1367.
254     RVOL=0.012
255     FVOL=0.0045
256     FGOLDU(1)=(RVOL+FVOL)/RVOL
257     print *,'Additional STAER=',FVOL,
258     & ' FGOLDU(1)=',FGOLDU(1)
259     C
260     endif
261     if(CO2FOR)then
262     print *,' Forcing due to ',CO2F,'xCO2 is calculated',
263     * 'in stead of cloud forcing'
264     endif
265     if(S0FOR)then
266     print *,' Forcing due to change in solar constant ',
267     * 'is calculated in stead of cloud forcing'
268     endif
269     if(FORBC)then
270     print *,' Forcing due to black carbon ',
271     * 'is calculated in stead of cloud forcing'
272     endif
273     if(STRARFOR)then
274     RVOL=0.012
275     FVOL=FVOLADD
276     FGOLDU1=(RVOL+FVOL)/RVOL
277     print *,'Additional STAER=',FVOL,
278     & ' FGOLDU1=',FGOLDU1
279     endif
280     if(GHSFALB)then
281     Print *,' Run with ghost albedo forcing'
282     print *,'ALBCF=',ALBCF
283     endif
284     if(GHSF)then
285     Print *,' Run with ghost forcing'
286     do l=1,lm+1
287     ghostfv(l)=0.0
288     do j=1,jm
289     ghostfv(l)=ghostfv(l)+ghostf(l,j)
290     enddo
291     ghostfv(l)=ghostfv(l)/float(jm)
292     enddo
293     print *,' Ghost forcing DFS=',ghostfv(1)
294     do l=LM+1,2,-1
295     print *,' Forsing at laeyr',l-1,'=',ghostfv(l)
296     enddo
297     c print *,' 2D ghost forcing'
298     endif
299     S0X=CFS0X*S0X
300     S0X0=S0X
301     print *,'Before CALL RADIA0'
302     CALL RADIA0 (IO,JM,CO2,READGHG) 5090.
303     print *,'After CALL RADIA0'
304     INCHM=NRAD/NDYN 5091.
305     C**** CLOUD LAYER INDICES USED FOR DIAGNOSTICS 5092.
306     DO 43 L=1,LM 5093.
307     LLOW=L 5094.
308     IF (.5*(PLE(L+1)+PLE(L+2)).LT.786.) GO TO 44 5095.
309     43 CONTINUE 5096.
310     44 LMID1=LLOW+1 5097.
311     DO 45 L=LMID1,LM 5098.
312     LMID=L 5099.
313     IF (.5*(PLE(L+1)+PLE(L+2)).LT.430.) GO TO 46 5100.
314     45 CONTINUE 5101.
315     46 LHI1=LMID+1 5102.
316     LHI=LM 5103.
317     IF (LHI1.GT.LHI) LHI=LHI1 5104.
318     WRITE (6,47) LLOW,LMID1,LMID,LHI1,LHI 5105.
319     47 FORMAT (' LOW CLOUDS IN LAYERS 1-',I2,' MID LEVEL CLOUDS IN',5106.
320     * ' LAYERS',I3,'-',I2,' HIGH CLOUDS IN LAYERS',I3,'-',I2) 5107.
321     C**** NO RADIATION AVERAGING IJRA=1 JRA=1 IRA=1 5108.
322     C**** RADIATION AVERAGING IN I 2 1 2 5109.
323     C**** RADIATION AVERAGING IN I AND J 4 2 2 5110.
324     JRA=(IJRA+2)/3 5111.
325     IRA=IJRA/JRA 5112.
326     50 JALTER=MOD(NSTEP,NRAD*JRA)/NRAD 5113.
327     JDAYR=JDAY
328     JYEARR=JYEAR
329     IALTER=MOD(NSTEP,NRAD*IJRA)/(NRAD*JRA) 5114.
330     TNOW=JYEAR+(JDAY-.5)/365. 5127.1
331     IF (JDAY.NE.JDLAST.AND.KSTREND.NE.0) then
332     IF (KSTREND.GT.0) then
333     if(abs(CFS0X-1.0).gt.1.e3)then
334     print *,'Wrong setting KSTREND=',KSTREND,' CFS0X=',
335     & CFS0X
336     endif
337     S0X=S0X0*(1.+S0RATE/70.*(TNOW-INYEAR))
338     c print *,'From radia S0X0=',S0X0,' S0X=',S0X
339     c print *,'S0X=',S0X,' TNOW=',TNOW
340     ELSE
341     call obssolar(S0X,TNOW)
342     ENDIF
343     ENDIF
344     S0=S0X*1367./RSDIST 5115.
345     S0R=S0
346     c print *,'S0=',S0
347     C**** CALCULATE AVERAGE COSINE OF ZENITH ANGLE FOR CURRENT COMP3 STEP 5116.
348     C**** AND RADIATION PERIOD 5117.
349     ROT1=TWOPI*TOFDAY/24. 5118.
350     if(DC25)then
351     ROT2=ROT1+TWOPI*DTCNDS/SDAY 5119.
352     CALL COSZT (IO,JM,SIND,COSD,ROT1,ROT2,COSZ1) 5120.
353     else
354     ROT2=ROT1+TWOPI
355     CALL COSZR (IO,JM,SIND,COSD,ROT1,ROT2,COSZ1)
356     end if
357     if(HPRNT)then
358     print *,' radia TAU=',TAU
359     print *,' CLDSS'
360     print *,(CLDSS(1,7,L),L=1,LM)
361     print *,' CLDMC'
362     print *,(CLDMC(1,7,L),L=1,LM)
363     endif
364     cprint *,' form radia TAU=',TAU,'MODRD=',MODRD
365     C
366     IF(MODRD.NE.0) GO TO 840 5121.
367     C
368     ROT2=ROT1+TWOPI*NRAD*DT/SDAY 5122.
369     CALL COSZS (IO,JM,SIND,COSD,ROT1,ROT2,COSZ2,COSZA) 5123.
370     C**** 5124.
371     C**** COMPUTE EARTH ALBEDOS AND OTHER PARAMETERS FOR BEGINNING OF DAY 5125.
372    
373     TNOW=JYEAR+(JDAY-.5)/365. 5127.1
374     KWRITE=0
375     if(JMONTH.ne.JMLAST) then
376     KWRITE=1
377     print *,' FROM radia READGHG=',READGHG
378     if(READGHG.eq.2) call tgases(CO2,JMONTH)
379     c print *,(zco2(1,j,1)
380     c & *28.97296245/44.0*1.e-3,j=1,jm)
381     !ppb(m) to kg per volume base
382     if(READGHG.eq.1) call rtgases(CO2,JMONTH)
383     #if ( defined PREDICTED_BC)
384     nrdbc=nrdbc+1
385     do l=1,lm
386     do j=1,jm
387     bcod(1,j,l)=bcodmn(j,l,nrdbc)
388     enddo
389     enddo
390     if(nrdbc.eq.12)then
391     nrdbc=0.
392     endif
393     #endif
394     endif
395     JMLAST=JMONTH
396     IF (JDAY.NE.JDLAST.AND.KTREND.GT.0)
397     & CALL FORGET(TNOW,KTREND,KWRITE)
398     IF (JDAY.NE.JDLAST)then
399     #ifdef PREDICTED_GASES
400     call chemglobal(P)
401     #endif
402     #if ( defined PREDICTED_BC)
403     call sulfr(BSO4LAND,BSO4OCEAN,1986.)
404     do j=1,jm
405     FLAND=FDATA(1,J,2)
406     BSO4TOTAL(j)=BSO4LAND(j)*FLAND+BSO4OCEAN(j)*(1.-FLAND)
407     enddo
408     #endif
409     CALL RCOMPT
410     if(CLDFEED)then
411     DTSURFAV=0.
412     C
413     do j=1,jm
414     DTSURFAV=DTSURFAV+DTSURF(J)*DXYP(j)
415     end do !j
416     DTSURFAV=DTSURFAV/AREAG
417     do j=1,jm
418     do k=1,3
419     cfcld(j,k)=1.+coefcl(k)*DTSURFAV
420     end do ! k
421     end do ! j
422     endif
423     ENDIF
424     JDLAST=JDAY 5129.
425     IHOUR=1.5+TOFDAY 5130.
426     if(CLDFEED)then
427     if (KWRITE.eq.1)then
428     print *,'cfcld'
429     print 9456,cfcld
430     print *,' DTSURF'
431     print 9456,DTSURF
432     print *,' DTSURFAV=',DTSURFAV
433     9456 format(12f6.2)
434     endif
435     C Save clouds predicted by model
436     do k=1,LM
437     do j=1,JM
438     CLDSSF(j,k)=CLDSS(1,j,k)
439     CLDMCF(j,k)=CLDMC(1,j,k)
440     enddo
441     enddo
442     do k=1,LM
443     if(k.le.5)then
444     k1=1
445     else
446     k1=2
447     endif
448     do j=1,JM
449     CLDSS(1,j,k)=cfcld(j,k1)*CLDSS(1,j,k)
450     CLDMC(1,j,k)=cfcld(j,3)*CLDMC(1,j,k)
451     enddo
452     enddo
453     endif
454     C
455     CB READING OF CLOUD
456     if(abs(PCLOUD-3.).lt.1.5)then
457     910 continue
458     if(nreadcld.eq.nrcldmax)go to 900
459     read(585,END=900)TFDAYF,JDATEF,JMNTHF,CLDSSF,CLDMCF,IRAND
460     nreadcld=nreadcld+1
461     if(IFIRST.eq.1)then
462     print *,' radia.f PCLOUD=',PCLOUD
463     if(PCLOUD.eq.2)print *,' FIXED MC and SS CLOUDS'
464     if(PCLOUD.eq.4)print *,' FIXED MC CLOUDS ONLY'
465     if(PCLOUD.eq.3)print *,' FIXED SS CLOUDS ONLY'
466     print *,TOFDAY,JDATE,JMONTH
467     print *,TFDAYF,JDATEF,JMNTHF
468     print *,' DTCNDS=',DTCNDS/3600.
469     print *,' DT*NRAD=',DT*NRAD/3600.
470     if(.not.WRCLD)IFIRST=0
471     endif
472     if(abs(TOFDAY-TFDAYF).gt.1.e-3.or.JDATE.ne.JDATEF.or.
473     * JMONTH.ne.JMNTHF)then
474     print *,' RADIA, disagrement in clouds'
475     print *,TOFDAY,JDATE,JMONTH
476     print *,TFDAYF,JDATEF,JMNTHF
477     stop
478     endif
479     go to 920
480     900 rewind 585
481     nreadcld=0
482     print *,' END OF file85'
483     print *,JYEAR
484     print *,TOFDAY,JDATE,JMONTH
485     print *,' REWIND 85'
486     go to 910
487     920 continue
488     c if (KWRITE.eq.1)then
489     c print *,'cfcld'
490     c print 9456,cfcld
491     c print *,' DTSURF'
492     c print 9456,DTSURF
493     c9456 format(12f6.2)
494     c endif
495     CALL RINIT (IRAND)
496     do 930 k=1,LM
497     do 930 j=1,JM
498     if(PCLOUD.ne.4)CLDSS(1,j,k)=CLDSSF(j,k)
499     if(PCLOUD.ne.3)CLDMC(1,j,k)=CLDMCF(j,k)
500     930 continue
501     endif
502     CE END OF READING OF CLOUD
503     if(WRCLD)then
504     if(NWRCLD.eq.1)then
505     CALL RFINAL(IRAND)
506     if(IFIRST.eq.1)print *,' SHORT CLOUDS RECORD'
507     write(581)TOFDAY,JDATE,JMONTH,CLDSS,CLDMC,IRAND
508     elseif(NWRCLD.eq.2)then
509     if(IFIRST.eq.1)print *,' LONG CLOUDS RECORD'
510     do 1150 k=1,14
511     do 1150 j=1,JM0
512     if(k.le.2)then
513     ODATA2(j,k)=ODATA(1,j,k)
514     BDATA2(j,k)=BLDATA(1,j,k)
515     FDATA2(j,k)=FDATA(1,j,k+1)
516     endif
517     if(k.le.3)RQT2(j,k)=RQT(1,j,k)
518     GDATA2(j,k)=GDATA(1,j,k)
519     1150 continue
520     CALL RFINAL(IRAND)
521     write(581)TOFDAY,JDATE,JMONTH,CLDSS,CLDMC,IRAND,
522     * JDAY,JYEAR,T,Q,P,
523     * ODATA2,BDATA2,FDATA2,GDATA2,RQT2
524     & ,CLDSSF,CLDMCF
525     else
526     print *,' NWRCLD=',NWRCLD
527     stop
528     endif
529     IFIRST=0
530     endif
531     C**** 5131.
532     C**** MAIN J LOOP 5132.
533     C**** 5133.
534     DO 600 J=1,JM 5134.
535     IF ((J-1)*(JM-J).NE.0) GO TO 140 5135.
536     C**** CONDITIONS AT THE POLES 5136.
537     POLE=.TRUE. 5137.
538     MODRJ=0 5138.
539     IMAX=1 5139.
540     GO TO 160 5140.
541     C**** CONDITIONS AT NON-POLAR POINTS 5141.
542     140 POLE=.FALSE. 5142.
543     MODRJ=MOD(J+JALTER,JRA) 5143.
544     IMAX=IM 5144.
545     160 XFRADJ=.2+1.2*COSP(J)*COSP(J) 5145.
546     JLAT=J 5145.1
547     if(GSOEQ)then
548     C
549     RVOL=0.012
550     FVOL=0.0045
551     FGOLDU(1)=(RVOL+FVOL)/RVOL
552     C
553     endif
554     IF(MODRJ.EQ.0) CALL RCOMPJ 5146.
555     C**** 5147.
556     C**** MAIN I LOOP 5148.
557     C**** 5149.
558     IM1=IM 5150.
559     DO 500 I=1,IMAX 5151.
560     MODRIJ=MODRJ+MOD(I+IALTER,IRA) 5152.
561     IF(POLE) MODRIJ=0 5153.
562     JR=J
563     C**** DETERMINE FRACTIONS FOR SURFACE TYPES AND COLUMN PRESSURE 5155.
564     PLAND=FDATA(I,J,2) 5156.
565     PWATER=1.-PLAND
566     POICE=ODATA(I,J,2)*(1.-PLAND) 5157.
567     POCEAN=(1.-PLAND)-POICE 5158.
568     if(POCEAN.LE.1.E-5)then
569     POCEAN=0.
570     POICE=PWATER
571     endif
572     PLICE=FDATA(I,J,3)*PLAND 5159.
573     PEARTH=PLAND-PLICE 5160.
574     SP=P(I,J) 5161.
575     C**** 5162.
576     C**** DETERMINE CLOUDS (AND THEIR OPTICAL DEPTHS) SEEN BY RADIATION 5163.
577     C**** 5164.
578     X=999999. 5164.1
579     c RANDSS=RANDU(X) 5165.
580     c RANDMC=RANDU(X) 5166.
581     CALL RANDUU(RANDSS,X)
582     CALL RANDUU(RANDMC,X)
583     if(HPRNT)then
584     print *,' radia TAU=',TAU
585     print *,' RANDSS=',RANDSS
586     print *,' RANDMC=',RANDMC
587     endif
588     C
589     CSS=0. 5167.
590     CMC=0. 5168.
591     DEPTH=0. 5169.
592     LTOP=0 5169.1
593     DO 210 L=1,LM 5170.
594     RTAU(L)=0. 5171.
595     210 TOTCLD(L)=0. 5172.
596     DO 240 L=1,LM 5173.
597     IF(CLDSS(I,J,L).LT.RANDSS) GO TO 220 5174.
598     RTAUSS=.013333*(PTOP-100.+SIG(L)*SP) 5175.
599     IF(RTAUSS.LT.0.) RTAUSS=0. 5176.
600     IF (T(I,J,L)*PK(I,J,L).LT.TCIR) RTAUSS=.3333333 5177.
601     RTAU(L)=RTAUSS 5178.
602     CSS=1. 5179.
603     AJL(J,L,28)=AJL(J,L,28)+CSS 5180.
604     TOTCLD(L)=1. 5181.
605     LTOP=L 5181.1
606     220 IF(CLDMC(I,J,L).LE.RANDMC) GO TO 240 5182.
607     RTAUMC=DSIG(L)*SP*.08 5183.
608     IF(RTAUMC.GT.RTAU(L)) RTAU(L)=RTAUMC 5184.
609     CMC=1. 5185.
610     AJL(J,L,29)=AJL(J,L,29)+CMC 5186.
611     TOTCLD(L)=1. 5187.
612     LTOP=L 5187.1
613     DEPTH=DEPTH+SP*DSIG(L) 5188.
614     AJL(J,L,19)=AJL(J,L,19)+TOTCLD(L) 5189.
615     240 CONTINUE
616     AJ(J,57)=AJ(J,57)+CSS*POCEAN 5190.
617     BJ(J,57)=BJ(J,57)+CSS*PLAND 5191.
618     CJ(J,57)=CJ(J,57)+CSS*POICE 5192.
619     DJ(JR,57)=DJ(JR,57)+CSS*DXYP(J) 5193.
620     AJ(J,58)=AJ(J,58)+CMC*POCEAN 5194.
621     BJ(J,58)=BJ(J,58)+CMC*PLAND 5195.
622     CJ(J,58)=CJ(J,58)+CMC*POICE 5196.
623     DJ(JR,58)=DJ(JR,58)+CMC*DXYP(J) 5197.
624     AIJ(I,J,17)=AIJ(I,J,17)+CMC 5198.
625     AJ(J,80)=AJ(J,80)+DEPTH*POCEAN 5199.
626     BJ(J,80)=BJ(J,80)+DEPTH*PLAND 5200.
627     CJ(J,80)=CJ(J,80)+DEPTH*POICE 5201.
628     DJ(JR,80)=DJ(JR,80)+DEPTH*DXYP(J) 5202.
629     CLDCV=CMC+CSS-CMC*CSS 5203.
630     cldd4tem(j)=cldd4tem(j)+CLDCV
631     ncldd4tem(j)=ncldd4tem(j)+1
632     AJ(J,59)=AJ(J,59)+CLDCV*POCEAN 5204.
633     BJ(J,59)=BJ(J,59)+CLDCV*PLAND 5205.
634     CJ(J,59)=CJ(J,59)+CLDCV*POICE 5206.
635     DJ(JR,59)=DJ(JR,59)+CLDCV*DXYP(J) 5207.
636     C**** clear sky condinion
637     CSS0=CSS
638     CMC0=CMC
639     do L=1,LM
640     RTAU0(L)=RTAU(L)
641     TOTCLD0(L)=TOTCLD(L)
642     enddo
643     if(CMC.le.0.and.CSS.le.0)then
644     CLEAR(J)=1
645     else
646     CLEAR(J)=0
647     endif
648     if(STRARFOR.or.CO2FOR.or.S0FOR.or.FORBC)then
649     CLEAR(J)=0
650     endif
651     C**** 5247.
652     C**** SET UP VERTICAL ARRAYS OMITTING THE I AND J INDICES 5248.
653     C**** 5249.
654     C**** EVEN PRESSURES 5250.
655     DO 340 L=1,LM 5251.
656     PLE(L)=SIGE(L)*SP+PTOP 5252.
657     C**** TEMPERATURES 5253.
658     TL(L)=T(I,J,L)*PK(I,J,L) 5254.
659     C**** MOISTURE VARIABLES 5255.
660     QL(L)=Q(I,J,L) 5256.
661     340 CONTINUE 5257.
662     C**** 5258.
663     C**** RADIATION, SOLAR AND THERMAL 5259.
664     C**** 5260.
665    
666     DO 420 K=1,3 5261.
667     420 TL(LM+K)=RQT(I,J,K) 5262.
668     COSZ=COSZA(I,J) 5263.
669     TGO=ODATA(I,J,1)+TF 5264.
670     TGOI=GDATA(I,J,3)+TF 5265.
671     TGLI=GDATA(I,J,13)+TF 5266.
672     TGE=GDATA(I,J,4)+TF 5267.
673     TS=BLDATA(I,J,2) 5268.
674     SNOWOI=GDATA(I,J,1) 5269.
675     SNOWLI=GDATA(I,J,12) 5270.
676     SNOWE=GDATA(I,J,2) 5271.
677     AGESN=GDATA(I,J,11) 5272.
678     c print *,'From radia J=',J,' TAU=',TAU,' TS=',TS
679     c print *,TGO,TGOI,TGLI,TGE
680     c print *,SNOWOI,SNOWLI,SNOWE,AGESN
681     WEARTH=(GDATA(I,J,5)+GDATA(I,J,6))/(VDATA(I,J,9)+1.E-20) 5273.
682     DO 430 K=1,8 5274.
683     430 PVT(K)=VDATA(I,J,K) 5275.
684     WS=BLDATA(I,J,1) 5276.
685     do 439 L=1,LM+1
686     SRHR(I,J,L)=0.
687     TRHR(I,J,L)=0.
688     if(L.le.4)then
689     SNFS(I,J,L)=0.
690     TNFS(I,J,L)=0.
691     if(L.le.3)then
692     SRHRS(I,J,L)=0.
693     TRHRS(I,J,L)=0.
694     endif
695     endif
696     439 continue
697     TRNFP0(J)=0.
698     TRNFP1(J)=0.
699     TRINCG(I,J)=0.
700     BTMPW(I,J)=0.
701     SRDAN=0.
702     SRNAN=0.
703     do 449 K=1,9
704     ALB(I,J,K)=0.
705     ALBJ(J,K)=0.
706     449 continue
707     do 499 ii=1,3
708     COSZ=COSZA(I,J)
709     PLAND=FDATA(I,J,2)
710     PWATER=1.-PLAND
711     POICE=ODATA(I,J,2)*(1.-PLAND)
712     POCEAN=(1.-PLAND)-POICE
713     if(POCEAN.LE.1.E-5)then
714     POCEAN=0.
715     POICE=PWATER
716     endif
717     PLICE=FDATA(I,J,3)*PLAND
718     PEARTH=PLAND-PLICE
719     if(ii.eq.1)then
720     BSO4BC=BSO4OCEAN(J)/BSO4TOTAL(J)
721     PTYPE=POCEAN
722     POICE=0.
723     POCEAN=1.
724     PLAND=0.
725     PEARTH=0.
726     PLICE=0.
727     TGAL=0.
728     else if(ii.eq.3)then
729     BSO4BC=BSO4OCEAN(J)/BSO4TOTAL(J)
730     PTYPE=POICE
731     POICE=1.
732     POCEAN=0.
733     PLAND=0.
734     PEARTH=0.
735     PLICE=0.
736     TGAL=TGOI
737     else
738     BSO4BC=BSO4LAND(J)/BSO4TOTAL(J)
739     PTYPE=PLAND
740     POCEAN=0.
741     POICE=0.
742     PWATER=0.
743     PLICE=FDATA(I,J,3)
744     PEARTH=1.-PLICE
745     TGAL=TGE*PEARTH+TGLI*PLICE
746     PLAND=1.
747     endif
748     if(PTYPE.lt.1.e-10)go to 499
749     if(ii.gt.1)then
750     if(TGAL.lt.263.)then
751     FRSNALB=0.30
752     elseif(TGAL.lt.273.)then
753     FRSNALB=0.30-0.015*(TGAL-263.)
754     else
755     FRSNALB=0.15
756     endif
757     endif !ii
758     c FGOLDU(2)=XFRADJ*(1.-PEARTH) 5277.
759     c FGOLDU(3)=XFRADJ*PEARTH 5278.
760     FGOLDU(2)=XFRADJ*(1.-PLAND)
761     FGOLDU(3)=XFRADJ*PLAND
762     if(CO2FOR)then
763     IF (KTREND.GT.0) THEN
764     CALL FORGET(TNOW,KTREND,0)
765     ELSE
766     FULGAS(2)=CO2
767     ENDIF
768     endif
769     !
770     #if ( defined FIXED_FOR )
771     FULGAS(2)=1.
772     #endif
773     !
774     if(S0FOR)then
775     S0=S0R
776     endif
777     CSS=CSS0
778     CMC=CMC0
779     do L=1,LM
780     RTAU(L)=RTAU0(L)
781     TOTCLD(L)=TOTCLD0(L)
782     enddo
783    
784     Chemstry Model Patch 092295
785     c
786     PTYPER=max(PTYPE,5.e-3)
787     FAERSOL=PLAND/PTYPER*3.0*CFAEROSOL
788     C Factor 3 is a "parameterization" for indirec aerosol effect
789     c For normal setting CFAEROSOL=1.0 (default)
790     c For hight forcing CFAEROSOL=0.5
791     c For low forcing CFAEROSOL=2.0
792     c
793     c ==========================
794     #if ( defined PREDICTED_BC)
795     FBC=BSO4BC*CFBC
796     #endif
797     ILON=I 5278.1
798     JLAT=J 5278.2
799     if(GHSFALB)then
800     BVSURFA=ALBCF
801     XVSURFA=ALBCF
802     BNSURFA=ALBCF
803     XNSURFA=ALBCF
804     else
805     BVSURFA=0.0
806     XVSURFA=0.0
807     BNSURFA=0.0
808     XNSURFA=0.0
809     endif
810     if(STRARFOR)then
811     c if(j.gt.JM/2)then
812     FGOLDU(1)=FGOLDU1
813     c else
814     c FGOLDU(1)=1.
815     c endif
816     endif
817     if(J.eq.-10)then
818     print *,'From radia GHSFALB=',GHSFALB,' ALBCF=',ALBCF
819     print *,'From radia BVSURFA=',BVSURFA
820     endif
821     c if(J.eq.1)then
822     c print *,ii,PTYPE,PWATER,POCEAN,POICE,PLAND,PEARTH,PLICE
823     c print *,(FGOLDU(if),if=1,5)
824     c endif
825     c print *,'S0 for real calculation=',S0
826     CALL RCOMPX 5279.
827     if (IRFIRST.eq.1.and.READGHG.eq.1)then
828     CALL WRITER(12)
829     if(ii.ge.2)IRFIRST=0
830     endif
831     IF(DMOD(TAU,365.*24.).EQ.0..and.J.eq.JM/2) then
832     print *,' tau=',TAU,' J=',J
833     CALL WRITER (1,0)
834     endif
835     c CALL WRITER (13) 5280.
836     SRHR(I,J,1)=SRHR(I,J,1)+SRNFLB(1)*PTYPE
837     c GHOST FORCINGS
838     if(GHSF)then
839    
840     do LFR=1,lm+1
841     if(LFR.gt.1)then
842     TRFCRL(LFR-1)=TRFCRL(LFR-1)-ghostf(LFR,J)
843     TRNFLB(LFR)=TRNFLB(LFR)+ghflux(LFR,J)
844     endif
845     enddo
846     endif
847     c
848     TRHR(I,J,1)=TRHR(I,J,1)+(STBO*(POCEAN*TGO**4+POICE*TGOI**4
849     * +PLICE*TGLI**4+PEARTH*TGE**4)-TRNFLB(1))*PTYPE
850     if(GHSF)then
851     TRHR(I,J,1)=TRHR(I,J,1)+ghostf(1,J)*PTYPE
852     TRDFLB(1)=TRDFLB(1)+ghostf(1,J)
853     endif
854     C *****
855     TRSURF(J,ii)=STBO*(POCEAN*TGO**4+POICE*TGOI**4
856     * +PLICE*TGLI**4+PEARTH*TGE**4)-TRNFLB(1)
857     if(GHSF)then
858     TRSURF(J,ii)=TRSURF(J,ii)+ghostf(1,J)
859     endif
860     SRSURF(J,ii)=SRNFLB(1)
861     DO 440 L=1,LM 5284.
862     SRHR(I,J,L+1)=SRHR(I,J,L+1)+SRFHRL(L)*PTYPE
863     440 TRHR(I,J,L+1)=TRHR(I,J,L+1)-TRFCRL(L)*PTYPE
864     DO 450 LR=1,3 5287.
865     SRHRS(I,J,LR)=SRHRS(I,J,LR)+SRFHRL(LM+LR)*PTYPE
866     450 TRHRS(I,J,LR)=TRHRS(I,J,LR)-TRFCRL(LM+LR)*PTYPE
867     DO 460 K=1,4 5290.
868     SNFS(I,J,K)=SNFS(I,J,K)+SRNFLB(K+LM)*PTYPE
869     460 TNFS(I,J,K)=TNFS(I,J,K)+(TRNFLB(K+LM)-TRNFLB(1))*PTYPE
870     TRNFP0(J)=TRNFP0(J)+TRNFLB(4+LM)*PTYPE
871     TRNFP1(J)=TRNFP1(J)+TRNFLB(1+LM)*PTYPE
872     TRINCG(I,J)=TRINCG(I,J)+TRDFLB(1)*PTYPE
873     BTMPW(I,J)=BTMPW(I,J)+(BTEMPW-TF)*PTYPE
874     c ALB(I,J,1)=SRNFLB(1)/(SRDFLB(1)+1.E-20) 5295.
875     SRDAN=SRDAN+SRDFLB(1)*PTYPE
876     SRNAN=SRNAN+SRNFLB(1)*PTYPE
877     ALB(I,J,2)=ALB(I,J,2)+PLAVIS*PTYPE
878     ALB(I,J,3)=ALB(I,J,3)+PLANIR*PTYPE
879     ALB(I,J,4)=ALB(I,J,4)+ALBVIS*PTYPE
880     ALB(I,J,5)=ALB(I,J,5)+ALBNIR*PTYPE
881     ALB(I,J,6)=ALB(I,J,6)+SRRVIS*PTYPE
882     ALB(I,J,7)=ALB(I,J,7)+SRRNIR*PTYPE
883     ALB(I,J,8)=ALB(I,J,8)+SRAVIS*PTYPE
884     ALB(I,J,9)=ALB(I,J,9)+SRANIR*PTYPE
885     ALB1=SRNFLB(1)/(SRDFLB(1)+1.E-20)
886     C **********
887     ALBJ(J,2)=PLAVIS
888     ALBJ(J,3)=PLANIR
889     ALBJ(J,4)=ALBVIS
890     ALBJ(J,5)=ALBNIR
891     ALBJ(J,6)=SRRVIS
892     ALBJ(J,7)=SRRNIR
893     ALBJ(J,8)=SRAVIS
894     ALBJ(J,9)=SRANIR
895     ALBJ(J,1)=SRNFLB(1)/(SRDFLB(1)+1.E-20)
896     C *********
897     COSZ=COSZ2(I,J)
898     do L=LM,1,-1
899     AJL(J,L,43)=AJL(J,L,43)-TRNFLB(L+1)*PTYPE
900     AJL(J,L,42)=AJL(J,L,42)+(SRNFLB(L+1)*COSZ)*PTYPE
901     enddo
902     if(ii.eq.2)then
903     #if ( defined CLM )
904     C for TEM CLM
905     DSWSRF(j)=SRDFLB(1)
906     DLWSRF(j)=TRSURF(J,2)
907     DSWVIS(j)=SRDVIS
908     DSWNIR(j)=SRDNIR
909     C for TEM CLM
910     #endif
911     PLAND=PTYPE
912     BJ(J,1)=BJ(J,1)+(S0*COSZ)*PLAND
913     BJ(J,2)=BJ(J,2)+(SRNFLB(4+LM)*COSZ)*PLAND
914     c BSNFS1=BSNFS1+(SRNFLB(1+LM)*COSZ)*PLAND
915     c BJ(J,5)=BJ(J,5)+(SRNFLB(1)*COSZ/(ALB1+1.E-20))*PLAND
916     BJ(J,5)=BJ(J,5)+(SRDFLB(1)*COSZ)*PLAND
917     BJ(J,6)=BJ(J,6)+(SRNFLB(1)*COSZ)*PLAND
918     BJ(J,55)=BJ(J,55)+(BTEMPW-TF)*PLAND
919     BJ(J,67)=BJ(J,67)+TRDFLB(1)*PLAND
920     BJ(J,70)=BJ(J,70)-(TRNFLB(4+LM)-TRNFLB(1))*PLAND
921     BJ(J,7)=BJ(J,7)-TRNFLB(4+LM)*PLAND
922     BJ(J,8)=BJ(J,8)-TRNFLB(1+LM)*PLAND
923     BJ(J,3)=BJ(J,3)+(SRNFLB(1+LM)*COSZ)*PLAND
924     BJ(J,71)=BJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*PLAND
925     DO 761 K=2,9
926     BJ(J,K+70)=BJ(J,K+70)+(S0*COSZ)*ALBJ(J,K)*PLAND
927     761 CONTINUE
928     else if(ii.eq.1)then
929     POCEAN=PTYPE
930     AJ(J,1)=AJ(J,1)+(S0*COSZ)*POCEAN
931     AJ(J,2)=AJ(J,2)+(SRNFLB(4+LM)*COSZ)*POCEAN
932     AJ(J,5)=AJ(J,5)+(SRDFLB(1)*COSZ)*POCEAN
933     AJ(J,6)=AJ(J,6)+(SRNFLB(1)*COSZ)*POCEAN
934     AJ(J,55)=AJ(J,55)+(BTEMPW-TF)*POCEAN
935     AJ(J,67)=AJ(J,67)+TRDFLB(1)*POCEAN
936     AJ(J,70)=AJ(J,70)-(TRNFLB(4+LM)-TRNFLB(1))*POCEAN
937     AJ(J,7)=AJ(J,7)-TRNFLB(4+LM)*POCEAN
938     AJ(J,8)=AJ(J,8)-TRNFLB(1+LM)*POCEAN
939     AJ(J,3)=AJ(J,3)+(SRNFLB(1+LM)*COSZ)*POCEAN
940     AJ(J,71)=AJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POCEAN
941     #if ( defined OCEAN_3D )
942     solarinc_ocean(J)=solarinc_ocean(J)+SRDFLB(1)*COSZ
943     solarnet_ocean(J)=solarnet_ocean(J)+SRNFLB(1)*COSZ
944     navrado(j)=navrado(j)+1
945     #endif
946     C
947     DO K=2,9
948     AJ(J,K+70)=AJ(J,K+70)+(S0*COSZ)*ALBJ(J,K)*POCEAN
949     END DO
950     C
951     else
952     POICE=PTYPE
953     CJ(J,1)=CJ(J,1)+(S0*COSZ)*POICE
954     CJ(J,2)=CJ(J,2)+(SRNFLB(4+LM)*COSZ)*POICE
955     CJ(J,5)=CJ(J,5)+(SRDFLB(1)*COSZ)*POICE
956     CJ(J,6)=CJ(J,6)+(SRNFLB(1)*COSZ)*POICE
957     CJ(J,55)=CJ(J,55)+(BTEMPW-TF)*POICE
958     CJ(J,67)=CJ(J,67)+TRDFLB(1)*POICE
959     CJ(J,70)=CJ(J,70)-(TRNFLB(4+LM)-TRNFLB(1))*POICE
960     CJ(J,7)=CJ(J,7)-TRNFLB(4+LM)*POICE
961     CJ(J,8)=CJ(J,8)-TRNFLB(1+LM)*POICE
962     CJ(J,3)=CJ(J,3)+(SRNFLB(1+LM)*COSZ)*POICE
963     CJ(J,71)=CJ(J,71)-(TRNFLB(1+LM)-TRNFLB(1))*POICE
964     #if ( defined OCEAN_3D )
965     solarinc_ice(J)=solarinc_ice(J)+SRDFLB(1)*COSZ
966     solarnet_ice(J)=solarnet_ice(J)+SRNFLB(1)*COSZ
967     navrad(j)=navrad(j)+1
968     #endif
969     C
970     DO K=2,9
971     CJ(J,K+70)=CJ(J,K+70)+(S0*COSZ)*ALBJ(J,K)*POICE
972     END DO
973     endif
974     if(CLEAR(J).eq.0)then
975     PLAND=FDATA(I,J,2)
976     PWATER=1.-PLAND
977     POICE=ODATA(I,J,2)*(1.-PLAND)
978     POCEAN=(1.-PLAND)-POICE
979     if(POCEAN.LE.1.E-5)then
980     POCEAN=0.
981     POICE=PWATER
982     endif
983     PLICE=FDATA(I,J,3)*PLAND
984     PEARTH=PLAND-PLICE
985     if(ii.eq.1)then
986     PTYPE=POCEAN
987     POICE=0.
988     POCEAN=1.
989     PLAND=0.
990     PEARTH=0.
991     PLICE=0.
992     else if(ii.eq.3)then
993     PTYPE=POICE
994     POICE=1.
995     POCEAN=0.
996     PLAND=0.
997     PEARTH=0.
998     PLICE=0.
999     else
1000     PTYPE=PLAND
1001     POCEAN=0.
1002     POICE=0.
1003     PWATER=0.
1004     PLICE=FDATA(I,J,3)
1005     PEARTH=1.-PLICE
1006     PLAND=1.
1007     endif
1008     COSZ=COSZA(I,J)
1009     if(STRARFOR)then
1010     FGOLDU(1)=1.0
1011     elseif(CO2FOR)then
1012     FULGAS(2)=CO2F
1013     elseif(S0FOR)then
1014     c print *,'For S0FOR ',S0X0,CFS0X,S0X
1015     S0=S0X0/CFS0X*1367./RSDIST
1016     c print *,'S0 for forcing calculation=',S0
1017     elseif(FORBC)then
1018     FBC=0.0
1019     else
1020     CSS=0.
1021     CMC=0.
1022     DEPTH=0.
1023     LTOP=0.
1024     do 1210 L=1,LM
1025     RTAU(L)=0.
1026     TOTCLD(L)=0.
1027     1210 continue
1028     endif
1029     568 continue
1030     CALL RCOMPX
1031     endif
1032     SRHRCL(J)=SRNFLB(1)
1033     TRHRCL(J)=-TRNFLB(1)
1034     ALBCL(J)=SRNFLB(1)/(SRDFLB(1)+1.e-20)
1035     SNP1CL(J)=SRNFLB(LM+1)
1036     SNP0CL(J)=SRNFLB(LM+4)
1037     TRINCL(J)=TRDFLB(1)
1038     TRP0CL(J)=TRNFLB(LM+4)
1039     TRP1CL(J)=TRNFLB(LM+1)
1040     COSZ=COSZ2(I,J)
1041     do L=LM,1,-1
1042     AJL(J,L,46)=AJL(J,L,46)-TRNFLB(L+1)*PTYPE
1043     AJL(J,L,45)=AJL(J,L,45)+(SRNFLB(L+1)*COSZ)*PTYPE
1044     enddo
1045     C *******
1046     NCLR(J)=NCLR(J)+1
1047     SRHE5=SRHRCL(J)*COSZ/(ALBCL(J)+1.E-20)
1048     if(ii.eq.1)then
1049     POCEAN=PTYPE
1050     AJCLR(J,1)=AJCLR(J,1)+(S0*COSZ)*POCEAN
1051     AJCLR(J,2)=AJCLR(J,2)+(SNP0CL(J)*COSZ)*POCEAN
1052     AJCLR(J,4)=AJCLR(J,4)+(SRHRCL(J)*COSZ)*POCEAN
1053     AJCLR(J,5)=AJCLR(J,5)+SRHE5*POCEAN
1054     AJCLR(J,6)=AJCLR(J,6)+TRINCL(J)*POCEAN
1055     AJCLR(J,8)=AJCLR(J,8)-TRP0CL(J)*POCEAN
1056     AJCLR(J,9)=AJCLR(J,9)-TRP1CL(J)*POCEAN
1057     AJCLR(J,3)=AJCLR(J,3)+(SNP1CL(J)*COSZ)*POCEAN
1058     AJCLR(J,7)=AJCLR(J,7)+TRHRCL(J)*POCEAN
1059     else if(ii.eq.2)then
1060     PLAND=PTYPE
1061     BJCLR(J,1)=BJCLR(J,1)+(S0*COSZ)*PLAND
1062     BJCLR(J,2)=BJCLR(J,2)+(SNP0CL(J)*COSZ)*PLAND
1063     BJCLR(J,4)=BJCLR(J,4)+(SRHRCL(J)*COSZ)*PLAND
1064     BJCLR(J,5)=BJCLR(J,5)+SRHE5*PLAND
1065     BJCLR(J,6)=BJCLR(J,6)+TRINCL(J)*PLAND
1066     BJCLR(J,8)=BJCLR(J,8)-TRP0CL(J)*PLAND
1067     BJCLR(J,9)=BJCLR(J,9)-TRP1CL(J)*PLAND
1068     BJCLR(J,3)=BJCLR(J,3)+(SNP1CL(J)*COSZ)*PLAND
1069     BJCLR(J,7)=BJCLR(J,7)+TRHRCL(J)*PLAND
1070     else
1071     POICE=PTYPE
1072     CJCLR(J,1)=CJCLR(J,1)+(S0*COSZ)*POICE
1073     CJCLR(J,2)=CJCLR(J,2)+(SNP0CL(J)*COSZ)*POICE
1074     CJCLR(J,4)=CJCLR(J,4)+(SRHRCL(J)*COSZ)*POICE
1075     CJCLR(J,5)=CJCLR(J,5)+SRHE5*POICE
1076     CJCLR(J,6)=CJCLR(J,6)+TRINCL(J)*POICE
1077     CJCLR(J,8)=CJCLR(J,8)-TRP0CL(J)*POICE
1078     CJCLR(J,9)=CJCLR(J,9)-TRP1CL(J)*POICE
1079     CJCLR(J,3)=CJCLR(J,3)+(SNP1CL(J)*COSZ)*POICE
1080     CJCLR(J,7)=CJCLR(J,7)+TRHRCL(J)*POICE
1081     endif
1082     C *********
1083     499 continue
1084     ALB(I,J,1)=SRNAN/(SRDAN+1.E-20)
1085     500 IM1=I 5304.
1086     C**** 5305.
1087     C**** END OF MAIN LOOP FOR I INDEX 5306.
1088     C**** 5307.
1089     600 CONTINUE 5345.
1090     C**** 5346.
1091     C**** END OF MAIN LOOP FOR J INDEX 5347.
1092     C**** 5348.
1093     C**** ACCUMULATE THE RADIATION DIAGNOSTICS 5394.
1094     C**** 5395.
1095     700 DO 780 J=1,JM 5396.
1096     DXYPJ=DXYP(J) 5397.
1097     IMAX=IM 5398.
1098     IF(J.EQ.1.OR.J.EQ.JM) IMAX=1 5399.
1099     DO 720 L=1,LM 5400.
1100     ASRHR=0. 5401.
1101     ATRHR=0. 5402.
1102     DO 710 I=1,IMAX 5403.
1103     ASRHR=ASRHR+SRHR(I,J,L+1)*COSZ2(I,J) 5404.
1104     710 ATRHR=ATRHR+TRHR(I,J,L+1) 5405.
1105     AJL(J,L,9)=AJL(J,L,9)+ASRHR 5406.
1106     720 AJL(J,L,10)=AJL(J,L,10)+ATRHR 5407.
1107     ASNFS1=0. 5408.
1108     BSNFS1=0. 5409.
1109     CSNFS1=0. 5410.
1110     ATNFS1=0. 5411.
1111     BTNFS1=0. 5412.
1112     CTNFS1=0. 5413.
1113     DO 770 I=1,IMAX 5414.
1114     SP=P(I,J) 5415.
1115     COSZ=COSZ2(I,J) 5416.
1116     PLAND=FDATA(I,J,2) 5417.
1117     PWATER=1.-PLAND
1118     POICE=ODATA(I,J,2)*(1.-PLAND) 5418.
1119     POCEAN=(1.-PLAND)-POICE 5419.
1120     if(POCEAN.LE.1.E-5)then
1121     POCEAN=0.
1122     POICE=PWATER
1123     endif
1124     JR=J
1125     DO 740 LR=1,3 5421.
1126     ASJL(J,LR,3)=ASJL(J,LR,3)+SRHRS(I,J,LR)*COSZ 5422.
1127     740 ASJL(J,LR,4)=ASJL(J,LR,4)+TRHRS(I,J,LR) 5423.
1128     DJ(JR,1)=DJ(JR,1)+(S0*COSZ)*DXYPJ 5440.
1129     DJ(JR,2)=DJ(JR,2)+(SNFS(I,J,4)*COSZ)*DXYPJ 5444.
1130     DJ(JR,3)=DJ(JR,3)+(SNFS(I,J,1)*COSZ)*DXYPJ 5448.
1131     DJ(JR,5)=DJ(JR,5)+(SRHR(I,J,1)*COSZ/(ALB(I,J,1)+1.E-20))*DXYPJ 5452.
1132     DJ(JR,6)=DJ(JR,6)+(SRHR(I,J,1)*COSZ)*DXYPJ 5456.
1133     DJ(JR,55)=DJ(JR,55)+BTMPW(I,J)*DXYPJ 5460.
1134     DJ(JR,67)=DJ(JR,67)+TRINCG(I,J)*DXYPJ 5464.
1135     DJ(JR,70)=DJ(JR,70)-TNFS(I,J,4)*DXYPJ 5468.
1136     DJ(JR,71)=DJ(JR,71)-TNFS(I,J,1)*DXYPJ 5472.
1137     770 CONTINUE 5485.
1138     780 CONTINUE 5492.
1139     C**** 5504.
1140     C**** UPDATE THE TEMPERATURES BY RADIATION 5505.
1141     C**** 5506.
1142     800 DO 820 J=1,JM 5507.
1143     IMAX=IM 5508.
1144     IF(J.EQ.1.OR.J.EQ.JM) IMAX=1 5509.
1145     DO 820 LR=1,3 5510.
1146     DO 820 I=1,IMAX 5511.
1147     820 RQT(I,J,LR)=RQT(I,J,LR)+(SRHRS(I,J,LR)*COSZ2(I,J) 5512.
1148     * +TRHRS(I,J,LR))*COE(LR+LM) 5513.
1149     840 DO 860 J=1,JM 5514.
1150     #if ( defined HR_DATA )
1151     dswhr(j)=DSWSRF(j)*COSZ1(1,j)
1152     dlwhr(j)=DLWSRF(j)
1153     #endif
1154     #if ( defined CLM )
1155     dsw4clm(j)=DSWSRF(j)*COSZ1(1,j)
1156     dlw4clm(j)=DLWSRF(j)
1157     swinr4clm(j)=DSWNIR(j)*COSZ1(1,j)
1158     swvis4clm(j)=DSWVIS(j)*COSZ1(1,j)
1159     c For TEM
1160     swtd4tem(j)=swtd4tem(j)+S0*COSZ1(1,j)
1161     swsd4tem(j)=swsd4tem(j)+DSWSRF(j)*COSZ1(1,j)
1162     nradd4tem(j)=nradd4tem(j)+1
1163     #endif
1164     IMAX=IM 5515.
1165     IF(J.EQ.1.OR.J.EQ.JM) IMAX=1 5516.
1166     if(HPRNT)then
1167     if(J.eq.1)then
1168     print *,' radia TAU=',TAU
1169     print *,' Before 860'
1170     print *,(T(1,1,L),L=1,LM)
1171     endif
1172     endif
1173     DO 860 L=1,LM 5517.
1174     DO 860 I=1,IMAX 5518.
1175     860 T(I,J,L)=T(I,J,L)+(SRHR(I,J,L+1)*COSZ1(I,J)+TRHR(I,J,L+1)) 5519.
1176     * *COE(L)/(P(I,J)*PK(I,J,L)) 5520.
1177     if(HPRNT)then
1178     print *,' radia TAU=',TAU
1179     print *,' after 860'
1180     print 'B66',T(1,1,1)
1181     print *,'COSZ1(1,1)',COSZ1(1,1),' COE(1)=',COE(1)
1182     print *,'SRHR'
1183     print *,(SRHR(1,1,L+1),L=1,LM)
1184     print *,'TRHR'
1185     print *,(TRHR(1,1,L+1),L=1,LM)
1186     print *,'T'
1187     print *,(T(1,1,L),L=1,LM)
1188     endif
1189     c if(ncallclm.ge.9)stop
1190     RETURN 5521.
1191     END 5522.

  ViewVC Help
Powered by ViewVC 1.1.22