/[MITgcm]/MITgcm/pkg/fizhi/fizhi_swrad.F
ViewVC logotype

Diff of /MITgcm/pkg/fizhi/fizhi_swrad.F

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

revision 1.1 by molod, Tue Jun 15 14:47:23 2004 UTC revision 1.17 by molod, Wed Aug 4 22:23:43 2004 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2  C $Name$  C $Name$
3                                                                                        
4        subroutine swrio (nymd,nhms,ndswr,myid,istrip,npcs,  #include "FIZHI_OPTIONS.h"
5       .                  pz,tz,qz,pkht,oz,co2,        subroutine swrio (nymd,nhms,bi,bj,ndswr,myid,istrip,npcs,
6       .                  albvisdr,albvisdf,albnirdr,albnirdf,       .        low_level,mid_level,im,jm,lm,
7       .                  dtradsw,dtswclr,radswg,swgclr,albedo,       .        pz,plz,plze,dpres,pkht,pkz,tz,qz,oz,co2,
8       .                  fdifpar,fdirpar,osr,osrclr,       .        albvisdr,albvisdf,albnirdr,albnirdf,
9       .                  im,jm,lm,sige,sig,dsig,ptop,       .        dtradsw,dtswclr,radswg,swgclr,
10       .                  nswcld,cldsw,cswmo,nswlz,swlz,       .        fdifpar,fdirpar,osr,osrclr,
11       .                  lpnt,qdiag,nd,       .        ptop,nswcld,cldsw,cswmo,nswlz,swlz,
12       .                  imstturb,qliqave,fccave,landtype,xlats,xlons)       .        lpnt,imstturb,qliqave,fccave,landtype,xlats,xlons)
13    
14        implicit none        implicit none
15        include 'diag.com'  #ifdef ALLOW_DIAGNOSTICS
16    #include "SIZE.h"
17    #include "diagnostics_SIZE.h"
18    #include "diagnostics.h"
19    #endif
20    
21  c Input Variables  c Input Variables
22  c ---------------  c ---------------
23        integer nymd,nhms,ndswr,istrip,npcs,nd        integer nymd,nhms,bi,bj,ndswr,myid,istrip,npcs
24          integer mid_level,low_level
25        integer im,jm,lm        ! Physics Grid        integer im,jm,lm        
26        real  ptop              ! Physics Grid        _RL  ptop
27        real  sige(lm+1)        ! Physics Grid        _RL pz(im,jm),plz(im,jm,lm),plze(im,jm,lm+1),dpres(im,jm,lm)
28        real   sig(lm)          ! Physics Grid        _RL pkht(im,jm,lm+1),pkz(im,jm,lm)
29        real  dsig(lm)          ! Physics Grid        _RL tz(im,jm,lm),qz(im,jm,lm)
30          _RL oz(im,jm,lm)
31        real    pz(im,jm)       ! Dynamics State        _RL co2
32        real    tz(im,jm,lm)    ! Dynamics State        _RL albvisdr(im,jm),albvisdf(im,jm),albnirdr(im,jm)
33        real  pkht(im,jm,lm)    ! Dynamics State        _RL albnirdf(im,jm)
34          _RL radswg(im,jm),swgclr(im,jm),fdifpar(im,jm),fdirpar(im,jm)
35        real    co2             ! Chemistry State        _RL osr(im,jm),osrclr(im,jm),dtradsw(im,jm,lm)
36        real    oz(im,jm,lm)    ! Chemistry Coupling        _RL dtswclr(im,jm,lm)
37        real    qz(im,jm,lm)    ! Chemistry Coupling + Dynamics State        integer nswcld,nswlz    
38          _RL cldsw(im,jm,lm),cswmo(im,jm,lm),swlz(im,jm,lm)  
39        real albvisdr(im,jm)    ! Land Coupling        logical lpnt            
40        real albvisdf(im,jm)    ! Land Coupling        integer imstturb        
41        real albnirdr(im,jm)    ! Land Coupling        _RL qliqave(im,jm,lm),fccave(im,jm,lm)  
42        real albnirdf(im,jm)    ! Land Coupling        integer landtype(im,jm)
43          _RL xlats(im,jm),xlons(im,jm)
       real   radswg(im,jm)    ! Shortwave Coupling  
       real   swgclr(im,jm)    ! Shortwave Coupling  
       real   albedo(im,jm)    ! Shortwave Coupling  
       real  fdifpar(im,jm)    ! Shortwave Coupling  
       real  fdirpar(im,jm)    ! Shortwave Coupling  
       real      osr(im,jm)    ! Shortwave Coupling  
       real   osrclr(im,jm)    ! Shortwave Coupling  
       real  dtradsw(im,jm,lm) ! Shortwave Tendency  
       real  dtswclr(im,jm,lm) ! Shortwave Tendency  
   
       integer nswcld,nswlz    ! Moist Coupling  
       real  cldsw(im,jm,lm)   ! Moist Coupling  
       real  cswmo(im,jm,lm)   ! Moist Coupling  
       real   swlz(im,jm,lm)   ! Moist Coupling  
   
       real  qdiag(im,jm,nd)   ! Diagnostics  
       logical lpnt            ! Point by Point Flag  
       integer imstturb        ! Turb Coupling  
       real qliqave(im,jm,lm)  ! Turb Coupling  
       real  fccave(im,jm,lm)  ! Turb Coupling  
   
       integer landtype(im,jm) ! Surface Land Type  
44    
45  c Local Variables  c Local Variables
46  c ---------------  c ---------------
47        integer   i,j,L,nn,nsecf,mid_level, low_level        integer   i,j,L,nn,nsecf
48        integer   nb2,ntmstp,nymd2,nhms2        integer   ntmstp,nymd2,nhms2
49        real      getcon,grav,cp,undef,pcheck        _RL      getcon,grav,cp,undef
50        real      ra,alf,reffw,reffi,tminv        _RL      ra,alf,reffw,reffi,tminv
51    
52        parameter ( reffw = 10.0 )   ! Effective radius for water droplets        parameter ( reffw = 10.0 )  
53        parameter ( reffi = 65.0 )   ! Effective radius for ice   particles        parameter ( reffi = 65.0 )  
54    
55        real      alat(im,jm)        _RL tdry(im,jm,lm)
56        real      alon(im,jm)        _RL TEMP1(im,jm)
57          _RL TEMP2(im,jm)
58        real          PKZ(im,jm,lm)        _RL zenith (im,jm)
59        real          PLZ(im,jm,lm)        _RL cldtot (im,jm,lm)
60        real         tdry(im,jm,lm)        _RL cldmxo (im,jm,lm)
61        real         PLZE(im,jm,lm+1)        _RL totcld (im,jm)
62        real        TEMP1(im,jm)        _RL cldlow (im,jm)
63        real        TEMP2(im,jm)        _RL cldmid (im,jm)
64        real      zenith (im,jm)        _RL cldhi  (im,jm)
65        real      cldtot (im,jm,lm)        _RL taulow (im,jm)
66        real      cldmxo (im,jm,lm)        _RL taumid (im,jm)
67        real      totcld (im,jm)        _RL tauhi  (im,jm)
68        real      cldlow (im,jm)        _RL tautype(im,jm,lm,3)
69        real      cldmid (im,jm)        _RL tau(im,jm,lm)
70        real      cldhi  (im,jm)        _RL albedo(im,jm)    
71        real      taulow (im,jm)  
72        real      taumid (im,jm)        _RL PK(ISTRIP,lm)
73        real      tauhi  (im,jm)        _RL qzl(ISTRIP,lm),CLRO(ISTRIP,lm)
74        real      tautype(im,jm,lm,3)        _RL TZL(ISTRIP,lm)
75        real      tau    (im,jm,lm)        _RL OZL(ISTRIP,lm)
76          _RL PLE(ISTRIP,lm+1)
77        real          PK(ISTRIP,lm)        _RL COSZ(ISTRIP)
78        real         qzl(ISTRIP,lm),  CLRO(ISTRIP,lm)        _RL dpstrip(ISTRIP,lm)
79        real         TZL(ISTRIP,lm)  
80        real         OZL(ISTRIP,lm)        _RL albuvdr(ISTRIP),albuvdf(ISTRIP)
81        real         PLE(ISTRIP,lm+1)        _RL albirdr(ISTRIP),albirdf(ISTRIP)
82        real        COSZ(ISTRIP)        _RL difpar (ISTRIP),dirpar (ISTRIP)
83    
84        real      albuvdr(ISTRIP),albuvdf(ISTRIP)        _RL fdirir(istrip),fdifir(istrip)
85        real      albirdr(ISTRIP),albirdf(ISTRIP)        _RL fdiruv(istrip),fdifuv(istrip)
86        real      difpar (ISTRIP),dirpar (ISTRIP)  
87          _RL flux(istrip,lm+1)
88        real      fdirir(istrip),fdifir(istrip)        _RL fluxclr(istrip,lm+1)
89        real      fdiruv(istrip),fdifuv(istrip)        _RL dtsw(istrip,lm)
90          _RL dtswc(istrip,lm)
91        real      flux   (istrip,lm+1)  
92        real      fluxclr(istrip,lm+1)        _RL taul   (istrip,lm)
93        real      dtsw   (istrip,lm)        _RL reff   (istrip,lm,2)
94        real      dtswc  (istrip,lm)        _RL tauc   (istrip,lm,2)
95          _RL taua   (istrip,lm)
96        real      taul   (istrip,lm)        _RL tstrip (istrip)
97        real      reff   (istrip,lm,2)  
98        real      tauc   (istrip,lm,2)        logical first
99        real      taua   (istrip,lm)        data first /.true./
       real      tstrip (istrip)  
   
       logical   first  
       data      first /.true./  
   
       integer   koz, kh2o  
       data      KOZ  /20/  
       data      kh2o /18/  
100    
101  C **********************************************************************  C **********************************************************************
102  C ****                       INITIALIZATION                         ****  C ****                       INITIALIZATION                         ****
# Line 135  C ************************************** Line 109  C **************************************
109        NTMSTP = nsecf(NDSWR)        NTMSTP = nsecf(NDSWR)
110        TMINV  = 1./float(ntmstp)        TMINV  = 1./float(ntmstp)
111    
       do j = 1,jm  
       do i = 1,im  
       PLZE(I,j,1) = SIGE(1)*PZ(I,j) + PTOP  
       enddo  
       enddo  
         
       DO L = 1,lm  
       do j = 1,jm  
       DO I = 1,im  
       PLZ (I,j,L  ) = SIG (L)  *PZ(I,j) + PTOP  
       PLZE(I,j,L+1) = SIGE(L+1)*PZ(I,j) + PTOP  
       ENDDO  
       ENDDO  
       ENDDO  
   
       call pkappa ( pz,pkht,pkz,ptop,sige,dsig,im,jm,lm )  
   
112  C Compute Temperature from Theta  C Compute Temperature from Theta
113  C ------------------------------  C ------------------------------
114        do L=1,lm        do L=1,lm
# Line 162  C ------------------------------ Line 119  C ------------------------------
119        enddo        enddo
120        enddo        enddo
121    
122  c Determine Level Indices for Low-Mid-High Cloud Regions        if (first .and. myid.eq.1 ) then
 c ------------------------------------------------------  
       low_level = lm  
       mid_level = lm  
       do L = lm-1,1,-1  
       pcheck = (1000.-ptop)*sig(l) + ptop  
       if (pcheck.gt.700.0) low_level = L  
       if (pcheck.gt.400.0) mid_level = L  
       enddo  
   
       if (first .and. myid.eq.0 ) then  
123        print *        print *
124        print *,'Low-Level Clouds are Grouped between levels: ',        print *,'Low-Level Clouds are Grouped between levels: ',
125       .         lm,' and ',low_level       .         lm,' and ',low_level
# Line 186  C ************************************** Line 133  C **************************************
133  C ****             CALCULATE COSINE OF THE ZENITH ANGLE             ****  C ****             CALCULATE COSINE OF THE ZENITH ANGLE             ****
134  C **********************************************************************  C **********************************************************************
135    
136        CALL ASTRO ( NYMD,   NHMS,  ALAT,ALON, im*jm, TEMP1,RA )        CALL ASTRO ( NYMD,   NHMS,  XLATS,XLONS, im*jm, TEMP1,RA )
137                     NYMD2 = NYMD                     NYMD2 = NYMD
138                     NHMS2 = NHMS                     NHMS2 = NHMS
139        CALL TICK  ( NYMD2,  NHMS2, NTMSTP )        CALL TICK  ( NYMD2,  NHMS2, NTMSTP )
140        CALL ASTRO ( NYMD2,  NHMS2, ALAT,ALON, im*jm, TEMP2,RA )        CALL ASTRO ( NYMD2,  NHMS2, XLATS,XLONS, im*jm, TEMP2,RA )
141    
142        do j = 1,jm        do j = 1,jm
143        do i = 1,im        do i = 1,im
# Line 202  C ************************************** Line 149  C **************************************
149        ENDDO        ENDDO
150        ENDDO        ENDDO
151    
   
152  C **********************************************************************  C **********************************************************************
153  c ****        Compute Two-Dimension Total Cloud Fraction (0-1)      ****  c ****        Compute Two-Dimension Total Cloud Fraction (0-1)      ****
154  C **********************************************************************  C **********************************************************************
# Line 223  c -------------------------------------- Line 169  c --------------------------------------
169          do L =1,lm          do L =1,lm
170          do j =1,jm          do j =1,jm
171          do i =1,im          do i =1,im
172           cldtot(i,j,L) =  min( 1.0,max(cldsw(i,j,L),fccave(i,j,L)/imstturb) )           cldtot(i,j,L)=min(1.0,max(cldsw(i,j,L),fccave(i,j,L)/imstturb))
173           cldmxo(i,j,L) =  min( 1.0,    cswmo(i,j,L) )           cldmxo(i,j,L)=min(1.0,cswmo(i,j,L))
174             swlz(i,j,L) =                swlz(i,j,L)+qliqave(i,j,L)/imstturb             swlz(i,j,L)=swlz(i,j,L)+qliqave(i,j,L)/imstturb
175          enddo          enddo
176          enddo          enddo
177          enddo          enddo
# Line 277  c ------------------------- Line 223  c -------------------------
223        if(icldfrc.gt.0) then        if(icldfrc.gt.0) then
224        do j=1,jm        do j=1,jm
225        do i=1,im        do i=1,im
226        qdiag(i,j,icldfrc) =  qdiag(i,j,icldfrc) + totcld(i,j)        qdiag(i,j,icldfrc,bi,bj) =  qdiag(i,j,icldfrc,bi,bj) + totcld(i,j)
227        enddo        enddo
228        enddo        enddo
229        ncldfrc = ncldfrc + 1        ncldfrc = ncldfrc + 1
# Line 287  c ------------------------- Line 233  c -------------------------
233        do L=1,lm        do L=1,lm
234        do j=1,jm        do j=1,jm
235        do i=1,im        do i=1,im
236        qdiag(i,j,icldras+L-1) = qdiag(i,j,icldras+L-1) + cswmo(i,j,L)        qdiag(i,j,icldras+L-1,bi,bj) = qdiag(i,j,icldras+L-1,bi,bj) +
237         .                                                     cswmo(i,j,L)
238        enddo        enddo
239        enddo        enddo
240        enddo        enddo
# Line 298  c ------------------------- Line 245  c -------------------------
245        do L=1,lm        do L=1,lm
246        do j=1,jm        do j=1,jm
247        do i=1,im        do i=1,im
248        qdiag(i,j,icldtot+L-1) = qdiag(i,j,icldtot+L-1) + cldtot(i,j,L)        qdiag(i,j,icldtot+L-1,bi,bj) = qdiag(i,j,icldtot+L-1,bi,bj) +
249         .                                                     cldtot(i,j,L)
250        enddo        enddo
251        enddo        enddo
252        enddo        enddo
# Line 308  c ------------------------- Line 256  c -------------------------
256        if( icldlow.gt.0 ) then        if( icldlow.gt.0 ) then
257        do j=1,jm        do j=1,jm
258        do i=1,im        do i=1,im
259        qdiag(i,j,icldlow) = qdiag(i,j,icldlow) + cldlow(i,j)        qdiag(i,j,icldlow,bi,bj) = qdiag(i,j,icldlow,bi,bj) + cldlow(i,j)
260        enddo        enddo
261        enddo        enddo
262        ncldlow = ncldlow + 1        ncldlow = ncldlow + 1
# Line 317  c ------------------------- Line 265  c -------------------------
265        if( icldmid.gt.0 ) then        if( icldmid.gt.0 ) then
266        do j=1,jm        do j=1,jm
267        do i=1,im        do i=1,im
268        qdiag(i,j,icldmid) = qdiag(i,j,icldmid) + cldmid(i,j)        qdiag(i,j,icldmid,bi,bj) = qdiag(i,j,icldmid,bi,bj) + cldmid(i,j)
269        enddo        enddo
270        enddo        enddo
271        ncldmid = ncldmid + 1        ncldmid = ncldmid + 1
# Line 326  c ------------------------- Line 274  c -------------------------
274        if( icldhi.gt.0 ) then        if( icldhi.gt.0 ) then
275        do j=1,jm        do j=1,jm
276        do i=1,im        do i=1,im
277        qdiag(i,j,icldhi) = qdiag(i,j,icldhi) + cldhi(i,j)        qdiag(i,j,icldhi,bi,bj) = qdiag(i,j,icldhi,bi,bj) + cldhi(i,j)
278        enddo        enddo
279        enddo        enddo
280        ncldhi = ncldhi + 1        ncldhi = ncldhi + 1
# Line 336  c ------------------------- Line 284  c -------------------------
284        do L=1,lm        do L=1,lm
285        do j=1,jm        do j=1,jm
286        do i=1,im        do i=1,im
287        qdiag(i,j,ilzrad+L-1) = qdiag(i,j,ilzrad+L-1) + swlz(i,j,L)*1.0e6        qdiag(i,j,ilzrad+L-1,bi,bj) = qdiag(i,j,ilzrad+L-1,bi,bj) +
288         .                                                     swlz(i,j,L)*1.0e6
289        enddo        enddo
290        enddo        enddo
291        enddo        enddo
# Line 348  c ------------------ Line 297  c ------------------
297        if( ialbvisdr.gt.0 ) then        if( ialbvisdr.gt.0 ) then
298        do j=1,jm        do j=1,jm
299        do i=1,im        do i=1,im
300        qdiag(i,j,ialbvisdr) = qdiag(i,j,ialbvisdr) + albvisdr(i,j)        qdiag(i,j,ialbvisdr,bi,bj) = qdiag(i,j,ialbvisdr,bi,bj) +
301         .                                                     albvisdr(i,j)
302        enddo        enddo
303        enddo        enddo
304        nalbvisdr = nalbvisdr + 1        nalbvisdr = nalbvisdr + 1
# Line 357  c ------------------ Line 307  c ------------------
307        if( ialbvisdf.gt.0 ) then        if( ialbvisdf.gt.0 ) then
308        do j=1,jm        do j=1,jm
309        do i=1,im        do i=1,im
310        qdiag(i,j,ialbvisdf) = qdiag(i,j,ialbvisdf) + albvisdf(i,j)        qdiag(i,j,ialbvisdf,bi,bj) = qdiag(i,j,ialbvisdf,bi,bj) +
311         .                                                     albvisdf(i,j)
312        enddo        enddo
313        enddo        enddo
314        nalbvisdf = nalbvisdf + 1        nalbvisdf = nalbvisdf + 1
# Line 366  c ------------------ Line 317  c ------------------
317        if( ialbnirdr.gt.0 ) then        if( ialbnirdr.gt.0 ) then
318        do j=1,jm        do j=1,jm
319        do i=1,im        do i=1,im
320        qdiag(i,j,ialbnirdr) = qdiag(i,j,ialbnirdr) + albnirdr(i,j)        qdiag(i,j,ialbnirdr,bi,bj) = qdiag(i,j,ialbnirdr,bi,bj) +
321         .                                                     albnirdr(i,j)
322        enddo        enddo
323        enddo        enddo
324        nalbnirdr = nalbnirdr + 1        nalbnirdr = nalbnirdr + 1
# Line 375  c ------------------ Line 327  c ------------------
327        if( ialbnirdf.gt.0 ) then        if( ialbnirdf.gt.0 ) then
328        do j=1,jm        do j=1,jm
329        do i=1,im        do i=1,im
330        qdiag(i,j,ialbnirdf) = qdiag(i,j,ialbnirdf) + albnirdf(i,j)        qdiag(i,j,ialbnirdf,bi,bj) = qdiag(i,j,ialbnirdf,bi,bj) +
331         .                                                     albnirdf(i,j)
332        enddo        enddo
333        enddo        enddo
334        nalbnirdf = nalbnirdf + 1        nalbnirdf = nalbnirdf + 1
# Line 383  c ------------------ Line 336  c ------------------
336    
337  C Compute Optical Thicknesses and Diagnostics  C Compute Optical Thicknesses and Diagnostics
338  C -------------------------------------------  C -------------------------------------------
339        call opthk ( tdry,plz,plze,swlz,cldtot,cldmxo,landtype,im,jm,lm,tautype )        call opthk(tdry,plz,plze,swlz,cldtot,cldmxo,landtype,im,jm,lm,
340         .                                                          tautype)
341    
342        do L = 1,lm        do L = 1,lm
343        do j = 1,jm        do j = 1,jm
344        do i = 1,im        do i = 1,im
345        tau(i,j,L) = tautype(i,j,L,1) + tautype(i,j,L,2) + tautype(i,j,L,3)        tau(i,j,L) = tautype(i,j,L,1)+tautype(i,j,L,2)+tautype(i,j,L,3)
346        enddo        enddo
347        enddo        enddo
348        enddo        enddo
# Line 397  C -------------------------------------- Line 351  C --------------------------------------
351        do L=1,lm        do L=1,lm
352        do j=1,jm        do j=1,jm
353        do i=1,im        do i=1,im
354        qdiag(i,j,itauave+L-1) = qdiag(i,j,itauave+L-1) + tau(i,j,L)*100/(plze(i,j,L+1)-plze(i,j,L))        qdiag(i,j,itauave+L-1,bi,bj) = qdiag(i,j,itauave+L-1,bi,bj) +
355         .                        tau(i,j,L)*100/(plze(i,j,L+1)-plze(i,j,L))
356        enddo        enddo
357        enddo        enddo
358        enddo        enddo
# Line 408  C -------------------------------------- Line 363  C --------------------------------------
363        do L=1,lm        do L=1,lm
364        do j=1,jm        do j=1,jm
365        do i=1,im        do i=1,im
366           if( cldtot(i,j,L).ne.0.0 ) then         if( cldtot(i,j,L).ne.0.0 ) then
367                qdiag(i,j,itaucld +L-1) = qdiag(i,j,itaucld +L-1) + tau(i,j,L)*100/(plze(i,j,L+1)-plze(i,j,L))          qdiag(i,j,itaucld +L-1,bi,bj) = qdiag(i,j,itaucld +L-1,bi,bj) +
368                qdiag(i,j,itaucldc+L-1) = qdiag(i,j,itaucldc+L-1) + 1.0       .                        tau(i,j,L)*100/(plze(i,j,L+1)-plze(i,j,L))
369           endif          qdiag(i,j,itaucldc+L-1,bi,bj) =
370         .                             qdiag(i,j,itaucldc+L-1,bi,bj) + 1.0
371           endif
372        enddo        enddo
373        enddo        enddo
374        enddo        enddo
# Line 420  C -------------------------------------- Line 377  C --------------------------------------
377  c Compute Low, Mid, and High Cloud Optical Depth Diagnostics  c Compute Low, Mid, and High Cloud Optical Depth Diagnostics
378  c ----------------------------------------------------------  c ----------------------------------------------------------
379        if( itaulow.ne.0 ) then        if( itaulow.ne.0 ) then
380            do j = 1,jm         do j = 1,jm
381            do i = 1,im         do i = 1,im
382            if( cldlow(i,j).ne.0.0 ) then          if( cldlow(i,j).ne.0.0 ) then
383                taulow(i,j) =  0.0           taulow(i,j) =  0.0
384                do L = low_level,lm           do L = low_level,lm
385                taulow(i,j) = taulow(i,j) + tau(i,j,L)            taulow(i,j) = taulow(i,j) + tau(i,j,L)
386                enddo           enddo
387                 qdiag(i,j,itaulow ) = qdiag(i,j,itaulow ) + taulow(i,j)           qdiag(i,j,itaulow,bi,bj ) = qdiag(i,j,itaulow,bi,bj ) +
388                 qdiag(i,j,itaulowc) = qdiag(i,j,itaulowc) + 1.0       .                                                    taulow(i,j)
389            endif           qdiag(i,j,itaulowc,bi,bj) = qdiag(i,j,itaulowc,bi,bj) + 1.0
390            enddo          endif
391            enddo         enddo
392           enddo
393        endif        endif
394    
395        if( itaumid.ne.0 ) then        if( itaumid.ne.0 ) then
396            do j = 1,jm         do j = 1,jm
397            do i = 1,im         do i = 1,im
398            if( cldmid(i,j).ne.0.0 ) then          if( cldmid(i,j).ne.0.0 ) then
399                taumid(i,j) =  0.0           taumid(i,j) =  0.0
400                do L = mid_level,low_level+1           do L = mid_level,low_level+1
401                taumid(i,j) = taumid(i,j) + tau(i,j,L)            taumid(i,j) = taumid(i,j) + tau(i,j,L)
402                enddo           enddo
403                 qdiag(i,j,itaumid ) = qdiag(i,j,itaumid ) + taumid(i,j)           qdiag(i,j,itaumid,bi,bj ) = qdiag(i,j,itaumid,bi,bj ) +
404                 qdiag(i,j,itaumidc) = qdiag(i,j,itaumidc) + 1.0       .                                                    taumid(i,j)
405            endif           qdiag(i,j,itaumidc,bi,bj) = qdiag(i,j,itaumidc,bi,bj) + 1.0
406            enddo          endif
407            enddo         enddo
408           enddo
409        endif        endif
410    
411        if( itauhi.ne.0 ) then        if( itauhi.ne.0 ) then
412            do j = 1,jm         do j = 1,jm
413            do i = 1,im         do i = 1,im
414            if( cldhi(i,j).ne.0.0 ) then          if( cldhi(i,j).ne.0.0 ) then
415                tauhi(i,j) =  0.0           tauhi(i,j) =  0.0
416                do L = 1,mid_level+1           do L = 1,mid_level+1
417                tauhi(i,j) = tauhi(i,j) + tau(i,j,L)            tauhi(i,j) = tauhi(i,j) + tau(i,j,L)
418                enddo           enddo
419                 qdiag(i,j,itauhi ) = qdiag(i,j,itauhi ) + tauhi(i,j)           qdiag(i,j,itauhi,bi,bj ) = qdiag(i,j,itauhi,bi,bj ) +
420                 qdiag(i,j,itauhic) = qdiag(i,j,itauhic) + 1.0       .                                                   tauhi(i,j)
421            endif           qdiag(i,j,itauhic,bi,bj) = qdiag(i,j,itauhic,bi,bj) + 1.0
422            enddo          endif
423            enddo         enddo
424           enddo
425        endif        endif
426    
427  C***********************************************************************  C***********************************************************************
# Line 476  C ************************************** Line 436  C **************************************
436    
437        CALL STRIP ( zenith,COSZ,im*jm,ISTRIP,1,NN )        CALL STRIP ( zenith,COSZ,im*jm,ISTRIP,1,NN )
438    
439        CALL STRIP ( plze, ple   ,im*jm,ISTRIP,lm+1,NN)        CALL STRIP ( plze,  ple   ,im*jm,ISTRIP,lm+1,NN)
440        CALL STRIP ( pkz , pk    ,im*jm,ISTRIP,lm  ,NN)        CALL STRIP ( pkz ,  pk    ,im*jm,ISTRIP,lm  ,NN)
441        CALL STRIP ( tdry, tzl   ,im*jm,ISTRIP,lm  ,NN)        CALL STRIP ( dpres,dpstrip,im*jm,ISTRIP,lm  ,NN)
442        CALL STRIP ( qz  , qzl   ,im*jm,ISTRIP,lm  ,NN)        CALL STRIP ( tdry,  tzl   ,im*jm,ISTRIP,lm  ,NN)
443        CALL STRIP ( oz  , ozl   ,im*jm,ISTRIP,lm  ,NN)        CALL STRIP ( qz  ,  qzl   ,im*jm,ISTRIP,lm  ,NN)
444        CALL STRIP ( tau , taul  ,im*jm,ISTRIP,lm  ,NN)        CALL STRIP ( oz  ,  ozl   ,im*jm,ISTRIP,lm  ,NN)
445          CALL STRIP ( tau ,  taul  ,im*jm,ISTRIP,lm  ,NN)
446    
447        CALL STRIP ( albvisdr,albuvdr,im*jm,ISTRIP,1,NN )        CALL STRIP ( albvisdr,albuvdr,im*jm,ISTRIP,1,NN )
448        CALL STRIP ( albvisdf,albuvdf,im*jm,ISTRIP,1,NN )        CALL STRIP ( albvisdf,albuvdf,im*jm,ISTRIP,1,NN )
# Line 528  C ****     Compute Mass-Weighted Theta T Line 489  C ****     Compute Mass-Weighted Theta T
489  C **********************************************************************  C **********************************************************************
490    
491        do l=1,lm        do l=1,lm
       alf = grav/(cp*dsig(L)*100)  
492        do i=1,istrip        do i=1,istrip
493          alf = grav*(ple(i,L+1)-ptop)/(cp*dpstrip(i,L)*100)
494        dtsw (i,L) = alf*( flux   (i,L)-flux   (i,L+1) )/pk(i,L)        dtsw (i,L) = alf*( flux   (i,L)-flux   (i,L+1) )/pk(i,L)
495        dtswc(i,L) = alf*( fluxclr(i,L)-fluxclr(i,L+1) )/pk(i,L)        dtswc(i,L) = alf*( fluxclr(i,L)-fluxclr(i,L+1) )/pk(i,L)
496        enddo        enddo
# Line 550  C ************************************** Line 511  C **************************************
511  c Calculate Mean Albedo  c Calculate Mean Albedo
512  c ---------------------  c ---------------------
513        do i=1,istrip        do i=1,istrip
514        if( cosz(i).gt.0.0 ) then         if( cosz(i).gt.0.0 ) then
515              tstrip(i) = 1.0 - flux(i,lm+1)/( fdirir(i)+fdifir(i)+dirpar(i)+difpar(i)          tstrip(i) = 1.0 - flux(i,lm+1)/
516       .                                     + fdiruv(i)+fdifuv(i) )       . ( fdirir(i)+fdifir(i)+dirpar(i)+difpar(i) + fdiruv(i)+fdifuv(i) )
517          if( tstrip(i).lt.0.0 ) tstrip(i) = undef          if( tstrip(i).lt.0.0 ) tstrip(i) = undef
518        else         else
519              tstrip(i) = undef          tstrip(i) = undef
520        endif         endif
521        enddo        enddo
522        call paste ( tstrip,albedo,istrip,im*jm,1,nn )        call paste ( tstrip,albedo,istrip,im*jm,1,nn )
523    
# Line 568  c ---------------------- Line 529  c ----------------------
529        do j=1,jm        do j=1,jm
530        do i=1,im        do i=1,im
531        if( albedo(i,j).ne.undef ) then        if( albedo(i,j).ne.undef ) then
532        qdiag(i,j,ialbedo ) = qdiag(i,j,ialbedo ) + albedo(i,j)        qdiag(i,j,ialbedo,bi,bj ) = qdiag(i,j,ialbedo,bi,bj )+albedo(i,j)
533        qdiag(i,j,ialbedoc) = qdiag(i,j,ialbedoc) + 1.0        qdiag(i,j,ialbedoc,bi,bj) = qdiag(i,j,ialbedoc,bi,bj) + 1.0
534        endif        endif
535        enddo        enddo
536        enddo        enddo
# Line 623  C                  tau(im,jm,lm,2):  Sus Line 584  C                  tau(im,jm,lm,2):  Sus
584  C                  tau(im,jm,lm,3):  Raindrops  C                  tau(im,jm,lm,3):  Raindrops
585  C  C
586  C***********************************************************************  C***********************************************************************
 C*                  GODDARD LABORATORY FOR ATMOSPHERES                 *  
 C***********************************************************************  
587    
588        implicit none        implicit none
589    
590        integer  im,jm,lm,i,j,L        integer  im,jm,lm,i,j,L
591    
592        real  tl(im,jm,lm)        _RL  tl(im,jm,lm)
593        real  pl(im,jm,lm)        _RL  pl(im,jm,lm)
594        real ple(im,jm,lm+1)        _RL ple(im,jm,lm+1)
595        real  lz(im,jm,lm)        _RL  lz(im,jm,lm)
596        real  cf(im,jm,lm)        _RL  cf(im,jm,lm)
597        real cfm(im,jm,lm)        _RL cfm(im,jm,lm)
598        real tau(im,jm,lm,3)        _RL tau(im,jm,lm,3)
599        integer lwi(im,jm)        integer lwi(im,jm)
600    
601        real dp, alf, fracls, fraccu        _RL dp, alf, fracls, fraccu
602        real tauice, tauh2o, tauras        _RL tauice, tauh2o, tauras
603    
604  c Compute Cloud Optical Depths  c Compute Cloud Optical Depths
605  c ----------------------------  c ----------------------------
# Line 672  c --------------- Line 631  c ---------------
631    
632  c Large-Scale Water  c Large-Scale Water
633  c -----------------  c -----------------
634    C Over Land
635           if( lwi(i,j).le.10 ) then           if( lwi(i,j).le.10 ) then
636                  tauh2o = max( 0.0020, 0.200*min(200*lz(i,j,L)*1000,1.0) )  ! Over Land            tauh2o = max( 0.0020, 0.200*min(200*lz(i,j,L)*1000,1.0) )  
637                  tau(i,j,L,3) = fracls*alf*tauh2o*dp            tau(i,j,L,3) = fracls*alf*tauh2o*dp
638           else           else
639              if( lz(i,j,L).eq.0.0 ) then  C Non-Precipitation Clouds Over Ocean
640                  tauh2o = .12                      ! Non-Precipitation Clouds Over Ocean            if( lz(i,j,L).eq.0.0 ) then
641                  tau(i,j,L,2) = fracls*alf*tauh2o*dp             tauh2o = .12                      
642              else             tau(i,j,L,2) = fracls*alf*tauh2o*dp
643                  tauh2o = max( 0.0020, 0.120*min( 20*lz(i,j,L)*1000,1.0) )  ! Over Ocean            else
644                  tau(i,j,L,3) = fracls*alf*tauh2o*dp  C Over Ocean
645              endif             tauh2o = max( 0.0020, 0.120*min( 20*lz(i,j,L)*1000,1.0) )  
646               tau(i,j,L,3) = fracls*alf*tauh2o*dp
647              endif
648           endif           endif
649    
650  c Sub-Grid Convective  c Sub-Grid Convective
# Line 805  c*************************************** Line 767  c***************************************
767    
768  c-----Explicit Inline Directives  c-----Explicit Inline Directives
769    
770  #if CRAY  #ifdef CRAY
771  #if f77  #ifdef f77
772  cfpp$ expand (expmn)  cfpp$ expand (expmn)
773  #endif  #endif
 #if f90  
 !DIR$ inline always expmn  
774  #endif  #endif
775  #endif        _RL expmn
       real expmn  
776    
777  c-----input parameters  c-----input parameters
778    
779        integer m,n,ndim,np,ict,icb        integer m,n,ndim,np,ict,icb
780        real pl(m,ndim,np+1),ta(m,ndim,np),wa(m,ndim,np),oa(m,ndim,np)        _RL pl(m,ndim,np+1),ta(m,ndim,np),wa(m,ndim,np),oa(m,ndim,np)
781        real  taucld(m,ndim,np,2),reff(m,ndim,np,2)        _RL  taucld(m,ndim,np,2),reff(m,ndim,np,2)
782        real  fcld(m,ndim,np),taual(m,ndim,np)        _RL  fcld(m,ndim,np),taual(m,ndim,np)
783        real  rsirbm(m,ndim),rsirdf(m,ndim),        _RL  rsirbm(m,ndim),rsirdf(m,ndim),
784       *     rsuvbm(m,ndim),rsuvdf(m,ndim),cosz(m,ndim),co2       *     rsuvbm(m,ndim),rsuvdf(m,ndim),cosz(m,ndim),co2
785    
786  c-----output parameters  c-----output parameters
787    
788        real  flx(m,ndim,np+1),flc(m,ndim,np+1)        _RL  flx(m,ndim,np+1),flc(m,ndim,np+1)
789        real  fdirir(m,ndim),fdifir(m,ndim)        _RL  fdirir(m,ndim),fdifir(m,ndim)
790        real  fdirpar(m,ndim),fdifpar(m,ndim)        _RL  fdirpar(m,ndim),fdifpar(m,ndim)
791        real  fdiruv(m,ndim),fdifuv(m,ndim)        _RL  fdiruv(m,ndim),fdifuv(m,ndim)
792    
793  c-----temporary array  c-----temporary array
794    
795        integer i,j,k,ik        integer i,j,k
796        real  cc(m,n,3),tauclb(m,n,np),tauclf(m,n,np)        _RL  cc(m,n,3),tauclb(m,n,np),tauclf(m,n,np)
797        real  dp(m,n,np),wh(m,n,np),oh(m,n,np),scal(m,n,np)        _RL  dp(m,n,np),wh(m,n,np),oh(m,n,np),scal(m,n,np)
798        real  swh(m,n,np+1),so2(m,n,np+1),df(m,n,np+1)        _RL  swh(m,n,np+1),so2(m,n,np+1),df(m,n,np+1)
799        real  sdf(m,n),sclr(m,n),csm(m,n),taux,x        _RL  sdf(m,n),sclr(m,n),csm(m,n),x
800    
801  c-----------------------------------------------------------------  c-----------------------------------------------------------------
802    
# Line 1050  c*************************************** Line 1009  c***************************************
1009  c-----input parameters  c-----input parameters
1010    
1011        integer m,n,ndim,np,ict,icb        integer m,n,ndim,np,ict,icb
1012        real  cosz(m,ndim),fcld(m,ndim,np),taucld(m,ndim,np,2)        _RL  cosz(m,ndim),fcld(m,ndim,np),taucld(m,ndim,np,2)
1013    
1014  c-----output parameters  c-----output parameters
1015    
1016        real  cc(m,n,3),tauclb(m,n,np),tauclf(m,n,np)        _RL  cc(m,n,3),tauclb(m,n,np),tauclf(m,n,np)
1017    
1018  c-----temporary variables  c-----temporary variables
1019    
1020        integer i,j,k,im,it,ia,kk        integer i,j,k,im,it,ia,kk
1021        real   fm,ft,fa,xai,taucl,taux        _RL   fm,ft,fa,xai,taux
1022    
1023  c-----pre-computed table  c-----pre-computed table
1024    
1025        integer   nm,nt,na        integer   nm,nt,na
1026        parameter (nm=11,nt=9,na=11)        parameter (nm=11,nt=9,na=11)
1027        real   dm,dt,da,t1,caib(nm,nt,na),caif(nt,na)        _RL   dm,dt,da,t1,caib(nm,nt,na),caif(nt,na)
1028        parameter (dm=0.1,dt=0.30103,da=0.1,t1=-0.9031)        parameter (dm=0.1,dt=0.30103,da=0.1,t1=-0.9031)
1029    
1030  c-----include the pre-computed table for cai  c-----include the pre-computed table for cai
1031    
1032        include 'cai.dat'  #include "cai-dat.h"
1033        save caib,caif        save caib,caif
1034    
1035    
# Line 1268  c*************************************** Line 1227  c***************************************
1227    
1228  c-----Explicit Inline Directives  c-----Explicit Inline Directives
1229    
1230  #if CRAY  #ifdef CRAY
1231  #if f77  #ifdef f77
1232  cfpp$ expand (deledd)  cfpp$ expand (deledd)
1233  cfpp$ expand (sagpol)  cfpp$ expand (sagpol)
1234  cfpp$ expand (expmn)  cfpp$ expand (expmn)
1235  #endif  #endif
 #if f90  
 !DIR$ inline always deledd  
 !DIR$ inline always sagpol  
 !DIR$ inline always expmn  
 #endif  
1236  #endif  #endif
1237        real expmn        _RL expmn
1238    
1239  c-----input parameters  c-----input parameters
1240    
1241        integer m,n,ndim,np,ict,icb        integer m,n,ndim,np,ict,icb
1242        real  taucld(m,ndim,np,2),reff(m,ndim,np,2),fcld(m,ndim,np)        _RL  taucld(m,ndim,np,2),reff(m,ndim,np,2),fcld(m,ndim,np)
1243        real  tauclb(m,n,np),tauclf(m,n,np),cc(m,n,3)        _RL  tauclb(m,n,np),tauclf(m,n,np),cc(m,n,3)
1244        real  rsirbm(m,ndim),rsirdf(m,ndim)        _RL  rsirbm(m,ndim),rsirdf(m,ndim)
1245        real  wh(m,n,np),taual(m,ndim,np),csm(m,n)        _RL  wh(m,n,np),taual(m,ndim,np),csm(m,n)
1246    
1247  c-----output (updated) parameters  c-----output (updated) parameters
1248    
1249        real  flx(m,ndim,np+1),flc(m,ndim,np+1)        _RL  flx(m,ndim,np+1),flc(m,ndim,np+1)
1250        real  fdirir(m,ndim),fdifir(m,ndim)        _RL  fdirir(m,ndim),fdifir(m,ndim)
1251    
1252  c-----static parameters  c-----static parameters
1253    
1254        integer nk,nband        integer nk,nband
1255        parameter (nk=10,nband=3)        parameter (nk=10,nband=3)
1256        real  xk(nk),hk(nband,nk),ssaal(nband),asyal(nband)        _RL  xk(nk),hk(nband,nk),ssaal(nband),asyal(nband)
1257        real  aia(nband,3),awa(nband,3),aig(nband,3),awg(nband,3)        _RL  aia(nband,3),awa(nband,3),aig(nband,3),awg(nband,3)
1258    
1259  c-----temporary array  c-----temporary array
1260    
1261        integer ib,ik,i,j,k        integer ib,ik,i,j,k
1262        real  ssacl(m,n,np),asycl(m,n,np)        _RL  ssacl(m,n,np),asycl(m,n,np)
1263        real  rr(m,n,np+1,2),tt(m,n,np+1,2),td(m,n,np+1,2),        _RL  rr(m,n,np+1,2),tt(m,n,np+1,2),td(m,n,np+1,2),
1264       *       rs(m,n,np+1,2),ts(m,n,np+1,2)       *       rs(m,n,np+1,2),ts(m,n,np+1,2)
1265        real  rssab(m,n,np+1),rabx(m,n,np+1),rsabx(m,n,np+1)        _RL  fall(m,n,np+1),fclr(m,n,np+1)
1266        real  fall(m,n,np+1),fclr(m,n,np+1)        _RL  fsdir(m,n),fsdif(m,n)
1267        real  fsdir(m,n),fsdif(m,n)  
1268          _RL  tauwv,tausto,ssatau,asysto,tauto,ssato,asyto
1269        real  tauwv,tausto,ssatau,asysto,tauto,ssato,asyto        _RL  taux,reff1,reff2,w1,w2,g1,g2
1270        real  taux,reff1,reff2,w1,w2,g1,g2        _RL  ssaclt(m,n),asyclt(m,n)
1271        real  ssaclt(m,n),asyclt(m,n)        _RL  rr1t(m,n),tt1t(m,n),td1t(m,n),rs1t(m,n),ts1t(m,n)
1272        real  rr1t(m,n),tt1t(m,n),td1t(m,n),rs1t(m,n),ts1t(m,n)        _RL  rr2t(m,n),tt2t(m,n),td2t(m,n),rs2t(m,n),ts2t(m,n)
       real  rr2t(m,n),tt2t(m,n),td2t(m,n),rs2t(m,n),ts2t(m,n)  
1273    
1274  c-----water vapor absorption coefficient for 10 k-intervals.  c-----water vapor absorption coefficient for 10 k-intervals.
1275  c     unit: cm^2/gm  c     unit: cm^2/gm
# Line 1449  c-----integration over the k-distributio Line 1402  c-----integration over the k-distributio
1402              do i= 1, m              do i= 1, m
1403    
1404               tauwv=xk(ik)*wh(i,j,k)               tauwv=xk(ik)*wh(i,j,k)
1405    
1406  c-----compute total optical thickness, single scattering albedo,  c-----compute total optical thickness, single scattering albedo,
1407  c     and asymmetry factor.  c     and asymmetry factor.
1408    
# Line 1575  c     in certain parallel processors. Line 1528  c     in certain parallel processors.
1528    
1529  c-----flux calculations  c-----flux calculations
1530    
1531           do k= 1, np+1
1532            do j= 1, n
1533             do i= 1, m
1534              fclr(i,j,k) = 0.
1535              fall(i,j,k) = 0.
1536             enddo
1537            enddo
1538           enddo
1539           do j= 1, n
1540            do i= 1, m
1541             fsdir(i,j) = 0.
1542             fsdif(i,j) = 0.
1543            enddo
1544           enddo
1545    
1546          call cldflx (m,n,np,ict,icb,cc,rr,tt,td,rs,ts,          call cldflx (m,n,np,ict,icb,cc,rr,tt,td,rs,ts,
1547       *               fclr,fall,fsdir,fsdif)       *               fclr,fall,fsdir,fsdif)
1548    
# Line 1683  c*************************************** Line 1651  c***************************************
1651    
1652  c-----Explicit Inline Directives    c-----Explicit Inline Directives  
1653        
1654  #if CRAY  #ifdef CRAY
1655  #if f77    #ifdef f77  
1656  cfpp$ expand (deledd)  cfpp$ expand (deledd)
1657  cfpp$ expand (sagpol)  cfpp$ expand (sagpol)
1658  #endif    #endif  
 #if f90    
 !DIR$ inline always deledd  
 !DIR$ inline always sagpol  
 #endif  
1659  #endif  #endif
1660    
1661  c-----input parameters  c-----input parameters
1662    
1663        integer m,n,ndim,np,ict,icb        integer m,n,ndim,np,ict,icb
1664        real  taucld(m,ndim,np,2),reff(m,ndim,np,2),fcld(m,ndim,np)        _RL  taucld(m,ndim,np,2),reff(m,ndim,np,2),fcld(m,ndim,np)
1665        real  tauclb(m,n,np),tauclf(m,n,np),cc(m,n,3)        _RL  tauclb(m,n,np),tauclf(m,n,np),cc(m,n,3)
1666        real  oh(m,n,np),dp(m,n,np),taual(m,ndim,np)        _RL  oh(m,n,np),dp(m,n,np),taual(m,ndim,np)
1667        real  rsuvbm(m,ndim),rsuvdf(m,ndim),csm(m,n)        _RL  rsuvbm(m,ndim),rsuvdf(m,ndim),csm(m,n)
1668    
1669  c-----output (updated) parameter  c-----output (updated) parameter
1670    
1671        real  flx(m,ndim,np+1),flc(m,ndim,np+1)        _RL  flx(m,ndim,np+1),flc(m,ndim,np+1)
1672        real  fdirpar(m,ndim),fdifpar(m,ndim)        _RL  fdirpar(m,ndim),fdifpar(m,ndim)
1673        real  fdiruv(m,ndim),fdifuv(m,ndim)        _RL  fdiruv(m,ndim),fdifuv(m,ndim)
1674    
1675  c-----static parameters  c-----static parameters
1676    
1677        integer nband        integer nband
1678        parameter (nband=8)        parameter (nband=8)
1679        real  hk(nband),xk(nband),ry(nband)        _RL  hk(nband),xk(nband),ry(nband)
1680        real  asyal(nband),ssaal(nband),aig(3),awg(3)        _RL  asyal(nband),ssaal(nband),aig(3),awg(3)
1681    
1682  c-----temporary array  c-----temporary array
1683    
1684        integer i,j,k,ib        integer i,j,k,ib
1685        real  taurs,tauoz,tausto,ssatau,asysto,tauto,ssato,asyto        _RL  taurs,tauoz,tausto,ssatau,asysto,tauto,ssato,asyto
1686        real  taux,reff1,reff2,g1,g2,asycl(m,n,np)        _RL  taux,reff1,reff2,g1,g2,asycl(m,n,np)
1687        real  td(m,n,np+1,2),rr(m,n,np+1,2),tt(m,n,np+1,2),        _RL  td(m,n,np+1,2),rr(m,n,np+1,2),tt(m,n,np+1,2),
1688       *       rs(m,n,np+1,2),ts(m,n,np+1,2)       *       rs(m,n,np+1,2),ts(m,n,np+1,2)
1689        real  upflux(m,n,np+1),dwflux(m,n,np+1),        _RL  fall(m,n,np+1),fclr(m,n,np+1),fsdir(m,n),fsdif(m,n)
1690       *     rssab(m,n,np+1),rabx(m,n,np+1),rsabx(m,n,np+1)        _RL  asyclt(m,n)
1691        real  fall(m,n,np+1),fclr(m,n,np+1),fsdir(m,n),fsdif(m,n)        _RL  rr1t(m,n),tt1t(m,n),td1t(m,n),rs1t(m,n),ts1t(m,n)
1692        real  asyclt(m,n)        _RL  rr2t(m,n),tt2t(m,n),td2t(m,n),rs2t(m,n),ts2t(m,n)
       real  rr1t(m,n),tt1t(m,n),td1t(m,n),rs1t(m,n),ts1t(m,n)  
       real  rr2t(m,n),tt2t(m,n),td2t(m,n),rs2t(m,n),ts2t(m,n)  
1693    
1694  c-----hk is the fractional extra-terrestrial solar flux.  c-----hk is the fractional extra-terrestrial solar flux.
1695  c     the sum of hk is 0.47074.  c     the sum of hk is 0.47074.
# Line 1938  c-----compute reflectance and transmitta Line 1900  c-----compute reflectance and transmitta
1900    
1901  c-----flux calculations  c-----flux calculations
1902    
1903           do k= 1, np+1
1904            do j= 1, n
1905             do i= 1, m
1906              fclr(i,j,k) = 0.
1907              fall(i,j,k) = 0.
1908             enddo
1909            enddo
1910           enddo
1911           do j= 1, n
1912            do i= 1, m
1913             fsdir(i,j) = 0.
1914             fsdif(i,j) = 0.
1915            enddo
1916           enddo
1917          call cldflx (m,n,np,ict,icb,cc,rr,tt,td,rs,ts,          call cldflx (m,n,np,ict,icb,cc,rr,tt,td,rs,ts,
1918       *               fclr,fall,fsdir,fsdif)       *               fclr,fall,fsdir,fsdif)
1919    
# Line 2004  c*************************************** Line 1980  c***************************************
1980    
1981  c-----Explicit Inline Directives    c-----Explicit Inline Directives  
1982        
1983  #if CRAY  #ifdef CRAY
1984  #if f77    #ifdef f77  
1985  cfpp$ expand (expmn)  cfpp$ expand (expmn)
1986  #endif    #endif  
 #if f90    
 !DIR$ inline always expmn  
1987  #endif  #endif
1988  #endif        _RL expmn
       real expmn  
1989    
1990        real  zero,one,two,three,four,fourth,seven,tumin        _RL  zero,one,two,three,four,fourth,seven,tumin
1991        parameter (one=1., three=3.)        parameter (one=1., three=3.)
1992        parameter (seven=7., two=2.)        parameter (seven=7., two=2.)
1993        parameter (four=4., fourth=.25)        parameter (four=4., fourth=.25)
1994        parameter (zero=0., tumin=1.e-20)        parameter (zero=0., tumin=1.e-20)
1995    
1996  c-----input parameters  c-----input parameters
1997        real  tau,ssc,g0,csm        _RL  tau,ssc,g0,csm
1998    
1999  c-----output parameters  c-----output parameters
2000        real  rr,tt,td        _RL  rr,tt,td
2001    
2002  c-----temporary parameters  c-----temporary parameters
2003    
2004        real  zth,ff,xx,taup,sscp,gp,gm1,gm2,gm3,akk,alf1,alf2,        _RL  zth,ff,xx,taup,sscp,gp,gm1,gm2,gm3,akk,alf1,alf2,
2005       *     all,bll,st7,st8,cll,dll,fll,ell,st1,st2,st3,st4       *     all,bll,st7,st8,cll,dll,fll,ell,st1,st2,st3,st4
2006  c  c
2007                  zth = one / csm                  zth = one / csm
# Line 2131  c*************************************** Line 2104  c***************************************
2104    
2105  c-----Explicit Inline Directives    c-----Explicit Inline Directives  
2106        
2107  #if CRAY  #ifdef CRAY
2108  #if f77    #ifdef f77  
2109  cfpp$ expand (expmn)  cfpp$ expand (expmn)
2110  #endif    #endif  
 #if f90    
 !DIR$ inline always expmn  
 #endif  
2111  #endif  #endif
2112        real expmn        _RL expmn
2113    
2114        real  one,three,four        _RL  one,three,four
2115        parameter (one=1., three=3., four=4.)        parameter (one=1., three=3., four=4.)
2116    
2117  c-----output parameters:  c-----output parameters:
2118    
2119        real  tau,ssc,g0        _RL  tau,ssc,g0
2120    
2121  c-----output parameters:  c-----output parameters:
2122    
2123        real  rll,tll        _RL  rll,tll
2124    
2125  c-----temporary arrays  c-----temporary arrays
2126    
2127        real  xx,uuu,ttt,emt,up1,um1,st1        _RL  xx,uuu,ttt,emt,up1,um1,st1
2128    
2129               xx  = one-ssc*g0               xx  = one-ssc*g0
2130               uuu = sqrt( xx/(one-ssc))               uuu = sqrt( xx/(one-ssc))
# Line 2176  c*************************************** Line 2146  c***************************************
2146    
2147  c*******************************************************************  c*******************************************************************
2148  c compute exponential for arguments in the range 0> fin > -10.  c compute exponential for arguments in the range 0> fin > -10.
2149    c*******************************************************************
2150          implicit none
2151          _RL  fin,expmn
2152    
2153          _RL one,expmin,e1,e2,e3,e4
2154        parameter (one=1.0, expmin=-10.0)        parameter (one=1.0, expmin=-10.0)
2155        parameter (e1=1.0,        e2=-2.507213e-1)        parameter (e1=1.0,        e2=-2.507213e-1)
2156        parameter (e3=2.92732e-2, e4=-3.827800e-3)        parameter (e3=2.92732e-2, e4=-3.827800e-3)
       real  fin,expmn  
2157    
2158        if (fin .lt. expmin) fin = expmin        if (fin .lt. expmin) fin = expmin
2159        expmn = ((e4*fin + e3)*fin+e2)*fin+e1        expmn = ((e4*fin + e3)*fin+e2)*fin+e1
# Line 2230  c-----input parameters Line 2203  c-----input parameters
2203    
2204        integer m,n,np,ict,icb        integer m,n,np,ict,icb
2205    
2206        real  rr(m,n,np+1,2),tt(m,n,np+1,2),td(m,n,np+1,2)        _RL  rr(m,n,np+1,2),tt(m,n,np+1,2),td(m,n,np+1,2)
2207        real  rs(m,n,np+1,2),ts(m,n,np+1,2)        _RL  rs(m,n,np+1,2),ts(m,n,np+1,2)
2208        real  cc(m,n,3)        _RL  cc(m,n,3)
2209    
2210  c-----temporary array  c-----temporary array
2211    
2212        integer i,j,k,ih,im,is        integer i,j,k,ih,im,is
2213        real  rra(m,n,np+1,2,2),tta(m,n,np+1,2,2),tda(m,n,np+1,2,2)        _RL  rra(m,n,np+1,2,2),tta(m,n,np+1,2,2),tda(m,n,np+1,2,2)
2214        real  rsa(m,n,np+1,2,2),rxa(m,n,np+1,2,2)        _RL  rsa(m,n,np+1,2,2),rxa(m,n,np+1,2,2)
2215        real  ch(m,n),cm(m,n),ct(m,n),flxdn(m,n,np+1)        _RL  ch(m,n),cm(m,n),ct(m,n),flxdn(m,n,np+1)
2216        real  fdndir(m,n),fdndif(m,n),fupdif        _RL  fdndir(m,n),fdndif(m,n),fupdif
2217        real  denm,xx        _RL  denm,xx
2218    
2219  c-----output parameters  c-----output parameters
2220    
2221        real  fclr(m,n,np+1),fall(m,n,np+1)        _RL  fclr(m,n,np+1),fall(m,n,np+1)
2222        real  fsdir(m,n),fsdif(m,n)        _RL  fsdir(m,n),fsdif(m,n)
2223    
2224  c-----initialize all-sky flux (fall) and surface downward fluxes  c-----initialize all-sky flux (fall) and surface downward fluxes
2225    
# Line 2546  c     due to co2 absorption. Line 2519  c     due to co2 absorption.
2519  c-----input parameters  c-----input parameters
2520    
2521        integer m,n,np        integer m,n,np
2522        real  csm(m,n),swc(m,n,np+1),swh(m,n,np+1),cah(22,19)        _RL  csm(m,n),swc(m,n,np+1),swh(m,n,np+1),cah(22,19)
2523    
2524  c-----output (undated) parameter  c-----output (undated) parameter
2525    
2526        real  df(m,n,np+1)        _RL  df(m,n,np+1)
2527    
2528  c-----temporary array  c-----temporary array
2529    
2530        integer i,j,k,ic,iw        integer i,j,k,ic,iw
2531        real  xx,clog,wlog,dc,dw,x1,x2,y2        _RL  xx,clog,wlog,dc,dw,x1,x2,y2
2532    
2533  c********************************************************************  c********************************************************************
2534  c-----include co2 look-up table  c-----include co2 look-up table
2535    
2536        include 'cah.dat'  #include "cah-dat.h"
2537        save cah  c     save cah
2538    
2539  c********************************************************************  c********************************************************************
2540  c-----table look-up for the reduction of clear-sky solar  c-----table look-up for the reduction of clear-sky solar
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.17
  ViewVC Help
Powered by ViewVC 1.1.22