/[MITgcm]/MITgcm/pkg/aim/phy_radiat.F
ViewVC logotype

Annotation of /MITgcm/pkg/aim/phy_radiat.F

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


Revision 1.5 - (hide annotations) (download)
Thu Sep 6 13:19:54 2001 UTC (22 years, 8 months ago) by adcroft
Branch: MAIN
CVS Tags: icebear2, checkpoint44h_pre, release1_p12, release1_p10, release1_p11, release1_p16, release1_p15, ecco_c44_e17, ecco_c44_e16, ecco_c44_e21, ecco_c44_e26, ecco_c44_e27, ecco_c44_e24, icebear5, icebear4, checkpoint44f_pre, icebear3, checkpoint46f_post, release1_p13_pre, checkpoint46d_pre, checkpoint46e_post, release1-branch_tutorials, release1_p14, checkpoint44g_post, checkpoint46h_pre, checkpoint44h_post, release1_p12_pre, checkpoint44e_post, checkpoint46e_pre, ecco-branch-mod4, checkpoint43a-release1mods, checkpoint45d_post, checkpoint45b_post, checkpoint46b_pre, chkpt44a_pre, release1-branch-end, release1_final_v1, ecco_c44_e19, checkpoint46, ecco_c44_e20, checkpoint44, release1_p13, ecco_c44_e18, checkpoint44f_post, release1_p17, release1_b1, checkpoint44b_post, chkpt44c_post, chkpt44d_post, checkpoint42, release1_p9, release1_p8, checkpoint43, checkpoint46g_pre, release1_p2, release1_p3, release1_p4, release1_p6, checkpoint46a_post, chkpt44a_post, checkpoint44b_pre, release1_p1, checkpoint46a_pre, ecco-branch-mod1, checkpoint45c_post, release1_p5, checkpoint44e_pre, chkpt44c_pre, checkpoint40pre9, release1_p7, ecco_ice2, ecco_ice1, checkpoint46b_post, checkpoint46d_post, ecco-branch-mod2, checkpoint46g_post, checkpoint45a_post, checkpoint46c_pre, ecco-branch-mod3, ecco-branch-mod5, ecco_c44_e22, release1_beta1, ecco_c44_e23, release1-branch_branchpoint, checkpoint46c_post, checkpoint40, checkpoint45, checkpoint46h_post, release1_chkpt44d_post, ecco_c44_e25, checkpoint41
Branch point for: c24_e25_ice, ecco-branch, release1_coupled, icebear, release1_final, release1-branch, release1, release1_50yr
Changes since 1.4: +2 -2 lines
Added declaration for i2. Caused failure on SUNs with -u.

1 adcroft 1.5 C $Header: /u/gcmpack/models/MITgcmUV/pkg/aim/phy_radiat.F,v 1.4 2001/06/18 17:39:58 cnh Exp $
2 cnh 1.3 C $Name: $
3 adcroft 1.2
4 cnh 1.3 SUBROUTINE SOL_OZ (SOLC,TYEAR,FSOL,OZONE,myThid)
5 adcroft 1.2
6     C--
7     C-- SUBROUTINE SOL_OZ (SOLC,TYEAR,FSOL,OZONE)
8     C--
9     C-- Purpose: Compute the flux of incoming solar radiation
10     C-- and a climatological ozone profile for SW absorption
11     C-- Input: SOLC = solar constant (area averaged)
12     C-- TYEAR = time as fraction of year (0-1, 0 = 1jan.h00)
13     C-- Output: FSOL = flux of incoming solar radiation (2-dim)
14     C-- OZONE = strat. ozone as fraction of global mean (2-dim)
15     C--
16    
17    
18     IMPLICIT rEAL*8 (A-H,O-Z)
19 cnh 1.3 INTEGER myThid
20 adcroft 1.2
21 cnh 1.3 #include "EEPARAMS.h"
22 adcroft 1.2
23     C Resolution parameters
24     C
25     #include "atparam.h"
26     #include "atparam1.h"
27     C
28 cnh 1.3 INTEGER NLON, NLAT, NLEV, NGP
29 cnh 1.4 INTEGER I, J, I2
30 adcroft 1.2 PARAMETER ( NLON=IX, NLAT=IL, NLEV=KX, NGP=NLON*NLAT )
31     C
32     C Constants + functions of sigma and latitude
33     C
34     #include "com_physcon.h"
35    
36     REAL FSOL(NLON,NLAT), OZONE(NLON,NLAT)
37     C
38     C ALPHA = year phase ( 0 - 2pi, 0 = winter solstice = 22dec.h00 )
39     ALPHA=4.*ASIN(1.)*(TYEAR+10./365.)
40    
41     CSR1=-0.796*COS(ALPHA)
42     CSR2= 0.147*COS(2*ALPHA)-0.477
43     COZ1= 0.0
44     C COZ1= 0.2*SIN(ALPHA)
45     COZ2= 0.3
46 cnh 1.3
47 adcroft 1.2 C
48     DO J=1,NLAT
49 cnh 1.4 DO I=1,NLON
50     I2=J
51     I2=NLON*(J-1)+I
52     FSOL(I,J)=
53     & SOLC*MAX(0.,1.0+CSR1*FMU(I2,1,myThid)+CSR2*FMU(I2,2,myThid))
54     OZONE(I,J)=1.0+COZ1*FMU(I2,1,myThid)+COZ2*FMU(I2,2,myThid)
55     ENDDO
56 adcroft 1.2 ENDDO
57 cnh 1.4 C DO J=1,NLAT
58     C FSOL(1,J)=
59     C & SOLC*MAX(0.,1.0+CSR1*FMU(J,1,myThid)+CSR2*FMU(J,2,myThid))
60     C OZONE(1,J)=1.0+COZ1*FMU(J,1,myThid)+COZ2*FMU(J,2,myThid)
61     C DO I=2,NLON
62     C FSOL(I,J)=FSOL(1,J)
63     C OZONE(I,J)=OZONE(1,J)
64     C ENDDO
65     C ENDDO
66 adcroft 1.2 C
67     RETURN
68     END
69    
70    
71     SUBROUTINE RADSW (PSA,QA,RH,
72 cnh 1.3 * FSOL,OZONE,ALB,TAU,
73     * CLOUDC,FTOP,FSFC,DFABS,
74     I myThid)
75 adcroft 1.2 C--
76     C-- SUBROUTINE RADSW (PSA,QA,RH,
77     C-- * FSOL,OZONE,ALB,
78     C-- * CLOUDC,FTOP,FSFC,DFABS)
79     C--
80     C-- Purpose: Compute the absorption of shortwave radiation and
81     C-- initialize arrays for longwave-radiation routines
82     C-- Input: PSA = norm. surface pressure [p/p0] (2-dim)
83     C-- QA = specific humidity [g/kg] (3-dim)
84     C-- RH = relative humidity (3-dim)
85     C-- FSOL = flux of incoming solar radiation (2-dim)
86     C-- OZONE = strat. ozone as fraction of global mean (2-dim)
87     C-- ALB = surface albedo (2-dim)
88     C-- Output: CLOUDC = total cloud cover (2-dim)
89     C-- FTOP = net downw. flux of sw rad. at the atm. top (2-dim)
90     C-- FSFC = net downw. flux of sw rad. at the surface (2-dim)
91     C-- DFABS = flux of sw rad. absorbed by each atm. layer (3-dim)
92     C--
93    
94    
95     IMPLICIT rEAL*8 (A-H,O-Z)
96 cnh 1.3 INTEGER myThid
97 adcroft 1.2
98    
99     C Resolution parameters
100     C
101     #include "atparam.h"
102     #include "atparam1.h"
103 cnh 1.3 #include "EEPARAMS.h"
104 adcroft 1.2 #include "Lev_def.h"
105     C
106 cnh 1.3 INTEGER NLON, NLAT, NLEV, NGP
107     INTEGER K, J
108 adcroft 1.2 PARAMETER ( NLON=IX, NLAT=IL, NLEV=KX, NGP=NLON*NLAT )
109     C
110     C Constants + functions of sigma and latitude
111     C
112     #include "com_physcon.h"
113     C
114     C Radiation parameters
115     C
116     #include "com_radcon.h"
117     C
118     REAL PSA(NGP), QA(NGP,NLEV), RH(NGP,NLEV),
119 cnh 1.3 * FSOL(NGP), OZONE(NGP), ALB(NGP), TAU(NGP,NLEV)
120 adcroft 1.2
121     REAL CLOUDC(NGP), FTOP(NGP), FSFC(NGP), DFABS(NGP,NLEV)
122    
123     REAL FLUX(NGP), FREFL(NGP), TAUOZ(NGP)
124     INTEGER NL1(NGP)
125     Cchdbg
126     INTEGER Npas
127     SAVE npas
128     LOGICAL Ifirst
129     SAVE Ifirst
130     DATA Ifirst /.TRUE./
131     REAL clsum(NGP)
132     SAVE clsum
133     REAL ABWLW1
134     cchdbg
135     C
136     DO J=1,NGP
137 cnh 1.3 NL1(J)=NLEVxy(J,myThid)-1
138 adcroft 1.2 ENDDO
139     C
140     C-- 1. Cloud cover:
141     C defined as a linear fun. of the maximum relative humidity
142     C in all tropospheric layers above PBL:
143     C CLOUDC = 0 for RHmax < RHCL1, = 1 for RHmax > RHCL2.
144     C This value is reduced by a factor (Qbase/QACL) if the
145     C cloud-base absolute humidity Qbase < QACL.
146     C
147     DRHCL=RHCL2-RHCL1
148     RCL=1./(DRHCL*QACL)
149     C
150     DO 122 J=1,NGP
151     CLOUDC(J)=0.
152     122 CONTINUE
153     C
154     DO 123 K=1,NLEV
155     DO 123 J=1,NGP
156     DFABS(J,K)=0.
157     123 CONTINUE
158 cnh 1.3
159 adcroft 1.2 C
160     DO 124 J=1,NGP
161     DO 124 K=2,NL1(J)
162     CLOUDC(J)=MAX(CLOUDC(J),(RH(J,K)-RHCL1))
163     124 CONTINUE
164     C
165     DO 126 J=1,NGP
166     IF ( NL1(J) .GT. 0 ) THEN
167     CLOUDC(J)=MIN(CLOUDC(J),DRHCL)*MIN(QA(J,NL1(J)),QACL)*RCL
168     ENDIF
169     cchdbg *******************************************
170     cchdbg CLOUDC(J)=MIN(CLOUDC(J),DRHCL)/DRHCL
171     cchdbg *******************************************
172     clear sky experiment
173 cnh 1.3 C cloudc(j) = 0.
174 adcroft 1.2 126 CONTINUE
175     C
176     C
177     C-- 2. Shortwave transmissivity:
178     C function of layer mass, ozone (in the statosphere),
179     C abs. humidity and cloud cover (in the troposphere)
180     C
181     DO 202 J=1,NGP
182     TAU(J,1)=EXP(-ABSSW*PSA(J)*DSIG(1))
183     TAUOZ(J)=EXP(-EPSSW*OZONE(J)*PSA(J))
184     202 CONTINUE
185     C
186     chhh WRITE(0,*) ' Hello from RADSW'
187     DO 204 J=1,NGP
188     DO 204 K=2,NL1(J)
189     TAU(J,K)=EXP(-(ABSSW+ABWSW*QA(J,K)
190     * +ABCSW*CLOUDC(J)*QA(J,NL1(J)))*PSA(J)*DSIG(K))
191     204 CONTINUE
192    
193     DO 206 J=1,NGP
194 cnh 1.3 IF ( NLEVxy(J,myThid) .GT. 0 ) THEN
195     TAU(J,NLEVxy(J,myThid))=EXP(-(ABSSW+ABWSW*QA(J,NLEVxy(J,myThid)))
196     * *PSA(J)*DSIG(NLEVxy(J,myThid)))
197 adcroft 1.2 ENDIF
198     206 CONTINUE
199     C
200     C--- 3. Shortwave downward flux
201     C
202     C 3.1 Absorption in the stratosphere
203     C
204     DO 312 J=1,NGP
205     FLUX(J)=TAU(J,1)*TAUOZ(J)*FSOL(J)
206     DFABS(J,1)=FSOL(J)-FLUX(J)
207     312 CONTINUE
208 cnh 1.3
209     C RETURN
210    
211 adcroft 1.2 C
212     C 3.2 Reflection at the top of the troposphere
213     C (proportional to cloud cover).
214     C
215     DO 322 J=1,NGP
216     FREFL(J)=ALBCL*CLOUDC(J)*FLUX(J)
217     FTOP(J) =FSOL(J)-FREFL(J)
218     FLUX(J) =FLUX(J)-FREFL(J)
219     322 CONTINUE
220     C
221     C 3.3 Absorption in the troposphere
222     C
223     DO 332 J=1,NGP
224 cnh 1.3 DO 332 K=2,NLEVxy(J,myThid)
225 adcroft 1.2 DFABS(J,K)=FLUX(J)
226     FLUX(J)=TAU(J,K)*FLUX(J)
227     DFABS(J,K)=DFABS(J,K)-FLUX(J)
228     332 CONTINUE
229 cnh 1.3
230     Cxx RETURN
231    
232 adcroft 1.2 C
233     C--- 4. Shortwave upward flux
234     C
235     C 4.1 Absorption and reflection at the surface
236     C
237     DO 412 J=1,NGP
238     FREFL(J)=ALB(J)*FLUX(J)
239     FSFC(J) =FLUX(J)-FREFL(J)
240     FLUX(J) =FREFL(J)
241     412 CONTINUE
242     C
243     C 4.2 Absorption in the atmosphere
244     C
245     DO 422 J=1,NGP
246 cnh 1.3 DO 422 K=NLEVxy(J,myThid),1,-1
247 adcroft 1.2 DFABS(J,K)=DFABS(J,K)+FLUX(J)
248     FLUX(J)=TAU(J,K)*FLUX(J)
249     DFABS(J,K)=DFABS(J,K)-FLUX(J)
250     422 CONTINUE
251 cnh 1.3
252     Cxx RETURN
253    
254 adcroft 1.2 C
255     C 4.3 Absorbed solar radiation = incoming - outgoing
256     C
257     DO 432 J=1,NGP
258     FTOP(J)=FTOP(J)-FLUX(J)
259     432 CONTINUE
260 cnh 1.3
261     C RETURN
262 adcroft 1.2 cdj
263     c write(0,*)'position j=20'
264     c j=20
265     c write(0,*)'ftop fsfc ftop-fsfc'
266     c write(0,*)ftop(j),fsfc(j),ftop(j)-fsfc(j)
267     c write(0,*)
268     c write(0,*)'k dfabs'
269     c do k = 1, nlevxy(j)
270     c write(0,*)k,dfabs(j,k)
271     c enddo
272     c write(0,*)'sum dfabs'
273     c write(0,*)sum(dfabs(j,:))
274     cdj
275     C
276     C--- 5. Initialization of longwave radiation model
277     C
278     C 5.1 Longwave transmissivity:
279     C function of layer mass, abs. humidity and cloud cover.
280     C
281     DO 512 J=1,NGP
282     TAU(J,1)=EXP(-ABSLW*PSA(J)*DSIG(1))
283     512 CONTINUE
284 cnh 1.3
285 adcroft 1.2 C
286     DO 514 J=1,NGP
287     DO 514 K=2,NL1(J)
288     TAU(J,K)=EXP(-(ABSLW+ABWLW*QA(J,K)
289     * +ABCLW*CLOUDC(J)*QA(J,NL1(J)))*PSA(J)*DSIG(K))
290     514 CONTINUE
291 cnh 1.3
292     C RETURN
293 adcroft 1.2 C
294     cchdbg ***************************************************
295     c ABCLW1=0.15
296     c DO 514 J=1,NGP
297     c DO 514 K=2,NL1(J)-1
298     c TAU(J,K)=EXP(-(ABSLW+ABWLW*QA(J,K)
299     c * +ABCLW1*CLOUDC(J)*QA(J,NL1(J)))*PSA(J)*DSIG(K))
300     c 514 CONTINUE
301     C
302     c DO 515 J=1,NGP
303     c DO 515 K=NL1(J),NL1(J)
304     c TAU(J,K)=EXP(-(ABSLW+ABWLW*QA(J,K)
305     c * +ABCLW*CLOUDC(J))*PSA(J)*DSIG(K))
306     c 515 CONTINUE
307     cchdbg ************************************************************
308     C
309     C *********************************************************************
310     C *********************************************************************
311     C *****************************************************************
312     cchdbg
313     c if(Ifirst) then
314     c npas=0
315     c do J=1,NGP
316     c clsum(J)=0.
317     c enddo
318     c ifirst=.FALSE.
319     c ENDIF
320     C
321     c npas=npas+1
322     c DO J=1,NGP
323     c clsum(J)=clsum(J)+ABCLW*CLOUDC(J)*QA(J,NL1(J))/5760.
324     c ENDDO
325     C
326     c IF(npas.eq.5760) then
327     c open(73,file='transmoy',form='unformatted')
328     c write(73) clsum
329     c close(73)
330     c ENDIF
331     Cchdbg
332     C
333     C *********************************************************************
334     C *********************************************************************
335    
336 cnh 1.3 C RETURN
337    
338 adcroft 1.2 ABWLW1=0.7
339     DO 516 J=1,NGP
340 cnh 1.3 IF ( NLEVxy(J,myThid) .GT. 0 ) THEN
341     TAU(J,NLEVxy(J,myThid))=EXP(-(ABSLW+ABWLW*QA(J,NLEVxy(J,myThid)))*PSA(J)
342     cchdbg TAU(J,NLEVxy(J,myThid))=EXP(-(ABSLW+ABWLW1*QA(J,NLEVxy(J,myThid)))*PSA(J)
343     * *DSIG(NLEVxy(J,myThid)))
344 adcroft 1.2 ENDIF
345     516 CONTINUE
346    
347     C
348     C---
349     RETURN
350     END
351    
352    
353     SUBROUTINE RADLW (IMODE,TA,TS,ST4S,
354 cnh 1.3 & TAU,ST4A,
355     * FTOP,FSFC,DFABS,FDOWN,
356     I myThid)
357 adcroft 1.2 C--
358     C-- SUBROUTINE RADLW (IMODE,TA,TS,ST4S,
359     C-- * FTOP,FSFC,DFABS)
360     C--
361     C-- Purpose: Compute the absorption of longwave radiation
362     C-- Input: IMODE = index for operation mode (see below)
363     C-- TA = absolute temperature (3-dim)
364     C-- TS = surface temperature (2-dim) [if IMODE=1]
365     C-- ST4S = surface blackbody emission (2-dim) [if IMODE=2]
366     C-- Output: ST4S = surface blackbody emission (2-dim) [if IMODE=1]
367     C-- FTOP = outgoing flux of lw rad. at the atm. top (2-dim)
368     C-- FSFC = net upw. flux of lw rad. at the surface (2-dim)
369     C-- DFABS = flux of lw rad. absorbed by each atm. layer (3-dim)
370     C-- FDOWN = downward flux of lw rad. at the surface (2-dim)
371     C--
372    
373    
374     IMPLICIT rEAL*8 (A-H,O-Z)
375 cnh 1.3 INTEGER myThid
376 adcroft 1.2
377     C Resolution parameters
378     C
379     #include "atparam.h"
380     #include "atparam1.h"
381 cnh 1.3 #include "EEPARAMS.h"
382 adcroft 1.2 #include "Lev_def.h"
383     C
384 cnh 1.3 INTEGER NLON, NLAT, NLEV, NGP
385     INTEGER K, J
386 adcroft 1.2 PARAMETER ( NLON=IX, NLAT=IL, NLEV=KX, NGP=NLON*NLAT )
387     C
388     C Constants + functions of sigma and latitude
389     C
390     #include "com_physcon.h"
391     C
392     C Radiation parameters
393     C
394     #include "com_radcon.h"
395     C
396 cnh 1.3 REAL TA(NGP,NLEV), TS(NGP), ST4S(NGP),
397     & TAU(NGP,NLEV), ST4A(NGP,NLEV,2)
398 adcroft 1.2
399     REAL FTOP(NGP), FSFC(NGP), DFABS(NGP,NLEV)
400     REAL FDOWN(NGP)
401    
402     REAL FLUX(NGP), BRAD(NGP), STCOR(NGP)
403     INTEGER NL1(NGP)
404 adcroft 1.5 INTEGER IMODE, J0, Jl, I2
405 adcroft 1.2 C
406     Cchdbg
407     INteger npas
408     SAVE npas
409     LOGICAL Ifirst
410     SAVE IFIRST
411     DATA Ifirst/.TRUE./
412     REAL FluxMoy(NGP)
413     REAL ST4SMoy(NGP)
414     SAVE FluxMoy, ST4SMoy
415     Cchdbg
416 cnh 1.3
417 adcroft 1.2 DO J=1,NGP
418 cnh 1.3 NL1(J)=NLEVxy(J,myThid)-1
419 adcroft 1.2 ENDDO
420 cnh 1.3
421 adcroft 1.2 C
422     DO K=1,NLEV
423     DO J=1,NGP
424     DFABS(J,K)=0.
425     ENDDO
426     ENDDO
427 cnh 1.3
428 adcroft 1.2 C
429     C--- 1. Blackbody emission from atmospheric full and half levels.
430     C Temperature is interpolated as a linear function of ln sigma.
431     C At the lower boundary, the emission is linearly extrapolated;
432     C at the upper boundary, the atmosphere is assumed isothermal.
433     C
434     DO 102 J=1,NGP
435 cnh 1.3 DO 102 K=1,NLEVxy(J,myThid)
436 adcroft 1.2 ST4A(J,K,1)=TA(J,K)*TA(J,K)
437     ST4A(J,K,1)=SBC*ST4A(J,K,1)*ST4A(J,K,1)
438     102 CONTINUE
439     C
440     DO 104 J=1,NGP
441     DO 104 K=1,NL1(J)
442     ST4A(J,K,2)=TA(J,K)+WVI(K,2)*(TA(J,K+1)-TA(J,K))
443     ST4A(J,K,2)=ST4A(J,K,2)*ST4A(J,K,2)
444     ST4A(J,K,2)=SBC*ST4A(J,K,2)*ST4A(J,K,2)
445     104 CONTINUE
446 cnh 1.3
447 adcroft 1.2 C
448     DO 106 J=1,NGP
449 cnh 1.3 IF ( NLEVxy(J,myThid) .GT. 0 ) THEN
450     ST4A(J,NLEVxy(J,myThid),2)=2.*ST4A(J,NLEVxy(J,myThid),1)-ST4A(J,NL1(J),2)
451 adcroft 1.2 ENDIF
452     106 CONTINUE
453     C
454     C--- 2. Empirical stratospheric correction
455     C
456     COR0= -13.
457     COR1= 0.
458     COR2= 24.
459     C
460     J0=0
461 cnh 1.4 DO JL=1,nlat
462     C CORR=COR0+COR1*FMU(JL,1,myThid)+COR2*FMU(JL,2,myThid)
463 adcroft 1.2 DO J=J0+1,J0+NLON
464 cnh 1.4 I2=JL
465     I2=J
466     STCOR(J)=COR0+COR1*FMU(I2,1,myThid)+COR2*FMU(I2,2,myThid)
467     C STCOR(J)=CORR
468 adcroft 1.2 ENDDO
469     J0=J0+NLON
470     ENDDO
471     C
472     C--- 3. Emission ad absorption of longwave downward flux.
473     C Downward emission is an average of the emission from the full level
474     C and the half-level below, weighted according to the transmissivity
475     C of the layer.
476     C
477     C 3.1 Stratosphere
478     C
479     DO 312 J=1,NGP
480     BRAD(J)=ST4A(J,1,2)+TAU(J,1)*(ST4A(J,1,1)-ST4A(J,1,2))
481     FLUX(J)=(1.-TAU(J,1))*BRAD(J)
482     DFABS(J,1)=STCOR(J)-FLUX(J)
483     312 CONTINUE
484     C
485     C 3.2 Troposphere
486     C
487     DO 322 J=1,NGP
488 cnh 1.3 DO 322 K=2,NLEVxy(J,myThid)
489 adcroft 1.2 DFABS(J,K)=FLUX(J)
490     BRAD(J)=ST4A(J,K,2)+TAU(J,K)*(ST4A(J,K,1)-ST4A(J,K,2))
491     FLUX(J)=TAU(J,K)*(FLUX(J)-BRAD(J))+BRAD(J)
492     DFABS(J,K)=DFABS(J,K)-FLUX(J)
493     322 CONTINUE
494     C
495     C--- 4. Emission ad absorption of longwave upward flux
496     C Upward emission is an average of the emission from the full level
497     C and the half-level above, weighted according to the transmissivity
498     C of the layer (for the top layer, full-level emission is used).
499     C Surface lw emission in the IR window goes directly into FTOP.
500     C
501     C 4.1 Surface
502     C
503     IF (IMODE.LE.1) THEN
504     DO 412 J=1,NGP
505     ST4S(J)=TS(J)*TS(J)
506     ST4S(J)=SBC*ST4S(J)*ST4S(J)
507     412 CONTINUE
508     ENDIF
509     C
510     C **************************************************************
511     Cchdbg
512     if(ifirst) then
513     DO J=1,NGP
514     ST4SMoy(J)=0.
515     FluxMoy(J)=0.
516     ENDDO
517     npas=0.
518     ifirst=.FALSE.
519     endif
520     C
521     npas=npas+1
522     DO 413 J=1,NGP
523     ST4SMoy(J)=ST4SMoy(J)+ ST4S(J)
524     FluxMoy(J)=FluxMoy(J)+ Flux(J)
525     413 CONTINUE
526     C
527     if(npas.eq.5760) then
528     DO J=1,NGP
529     ST4SMoy(J)=ST4SMoy(J)/float(npas)
530     FluxMoy(J)=FluxMoy(J)/float(npas)
531     ENDDO
532     open(73,file='ST4Smoy',form='unformatted')
533     write(73) ST4SMoy
534     close(73)
535     open(74,file='FluxMoy',form='unformatted')
536     write(74) FluxMoy
537     close(74)
538     ENDIF
539     Cchdbg
540     C ****************************************************************
541     C
542     C
543     DO 414 J=1,NGP
544     FSFC(J)=ST4S(J)-FLUX(J)
545     FDOWN(J)=FLUX(J)
546     FTOP(J)=EPSLW*ST4S(J)
547     FLUX(J)=ST4S(J)-FTOP(J)
548     414 CONTINUE
549     C
550     C 4.2 Troposphere
551     C
552     DO 422 J=1,NGP
553 cnh 1.3 DO 422 K=NLEVxy(J,myThid),2,-1
554 adcroft 1.2 DFABS(J,K)=DFABS(J,K)+FLUX(J)
555     BRAD(J)=ST4A(J,K-1,2)+TAU(J,K)*(ST4A(J,K,1)-ST4A(J,K-1,2))
556     FLUX(J)=TAU(J,K)*(FLUX(J)-BRAD(J))+BRAD(J)
557     DFABS(J,K)=DFABS(J,K)-FLUX(J)
558     422 CONTINUE
559     C
560     C 4.3 Stratosphere
561     C
562     DO 432 J=1,NGP
563     DFABS(J,1)=DFABS(J,1)+FLUX(J)
564     FLUX(J)=TAU(J,1)*(FLUX(J)-ST4A(J,1,1))+ST4A(J,1,1)
565     DFABS(J,1)=DFABS(J,1)-FLUX(J)
566     432 CONTINUE
567     C
568     C 4.4 Outgoing longwave radiation
569     C
570     DO 442 J=1,NGP
571     cdj FTOP(J)=FTOP(J)+FLUX(J)
572     FTOP(J)=FTOP(J)+FLUX(J)-STCOR(J)
573     442 CONTINUE
574     cdj
575     c write(0,*)'position j=20'
576     c j=20
577     c write(0,*)'ftop fsfc ftop-fsfc'
578     c write(0,*)ftop(j),fsfc(j),ftop(j)-fsfc(j)
579     c write(0,*)
580     c write(0,*)'k dfabs'
581     c do k = 1, nlevxy(j)
582     c write(0,*)k,dfabs(j,k)
583     c enddo
584     c write(0,*)'sum dfabs'
585     c write(0,*)sum(dfabs(j,:))
586     c open(74,file='ftop0',form='unformatted',status='unknown')
587     c write(74) ftop
588     c open(75,file='stcor',form='unformatted',status='unknown')
589     c write(75) stcor
590     c stop
591     cdj
592     C
593     C---
594     RETURN
595     END

  ViewVC Help
Powered by ViewVC 1.1.22