source: LMDZ5/trunk/libf/phylmd/radlwsw_aero.F90 @ 2093

Last change on this file since 2093 was 2009, checked in by Laurent Fairhead, 11 years ago

Rajout de la clé CPP_RRTM pour la compilation avec RRTM


Added CPP_RRTM key for RRTM compilation

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:keywords set to Author Date Id Revi
File size: 33.2 KB
Line 
1!
2! $Id: radlwsw_aero.F90 2009 2014-04-08 08:48:17Z fhourdin $
3!
4SUBROUTINE radlwsw_aero( &
5   dist, rmu0, fract, &
6   paprs, pplay,tsol,alb1, alb2,&
7   t,q,wo,&
8   cldfra, cldemi, cldtaupd,&
9   ok_ade, ok_aie,&
10   tau_aero, piz_aero, cg_aero,&
11   cldtaupi, new_aod, &
12   heat,heat0,cool,cool0,radsol,albpla,&
13   topsw,toplw,solsw,sollw,&
14   sollwdown,&
15   topsw0,toplw0,solsw0,sollw0,&
16   lwdn0, lwdn, lwup0, lwup,&
17   swdn0, swdn, swup0, swup,&
18   topswad_aero, solswad_aero,&
19   topswai_aero, solswai_aero, &
20   topswad0_aero, solswad0_aero,&
21   topsw_aero, topsw0_aero,&
22   solsw_aero, solsw0_aero,qsat,flwc,fiwc)
23
24
25
26  USE DIMPHY
27  USE comgeomphy
28  USE write_field_phy
29!    modules necessaires au rayonnement
30!    -----------------------------------------
31! USE YOMCST   , ONLY : RG       ,RD       ,RTT      ,RPI
32! USE YOERAD   , ONLY : NSW      ,LRRTM    ,LINHOM   , LCCNL,LCCNO,
33! USE YOERAD   , ONLY : NSW      ,LRRTM    ,LCCNL    ,LCCNO ,&
34! NSW mis dans .def MPL 20140211
35#ifdef CPP_RRTM
36!  USE YOERAD   , ONLY : LRRTM    ,LCCNL    ,LCCNO ,&
37!      NRADIP   , NRADLP , NICEOPT, NLIQOPT ,RCCNLND  , RCCNSEA
38!  USE YOELW    , ONLY : NSIL     ,NTRA     ,NUA      ,TSTAND   ,XP
39!  USE YOESW    , ONLY : RYFWCA   ,RYFWCB   ,RYFWCC   ,RYFWCD,&   
40!      RYFWCE   ,RYFWCF   ,REBCUA   ,REBCUB   ,REBCUC,&   
41!      REBCUD   ,REBCUE   ,REBCUF   ,REBCUI   ,REBCUJ,& 
42!      REBCUG   ,REBCUH   ,RHSAVI   ,RFULIO   ,RFLAA0,& 
43!      RFLAA1   ,RFLBB0   ,RFLBB1   ,RFLBB2   ,RFLBB3,& 
44!      RFLCC0   ,RFLCC1   ,RFLCC2   ,RFLCC3   ,RFLDD0,& 
45!      RFLDD1   ,RFLDD2   ,RFLDD3   ,RFUETA   ,RASWCA,&
46!      RASWCB   ,RASWCC   ,RASWCD   ,RASWCE   ,RASWCF
47!&    RASWCB   ,RASWCC   ,RASWCD   ,RASWCE   ,RASWCF, RLINLI
48!  USE YOERDU   , ONLY : NUAER  ,NTRAER ,REPLOG ,REPSC  ,REPSCW ,DIFF
49!  USE YOETHF   , ONLY : RTICE
50!  USE YOERRTWN , ONLY : DELWAVE   ,TOTPLNK     
51  USE YOMPHY3  , ONLY : RII0
52#endif
53     
54  IMPLICIT NONE
55
56  !======================================================================
57  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
58  ! Objet: interface entre le modele et les rayonnements
59  ! Arguments:
60  ! dist-----input-R- distance astronomique terre-soleil
61  ! rmu0-----input-R- cosinus de l'angle zenithal
62  ! fract----input-R- duree d'ensoleillement normalisee
63  ! co2_ppm--input-R- concentration du gaz carbonique (en ppm)
64  ! solaire--input-R- constante solaire (W/m**2)
65  ! paprs----input-R- pression a inter-couche (Pa)
66  ! pplay----input-R- pression au milieu de couche (Pa)
67  ! tsol-----input-R- temperature du sol (en K)
68  ! alb1-----input-R- albedo du sol(entre 0 et 1) dans l'interval visible
69  ! alb2-----input-R- albedo du sol(entre 0 et 1) dans l'interval proche infra-rouge
70  ! t--------input-R- temperature (K)
71  ! q--------input-R- vapeur d'eau (en kg/kg)
72  ! wo-------input-R- contenu en ozone (en kg/kg) correction MPL 100505
73  ! cldfra---input-R- fraction nuageuse (entre 0 et 1)
74  ! cldtaupd---input-R- epaisseur optique des nuages dans le visible (present-day value)
75  ! cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
76  ! ok_ade---input-L- apply the Aerosol Direct Effect or not?
77  ! ok_aie---input-L- apply the Aerosol Indirect Effect or not?
78  ! tau_ae, piz_ae, cg_ae-input-R- aerosol optical properties (calculated in aeropt.F)
79  ! cldtaupi-input-R- epaisseur optique des nuages dans le visible
80  !                   calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
81  !                   droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
82  !                   it is needed for the diagnostics of the aerosol indirect radiative forcing     
83  !
84  ! heat-----output-R- echauffement atmospherique (visible) (K/jour)
85  ! cool-----output-R- refroidissement dans l'IR (K/jour)
86  ! radsol---output-R- bilan radiatif net au sol (W/m**2) (+ vers le bas)
87  ! albpla---output-R- albedo planetaire (entre 0 et 1)
88  ! topsw----output-R- flux solaire net au sommet de l'atm.
89  ! toplw----output-R- ray. IR montant au sommet de l'atmosphere
90  ! solsw----output-R- flux solaire net a la surface
91  ! sollw----output-R- ray. IR montant a la surface
92  ! solswad---output-R- ray. solaire net absorbe a la surface (aerosol dir)
93  ! topswad---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol dir)
94  ! solswai---output-R- ray. solaire net absorbe a la surface (aerosol ind)
95  ! topswai---output-R- ray. solaire absorbe au sommet de l'atm. (aerosol ind)
96  !
97  ! ATTENTION: swai and swad have to be interpreted in the following manner:
98  ! ---------
99  ! ok_ade=F & ok_aie=F -both are zero
100  ! ok_ade=T & ok_aie=F -aerosol direct forcing is F_{AD} = topsw-topswad
101  !                        indirect is zero
102  ! ok_ade=F & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
103  !                        direct is zero
104  ! ok_ade=T & ok_aie=T -aerosol indirect forcing is F_{AI} = topsw-topswai
105  !                        aerosol direct forcing is F_{AD} = topswai-topswad
106  !
107 
108  !======================================================================
109 
110  ! ====================================================================
111  ! Adapte au modele de chimie INCA par Celine Deandreis & Anne Cozic -- 2009
112  ! 1 = ZERO   
113  ! 2 = AER total   
114  ! 3 = NAT   
115  ! 4 = BC   
116  ! 5 = SO4   
117  ! 6 = POM   
118  ! 7 = DUST   
119  ! 8 = SS   
120  ! 9 = NO3   
121  !
122  ! ====================================================================
123  include "YOETHF.h"
124  include "YOMCST.h"
125  include "clesphys.h"
126
127! Input arguments
128! REAL,    INTENT(in)  :: solaire
129  REAL,    INTENT(in)  :: dist
130  REAL,    INTENT(in)  :: rmu0(KLON), fract(KLON)
131  REAL,    INTENT(in)  :: paprs(KLON,KLEV+1), pplay(KLON,KLEV)
132  REAL,    INTENT(in)  :: alb1(KLON), alb2(KLON),tsol(KLON)
133  REAL,    INTENT(in)  :: t(KLON,KLEV), q(KLON,KLEV), wo(KLON,KLEV)
134  LOGICAL, INTENT(in)  :: ok_ade, ok_aie                                 ! switches whether to use aerosol direct (indirect) effects or not
135  REAL,    INTENT(in)  :: cldfra(KLON,KLEV), cldemi(KLON,KLEV), cldtaupd(KLON,KLEV)
136  REAL,    INTENT(in)  :: tau_aero(KLON,KLEV,9,2)                        ! aerosol optical properties (see aeropt.F)
137  REAL,    INTENT(in)  :: piz_aero(KLON,KLEV,9,2)                        ! aerosol optical properties (see aeropt.F)
138  REAL,    INTENT(in)  :: cg_aero(KLON,KLEV,9,2)                         ! aerosol optical properties (see aeropt.F)
139  REAL,    INTENT(in)  :: cldtaupi(KLON,KLEV)                            ! cloud optical thickness for pre-industrial aerosol concentrations
140!MPL input supplementaires pour RECMWFL
141! flwc, fiwc = Liquid Water Content & Ice Water Content (kg/kg)
142  REAL*8  GEMU(klon)
143  REAL*8  qsat(klon,klev),flwc(klon,klev),fiwc(klon,klev)
144!MPL input RECMWFL:
145!Tableaux aux niveaux inverses pour respecter convention Arpege
146  REAL*8  paprs_i(klon,klev+1)
147  REAL*8  pplay_i(klon,klev)
148  REAL*8  cldfra_i(klon,klev)
149  REAL*8  POZON_i(kdlon,kflev)
150!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6     
151  REAL*8  PAER_i(kdlon,kflev,6)
152  REAL*8  PDP_i(klon,klev)
153  REAL*8  t_i(klon,klev),q_i(klon,klev),qsat_i(klon,klev)
154  REAL*8  flwc_i(klon,klev),fiwc_i(klon,klev)
155! new_aod: flag pour retrouver les resultats exacts de l'AR4 dans le cas ou l'on ne travaille qu'avec les sulfates
156  LOGICAL, INTENT(in)  :: new_aod
157  LOGICAL lldebug
158
159! Output arguments
160  REAL,    INTENT(out) :: heat(KLON,KLEV), cool(KLON,KLEV)
161  REAL,    INTENT(out) :: heat0(KLON,KLEV), cool0(KLON,KLEV)
162  REAL,    INTENT(out) :: radsol(KLON), topsw(KLON), toplw(KLON)
163  REAL,    INTENT(out) :: solsw(KLON), sollw(KLON), albpla(KLON)
164  REAL,    INTENT(out) :: topsw0(KLON), toplw0(KLON), solsw0(KLON), sollw0(KLON)
165  REAL,    INTENT(out) :: sollwdown(KLON)
166  REAL,    INTENT(out) :: swdn(KLON,kflev+1),swdn0(KLON,kflev+1)
167  REAL,    INTENT(out) :: swup(KLON,kflev+1),swup0(KLON,kflev+1)
168  REAL,    INTENT(out) :: lwdn(KLON,kflev+1),lwdn0(KLON,kflev+1)
169  REAL,    INTENT(out) :: lwup(KLON,kflev+1),lwup0(KLON,kflev+1)
170  REAL,    INTENT(out) :: topswad_aero(KLON), solswad_aero(KLON)         ! output: aerosol direct forcing at TOA and surface
171  REAL,    INTENT(out) :: topswai_aero(KLON), solswai_aero(KLON)         ! output: aerosol indirect forcing atTOA and surface
172  REAL, DIMENSION(klon), INTENT(out)    :: topswad0_aero
173  REAL, DIMENSION(klon), INTENT(out)    :: solswad0_aero
174  REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw_aero
175  REAL, DIMENSION(kdlon,9), INTENT(out) :: topsw0_aero
176  REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw_aero
177  REAL, DIMENSION(kdlon,9), INTENT(out) :: solsw0_aero
178! --------- output RECMWFL
179!  ZEMTD (KPROMA,KLEV+1)         ; TOTAL DOWNWARD LONGWAVE EMISSIVITY
180!  ZEMTU (KPROMA,KLEV+1)         ; TOTAL UPWARD   LONGWAVE EMISSIVITY
181!  ZTRSO (KPROMA,KLEV+1)         ; TOTAL SHORTWAVE TRANSMISSIVITY
182!  ZTH   (KPROMA,KLEV+1)         ; HALF LEVEL TEMPERATURE
183!  ZCTRSO(KPROMA,2)              ; CLEAR-SKY SHORTWAVE TRANSMISSIVITY
184!  ZCEMTR(KPROMA,2)              ; CLEAR-SKY NET LONGWAVE EMISSIVITY
185!  ZTRSOD(KPROMA)                ; TOTAL-SKY SURFACE SW TRANSMISSITY
186!  ZLWFC (KPROMA,2)              ; CLEAR-SKY LONGWAVE FLUXES
187!  ZLWFT (KPROMA,KLEV+1)         ; TOTAL-SKY LONGWAVE FLUXES
188!  ZLWFT0(KPROMA,KLEV+1)         ; CLEAR-SKY LONGWAVE FLUXES      ! added by MPL 090109
189!  ZSWFC (KPROMA,2)              ; CLEAR-SKY SHORTWAVE FLUXES
190!  ZSWFT (KPROMA,KLEV+1)         ; TOTAL-SKY SHORTWAVE FLUXES
191!  ZSWFT0(KPROMA,KLEV+1)         ; CLEAR-SKY SHORTWAVE FLUXES     ! added by MPL 090109
192!  ZFLUX (KLON,2,KLEV+1)         ; TOTAL LW FLUXES  1=up, 2=DWN   ! added by MPL 080411
193!  ZFLUC (KLON,2,KLEV+1)         ; CLEAR SKY LW FLUXES            ! added by MPL 080411
194!  ZFSDWN(klon,KLEV+1)           ; TOTAL SW  DWN FLUXES           ! added by MPL 080411
195!  ZFCDWN(klon,KLEV+1)           ; CLEAR SKY SW  DWN FLUXES       ! added by MPL 080411
196!  ZFSUP (klon,KLEV+1)           ; TOTAL SW  UP  FLUXES           ! added by MPL 080411
197!  ZFCUP (klon,KLEV+1)           ; CLEAR SKY SW  UP  FLUXES       ! added by MPL 080411
198!MPL output RECMWFL:
199  REAL*8  ZEMTD (klon,klev+1),ZEMTD_i (klon,klev+1)       
200  REAL*8  ZEMTU (klon,klev+1),ZEMTU_i (klon,klev+1)     
201  REAL*8  ZTRSO (klon,klev+1),ZTRSO_i (klon,klev+1)   
202  REAL*8  ZTH   (klon,klev+1),ZTH_i   (klon,klev+1)   
203  REAL*8  ZCTRSO(klon,2)       
204  REAL*8  ZCEMTR(klon,2)     
205  REAL*8  ZTRSOD(klon)       
206  REAL*8  ZLWFC (klon,2)     
207  REAL*8  ZLWFT (klon,klev+1),ZLWFT_i (klon,klev+1)   
208  REAL*8  ZSWFC (klon,2)     
209  REAL*8  ZSWFT (klon,klev+1),ZSWFT_i (klon,klev+1)
210  REAL*8  ZSWFT0(klon,klev+1),ZLWFT0 (klon,klev+1)
211  REAL*8  PPIZA_DST(klon,klev,NSW)
212  REAL*8  PCGA_DST(klon,klev,NSW)
213  REAL*8  PTAUREL_DST(klon,klev,NSW)
214  REAL*8  PSFSWDIR(klon,NSW)
215  REAL*8  PSFSWDIF(klon,NSW)
216  REAL*8  PFSDNN(klon)
217  REAL*8  PFSDNV(klon)
218!MPL On ne redefinit pas les tableaux ZFLUX,ZFLUC,
219!MPL ZFSDWN,ZFCDWN,ZFSUP,ZFCUP car ils existent deja
220!MPL sous les noms de ZFLDN,ZFLDN0,ZFLUP,ZFLUP0,
221!MPL ZFSDN,ZFSDN0,ZFSUP,ZFSUP0
222  REAL*8  ZFLUX_i (klon,2,klev+1)
223  REAL*8  ZFLUC_i (klon,2,klev+1)
224  REAL*8  ZFSDWN_i (klon,klev+1)
225  REAL*8  ZFCDWN_i (klon,klev+1)
226  REAL*8  ZFSUP_i (klon,klev+1)
227  REAL*8  ZFCUP_i (klon,klev+1)
228
229! Local variables
230  REAL*8 ZFSUP(KDLON,KFLEV+1)
231  REAL*8 ZFSDN(KDLON,KFLEV+1)
232  REAL*8 ZFSUP0(KDLON,KFLEV+1)
233  REAL*8 ZFSDN0(KDLON,KFLEV+1)
234  REAL*8 ZFLUP(KDLON,KFLEV+1)
235  REAL*8 ZFLDN(KDLON,KFLEV+1)
236  REAL*8 ZFLUP0(KDLON,KFLEV+1)
237  REAL*8 ZFLDN0(KDLON,KFLEV+1)
238  REAL*8 zx_alpha1, zx_alpha2
239  INTEGER k, kk, i, j, iof, nb_gr
240  INTEGER ist,iend,ktdia,kmode
241  INTEGER , SAVE :: iprint=0
242  REAL*8 PSCT
243  REAL*8 PALBD(kdlon,2), PALBP(kdlon,2)
244  REAL*8 PALBD_NEW(kdlon,NSW), PALBP_NEW(kdlon,NSW)
245  REAL*8 PEMIS(kdlon), PDT0(kdlon), PVIEW(kdlon)
246  REAL*8 PPSOL(kdlon), PDP(kdlon,KLEV)
247  REAL*8 PTL(kdlon,kflev+1), PPMB(kdlon,kflev+1)
248  REAL*8 PTAVE(kdlon,kflev)
249  REAL*8 PWV(kdlon,kflev), PQS(kdlon,kflev), POZON(kdlon,kflev)
250!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6
251  REAL*8 PAER(kdlon,kflev,6)
252  REAL*8 PCLDLD(kdlon,kflev)
253  REAL*8 PCLDLU(kdlon,kflev)
254  REAL*8 PCLDSW(kdlon,kflev)
255  REAL*8 PTAU(kdlon,2,kflev)
256  REAL*8 POMEGA(kdlon,2,kflev)
257  REAL*8 PCG(kdlon,2,kflev)
258  REAL*8 zfract(kdlon), zrmu0(kdlon), zdist
259  REAL*8 zheat(kdlon,kflev), zcool(kdlon,kflev)
260  REAL*8 zheat0(kdlon,kflev), zcool0(kdlon,kflev)
261  REAL*8 ztopsw(kdlon), ztoplw(kdlon)
262  REAL*8 zsolsw(kdlon), zsollw(kdlon), zalbpla(kdlon)
263  REAL*8 zsollwdown(kdlon)
264  REAL*8 ztopsw0(kdlon), ztoplw0(kdlon)
265  REAL*8 zsolsw0(kdlon), zsollw0(kdlon)
266  REAL*8 zznormcp
267  REAL*8 tauaero(kdlon,kflev,9,2)                     ! aer opt properties
268  REAL*8 pizaero(kdlon,kflev,9,2)
269  REAL*8 cgaero(kdlon,kflev,9,2)
270  REAL*8 PTAUA(kdlon,2,kflev)                         ! present-day value of cloud opt thickness (PTAU is pre-industrial value), local use
271  REAL*8 POMEGAA(kdlon,2,kflev)                       ! dito for single scatt albedo
272  REAL*8 ztopswadaero(kdlon), zsolswadaero(kdlon)     ! Aerosol direct forcing at TOAand surface
273  REAL*8 ztopswad0aero(kdlon), zsolswad0aero(kdlon)   ! Aerosol direct forcing at TOAand surface
274  REAL*8 ztopswaiaero(kdlon), zsolswaiaero(kdlon)     ! dito, indirect
275  REAL*8 ztopsw_aero(kdlon,9), ztopsw0_aero(kdlon,9)
276  REAL*8 zsolsw_aero(kdlon,9), zsolsw0_aero(kdlon,9)
277
278  CHARACTER (LEN=20) :: modname
279  CHARACTER (LEN=80) :: abort_message
280
281  ! initialisation
282  ist=1
283  iend=klon
284  ktdia=1
285  kmode=ist 
286  lldebug=.FALSE.
287  ! initialisation
288  tauaero(:,:,:,:)=0.
289  pizaero(:,:,:,:)=0.
290  cgaero(:,:,:,:)=0.
291 
292  !
293  !-------------------------------------------
294  nb_gr = KLON / kdlon
295  IF (nb_gr*kdlon .NE. KLON) THEN
296      PRINT*, "kdlon mauvais:", KLON, kdlon, nb_gr
297      CALL abort
298  ENDIF
299  IF (kflev .NE. KLEV) THEN
300      PRINT*, "kflev differe de KLEV, kflev, KLEV"
301      CALL abort
302  ENDIF
303  !-------------------------------------------
304!     print *,'Entree de radlwsw, iflag_rrtm tsol=',iflag_rrtm,tsol
305  IF (iprint>10) THEN
306  DO k = 1, KLEV
307    DO i = 1, KLON
308!        print *,'En entree de radlwsw: k tsol temp',k,tsol,t(1,k)
309      heat(i,k)=0.
310      cool(i,k)=0.
311      heat0(i,k)=0.
312      cool0(i,k)=0.
313    ENDDO
314  ENDDO
315  ENDIF
316  !
317  zdist = dist
318  !
319  PSCT = solaire/zdist/zdist
320  DO j = 1, nb_gr
321    iof = kdlon*(j-1)
322    DO i = 1, kdlon
323      zfract(i) = fract(iof+i)
324      zrmu0(i) = rmu0(iof+i)
325      PALBD(i,1) = alb1(iof+i)
326!     PALBD(i,2) = alb1(iof+i)
327      PALBD(i,2) = alb2(iof+i)
328!
329      PALBD_NEW(i,1) = alb1(iof+i)
330      DO kk=2,NSW
331        PALBD_NEW(i,kk) = alb2(iof+i)
332      ENDDO
333!     PALBD_NEW(i,2) = alb2(iof+i)
334!     PALBD_NEW(i,3) = alb2(iof+i)
335!     PALBD_NEW(i,4) = alb2(iof+i)
336!     PALBD_NEW(i,5) = alb2(iof+i)
337!     PALBD_NEW(i,6) = alb2(iof+i)
338!
339      PALBP(i,1) = alb1(iof+i)
340!     PALBP(i,2) = alb1(iof+i)
341      PALBP(i,2) = alb2(iof+i)
342!
343      PALBP_NEW(i,1) = alb1(iof+i)
344      DO kk=2,NSW
345        PALBP_NEW(i,kk) = alb2(iof+i)
346      ENDDO
347!     PALBP_NEW(i,2) = alb2(iof+i)
348!     PALBP_NEW(i,3) = alb2(iof+i)
349!     PALBP_NEW(i,4) = alb2(iof+i)
350!     PALBP_NEW(i,5) = alb2(iof+i)
351!     PALBP_NEW(i,6) = alb2(iof+i)
352      PEMIS(i) = 1.0
353      PVIEW(i) = 1.66
354      PPSOL(i) = paprs(iof+i,1)
355      zx_alpha1 = (paprs(iof+i,1)-pplay(iof+i,2))/(pplay(iof+i,1)-pplay(iof+i,2))
356      zx_alpha2 = 1.0 - zx_alpha1
357      PTL(i,1) = t(iof+i,1) * zx_alpha1 + t(iof+i,2) * zx_alpha2
358      PTL(i,KLEV+1) = t(iof+i,KLEV)
359      PDT0(i) = tsol(iof+i) - PTL(i,1)
360    ENDDO
361    DO k = 2, kflev
362      DO i = 1, kdlon
363        PTL(i,k) = (t(iof+i,k)+t(iof+i,k-1))*0.5
364      ENDDO
365    ENDDO
366    DO k = 1, kflev
367      DO i = 1, kdlon
368        PDP(i,k) = paprs(iof+i,k)-paprs(iof+i,k+1)
369        PTAVE(i,k) = t(iof+i,k)
370        PWV(i,k) = MAX (q(iof+i,k), 1.0e-12)
371        PQS(i,k) = PWV(i,k)
372        ! wo:    cm.atm (epaisseur en cm dans la situation standard)
373        ! POZON: kg/kg
374        POZON(i,k) = MAX(wo(iof+i,k),1.0e-12)*RG/46.6968 &
375           /(paprs(iof+i,k)-paprs(iof+i,k+1))&
376           *(paprs(iof+i,1)/101325.0)
377        PCLDLD(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
378        PCLDLU(i,k) = cldfra(iof+i,k)*cldemi(iof+i,k)
379        PCLDSW(i,k) = cldfra(iof+i,k)
380        PTAU(i,1,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! 1e-12 serait instable
381        PTAU(i,2,k) = MAX(cldtaupi(iof+i,k), 1.0e-05)! pour 32-bit machines
382        POMEGA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAU(i,1,k))
383        POMEGA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAU(i,2,k))
384        PCG(i,1,k) = 0.865
385        PCG(i,2,k) = 0.910
386        !-
387        ! Introduced for aerosol indirect forcings.
388        ! The following values use the cloud optical thickness calculated from
389        ! present-day aerosol concentrations whereas the quantities without the
390        ! "A" at the end are for pre-industial (natural-only) aerosol concentrations
391        !
392        PTAUA(i,1,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! 1e-12 serait instable
393        PTAUA(i,2,k) = MAX(cldtaupd(iof+i,k), 1.0e-05)! pour 32-bit machines
394        POMEGAA(i,1,k) = 0.9999 - 5.0e-04 * EXP(-0.5 * PTAUA(i,1,k))
395        POMEGAA(i,2,k) = 0.9988 - 2.5e-03 * EXP(-0.05 * PTAUA(i,2,k))
396      ENDDO
397    ENDDO
398    !
399    DO k = 1, kflev+1
400      DO i = 1, kdlon
401        PPMB(i,k) = paprs(iof+i,k)/100.0
402      ENDDO
403    ENDDO
404    !
405!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6 
406    DO kk = 1, 6
407      DO k = 1, kflev
408        DO i = 1, kdlon
409          PAER(i,k,kk) = 1.0E-15
410        ENDDO
411      ENDDO
412    ENDDO
413    DO k = 1, kflev
414      DO i = 1, kdlon
415        tauaero(i,k,:,1)=tau_aero(iof+i,k,:,1)
416        pizaero(i,k,:,1)=piz_aero(iof+i,k,:,1)
417        cgaero(i,k,:,1) =cg_aero(iof+i,k,:,1)
418        tauaero(i,k,:,2)=tau_aero(iof+i,k,:,2)
419        pizaero(i,k,:,2)=piz_aero(iof+i,k,:,2)
420        cgaero(i,k,:,2) =cg_aero(iof+i,k,:,2)
421      ENDDO
422    ENDDO
423!
424!======================================================================
425!===== si iflag_rrtm=0 ================================================
426!IM ctes ds clesphys.h   CALL LW(RCO2,RCH4,RN2O,RCFC11,RCFC12,
427!IM ctes ds clesphys.h   CALL SW(PSCT, RCO2, zrmu0, zfract,
428!     
429      IF (iflag_rrtm.eq.0) then
430
431!----- Mise a zero des tableaux output du rayonnement LW-AR4 ----------             
432      DO k = 1, kflev+1
433      DO i = 1, kdlon
434      ZFLUP(i,k)=0.
435      ZFLDN(i,k)=0.
436      ZFLUP0(i,k)=0.
437      ZFLDN0(i,k)=0.
438      ENDDO
439      ENDDO
440      DO k = 1, kflev
441      DO i = 1, kdlon
442      zcool(i,k)=0.
443      zcool0(i,k)=0.
444      ENDDO
445      ENDDO
446      DO i = 1, kdlon
447      ztoplw(i)=0.
448      zsollw(i)=0.
449      ztoplw0(i)=0.
450      zsollw0(i)=0.
451      zsollwdown(i)=0.
452      ENDDO
453!
454    CALL LW_LMDAR4(&
455       PPMB, PDP,&
456       PPSOL,PDT0,PEMIS,&
457       PTL, PTAVE, PWV, POZON, PAER,&
458       PCLDLD,PCLDLU,&
459       PVIEW,&
460       zcool, zcool0,&
461       ztoplw,zsollw,ztoplw0,zsollw0,&
462       zsollwdown,&
463       ZFLUP, ZFLDN, ZFLUP0,ZFLDN0)
464
465!----- Mise a zero des tableaux output du rayonnement SW-AR4 ----------             
466      DO k = 1, kflev+1
467      DO i = 1, kdlon
468      ZFSUP(i,k)=0.
469      ZFSDN(i,k)=0.
470      ZFSUP0(i,k)=0.
471      ZFSDN0(i,k)=0.
472      ENDDO
473      ENDDO
474      DO k = 1, kflev
475      DO i = 1, kdlon
476      zheat(i,k)=0.
477      zheat0(i,k)=0.
478      ENDDO
479      ENDDO
480      DO i = 1, kdlon
481      zalbpla(i)=0.
482      ztopsw(i)=0.
483      zsolsw(i)=0.
484      ztopsw0(i)=0.
485      zsolsw0(i)=0.
486      ztopswadaero(i)=0.
487      zsolswadaero(i)=0.
488      ztopswaiaero(i)=0.
489      zsolswaiaero(i)=0.
490      ENDDO
491    IF (.NOT. new_aod) THEN
492       ! use old version
493        CALL SW_LMDAR4(PSCT, zrmu0, zfract,&
494           PPMB, PDP, &
495           PPSOL, PALBD, PALBP,&
496           PTAVE, PWV, PQS, POZON, PAER,&
497           PCLDSW, PTAU, POMEGA, PCG,&
498           zheat, zheat0,&
499           zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,&
500           ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
501           tau_aero(:,:,5,:), piz_aero(:,:,5,:), cg_aero(:,:,5,:),&
502           PTAUA, POMEGAA,&
503           ztopswadaero,zsolswadaero,&
504           ztopswaiaero,zsolswaiaero,&
505           ok_ade, ok_aie)
506    ELSE
507
508        CALL SW_AERO(PSCT, zrmu0, zfract,&
509           PPMB, PDP,&
510           PPSOL, PALBD, PALBP,&
511           PTAVE, PWV, PQS, POZON, PAER,&
512           PCLDSW, PTAU, POMEGA, PCG,&
513           zheat, zheat0,&
514           zalbpla,ztopsw,zsolsw,ztopsw0,zsolsw0,&
515           ZFSUP,ZFSDN,ZFSUP0,ZFSDN0,&
516           tauaero, pizaero, cgaero, &
517           PTAUA, POMEGAA,&
518           ztopswadaero,zsolswadaero,&
519           ztopswad0aero,zsolswad0aero,&
520           ztopswaiaero,zsolswaiaero, &
521           ztopsw_aero,ztopsw0_aero,&
522           zsolsw_aero,zsolsw0_aero,&
523           ok_ade, ok_aie)
524   
525    ENDIF
526
527!===== si iflag_rrtm=1, on passe dans SW via RECMWFL ===============
528!----- Mise a zero des tableaux output de RECMWF -------------------             
529
530      else
531#ifdef CPP_RRTM
532      DO k = 1, kflev+1
533      DO i = 1, kdlon
534      ZEMTD_i(i,k)=0.
535      ZEMTU_i(i,k)=0.
536      ZTRSO_i(i,k)=0.
537      ZTH_i(i,k)=0.
538      ZLWFT_i(i,k)=0.
539      ZSWFT_i(i,k)=0.
540      ZFLUX_i(i,1,k)=0.
541      ZFLUX_i(i,2,k)=0.
542      ZFLUC_i(i,1,k)=0.
543      ZFLUC_i(i,2,k)=0.
544      ZFSDWN_i(i,k)=0.
545      ZFCDWN_i(i,k)=0.
546      ZFSUP_i(i,k)=0.
547      ZFCUP_i(i,k)=0.
548      ENDDO
549      ENDDO
550!     
551      DO k = 1, kflev
552      DO i = 1, kdlon
553      DO kk = 1, NSW
554      PPIZA_DST(i,k,kk)=0.
555      PCGA_DST(i,k,kk)=0.
556      PTAUREL_DST(i,k,kk)=0.
557      ENDDO
558      ENDDO
559      ENDDO
560!     
561      DO i = 1, kdlon
562      ZCTRSO(i,1)=0.
563      ZCTRSO(i,2)=0.
564      ZCEMTR(i,1)=0.
565      ZCEMTR(i,2)=0.
566      ZTRSOD(i)=0.
567      ZLWFC(i,1)=0.
568      ZLWFC(i,2)=0.
569      ZSWFC(i,1)=0.
570      ZSWFC(i,2)=0.
571      PFSDNN(i)=0.
572      PFSDNV(i)=0.
573      DO kk = 1, NSW
574      PSFSWDIR(i,kk)=0.
575      PSFSWDIF(i,kk)=0.
576      ENDDO
577      ENDDO
578!----- Fin des mises a zero des tableaux output de RECMWF -------------------             
579         GEMU(1:klon)=sin(rlatd(1:klon))
580! On met les donnees dans l'ordre des niveaux arpege
581         paprs_i(:,1)=paprs(:,klev+1)
582         DO k=1,klev
583            paprs_i(1:klon,k+1) =paprs(1:klon,klev+1-k)
584            pplay_i(1:klon,k)   =pplay(1:klon,klev+1-k)
585            cldfra_i(1:klon,k)  =cldfra(1:klon,klev+1-k)
586            PDP_i(1:klon,k)     =PDP(1:klon,klev+1-k)
587            t_i(1:klon,k)       =t(1:klon,klev+1-k)
588            q_i(1:klon,k)       =q(1:klon,klev+1-k)
589            qsat_i(1:klon,k)    =qsat(1:klon,klev+1-k)
590            flwc_i(1:klon,k)    =flwc(1:klon,klev+1-k)
591            fiwc_i(1:klon,k)    =fiwc(1:klon,klev+1-k)
592         ENDDO
593         DO k=1,kflev
594            POZON_i(1:klon,k)=POZON(1:klon,kflev+1-k)
595!           print *,'Juste avant RECMWFL: k tsol temp',k,tsol,t(1,k)
596!!!!! Modif MPL 6.01.09 avec RRTM, on passe de 5 a 6     
597            DO i=1,6
598            PAER_i(1:klon,k,i)=PAER(1:klon,kflev+1-k,i)
599            ENDDO
600         ENDDO
601!        print *,'RADLWSW: avant RECMWFL, RI0,rmu0 =',solaire,rmu0
602
603!  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
604! La version ARPEGE1D utilise differentes valeurs de la constante
605! solaire suivant le rayonnement utilise.
606! A controler ...
607! SOLAR FLUX AT THE TOP (/YOMPHY3/)
608! introduce season correction
609!--------------------------------------
610! RII0 = RIP0
611! IF(LRAYFM)
612! RII0 = RIP0M   ! =rip0m if Morcrette non-each time step call.
613! IF(LRAYFM15)
614! RII0 = RIP0M15 ! =rip0m if Morcrette non-each time step call.
615         RII0=solaire
616!  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
617! Ancien appel a RECMWF (celui du cy25)
618!        CALL RECMWF (ist , iend, klon , ktdia , klev   , kmode ,
619!    s   PALBD    , PALBP   , paprs_i , pplay_i , RCO2   , cldfra_i,
620!    s   POZON_i  , PAER_i  , PDP_i   , PEMIS   , GEMU   , rmu0,
621!    s    q_i     , qsat_i  , fiwc_i  , flwc_i  , zmasq  , t_i  ,tsol,
622!    s   ZEMTD_i  , ZEMTU_i , ZTRSO_i ,
623!    s   ZTH_i    , ZCTRSO  , ZCEMTR  , ZTRSOD  ,
624!    s   ZLWFC    , ZLWFT_i , ZSWFC   , ZSWFT_i ,
625!    s   ZFLUX_i  , ZFLUC_i , ZFSDWN_i, ZFSUP_i , ZFCDWN_i,ZFCUP_i)
626!    s   'RECMWF ')
627!
628      IF(lldebug) then
629        call writefield_phy('paprs_i',paprs_i,klev+1)
630        call writefield_phy('pplay_i',pplay_i,klev)
631        call writefield_phy('cldfra_i',cldfra_i,klev)
632        call writefield_phy('pozon_i',POZON_i,klev)
633        call writefield_phy('paer_i',PAER_i,klev)
634        call writefield_phy('pdp_i',PDP_i,klev)
635        call writefield_phy('q_i',q_i,klev)
636        call writefield_phy('qsat_i',qsat_i,klev)
637        call writefield_phy('fiwc_i',fiwc_i,klev)
638        call writefield_phy('flwc_i',flwc_i,klev)
639        call writefield_phy('t_i',t_i,klev)
640        call writefield_phy('palbd_new',PALBD_NEW,NSW)
641        call writefield_phy('palbp_new',PALBP_NEW,NSW)
642      ENDIF
643
644! Nouvel appel a RECMWF (celui du cy32t0)
645         CALL RECMWF (ist , iend, klon , ktdia  , klev   , kmode ,&
646         PALBD_NEW,PALBP_NEW, paprs_i , pplay_i , RCO2   , cldfra_i,&
647         POZON_i  , PAER_i  , PDP_i   , PEMIS   , rmu0   ,&
648          q_i     , qsat_i  , fiwc_i  , flwc_i  , zmasq  , t_i  ,tsol,&
649         ZEMTD_i  , ZEMTU_i , ZTRSO_i ,&
650         ZTH_i    , ZCTRSO  , ZCEMTR  , ZTRSOD  ,&
651         ZLWFC    , ZLWFT_i , ZSWFC   , ZSWFT_i ,&
652         PSFSWDIR , PSFSWDIF, PFSDNN  , PFSDNV  ,&
653         PPIZA_DST, PCGA_DST,PTAUREL_DST,ZFLUX_i  , ZFLUC_i ,&
654         ZFSDWN_i , ZFSUP_i , ZFCDWN_i, ZFCUP_i)
655           
656         print *,'RADLWSW: apres RECMWF'
657      IF(lldebug) THEN
658        call writefield_phy('zemtd_i',ZEMTD_i,klev+1)
659        call writefield_phy('zemtu_i',ZEMTU_i,klev+1)
660        call writefield_phy('ztrso_i',ZTRSO_i,klev+1)
661        call writefield_phy('zth_i',ZTH_i,klev+1)
662        call writefield_phy('zctrso',ZCTRSO,2)
663        call writefield_phy('zcemtr',ZCEMTR,2)
664        call writefield_phy('ztrsod',ZTRSOD,1)
665        call writefield_phy('zlwfc',ZLWFC,2)
666        call writefield_phy('zlwft_i',ZLWFT_i,klev+1)
667        call writefield_phy('zswfc',ZSWFC,2)
668        call writefield_phy('zswft_i',ZSWFT_i,klev+1)
669        call writefield_phy('psfswdir',PSFSWDIR,6)
670        call writefield_phy('psfswdif',PSFSWDIF,6)
671        call writefield_phy('pfsdnn',PFSDNN,1)
672        call writefield_phy('pfsdnv',PFSDNV,1)
673        call writefield_phy('ppiza_dst',PPIZA_DST,klev)
674        call writefield_phy('pcga_dst',PCGA_DST,klev)
675        call writefield_phy('ptaurel_dst',PTAUREL_DST,klev)
676        call writefield_phy('zflux_i',ZFLUX_i,klev+1)
677        call writefield_phy('zfluc_i',ZFLUC_i,klev+1)
678        call writefield_phy('zfsdwn_i',ZFSDWN_i,klev+1)
679        call writefield_phy('zfsup_i',ZFSUP_i,klev+1)
680        call writefield_phy('zfcdwn_i',ZFCDWN_i,klev+1)
681        call writefield_phy('zfcup_i',ZFCUP_i,klev+1)
682      ENDIF
683! --------- output RECMWFL
684!  ZEMTD        (KPROMA,KLEV+1)  ; TOTAL DOWNWARD LONGWAVE EMISSIVITY
685!  ZEMTU        (KPROMA,KLEV+1)  ; TOTAL UPWARD   LONGWAVE EMISSIVITY
686!  ZTRSO        (KPROMA,KLEV+1)  ; TOTAL SHORTWAVE TRANSMISSIVITY
687!  ZTH          (KPROMA,KLEV+1)  ; HALF LEVEL TEMPERATURE
688!  ZCTRSO       (KPROMA,2)       ; CLEAR-SKY SHORTWAVE TRANSMISSIVITY
689!  ZCEMTR       (KPROMA,2)       ; CLEAR-SKY NET LONGWAVE EMISSIVITY
690!  ZTRSOD       (KPROMA)         ; TOTAL-SKY SURFACE SW TRANSMISSITY
691!  ZLWFC        (KPROMA,2)       ; CLEAR-SKY LONGWAVE FLUXES
692!  ZLWFT        (KPROMA,KLEV+1)  ; TOTAL-SKY LONGWAVE FLUXES
693!  ZSWFC        (KPROMA,2)       ; CLEAR-SKY SHORTWAVE FLUXES
694!  ZSWFT        (KPROMA,KLEV+1)  ; TOTAL-SKY SHORTWAVE FLUXES
695!  PPIZA_DST    (KPROMA,KLEV,NSW); Single scattering albedo of dust
696!  PCGA_DST     (KPROMA,KLEV,NSW); Assymetry factor for dust
697!  PTAUREL_DST  (KPROMA,KLEV,NSW); Optical depth of dust relative to at 550nm
698!  PSFSWDIR     (KPROMA,NSW)     ;
699!  PSFSWDIF     (KPROMA,NSW)     ;
700!  PFSDNN       (KPROMA)         ;
701!  PFSDNV       (KPROMA)         ;
702! ---------
703! On retablit l'ordre des niveaux lmd pour les tableaux de sortie
704      DO k=0,klev
705         ZEMTD(1:klon,k+1)  = ZEMTD_i(1:klon,k+1)
706         ZEMTU(1:klon,k+1)  = ZEMTU_i(1:klon,k+1)
707         ZTRSO(1:klon,k+1)  = ZTRSO_i(1:klon,k+1)
708         ZTH(1:klon,k+1)    = ZTH_i(1:klon,k+1)
709!        ZLWFT(1:klon,k+1)  = ZLWFT_i(1:klon,klev+1-k)
710!        ZSWFT(1:klon,k+1)  = ZSWFT_i(1:klon,klev+1-k)
711         ZFLUP(1:klon,k+1)  = ZFLUX_i(1:klon,1,k+1)
712         ZFLDN(1:klon,k+1)  = ZFLUX_i(1:klon,2,k+1)
713         ZFLUP0(1:klon,k+1) = ZFLUC_i(1:klon,1,k+1)
714         ZFLDN0(1:klon,k+1) = ZFLUC_i(1:klon,2,k+1)
715         ZFSDN(1:klon,k+1)  = ZFSDWN_i(1:klon,k+1)
716         ZFSDN0(1:klon,k+1) = ZFCDWN_i(1:klon,k+1)
717         ZFSUP (1:klon,k+1) = ZFSUP_i(1:klon,k+1)
718         ZFSUP0(1:klon,k+1) = ZFCUP_i(1:klon,k+1)
719!   Nouveau calcul car visiblement ZSWFT et ZSWFC sont nuls dans RRTM cy32
720!   en sortie de radlsw.F90 - MPL 7.01.09
721         ZSWFT(1:klon,k+1)  = ZFSDWN_i(1:klon,k+1)-ZFSUP_i(1:klon,k+1)
722         ZSWFT0(1:klon,k+1) = ZFCDWN_i(1:klon,k+1)-ZFCUP_i(1:klon,k+1)
723!        WRITE(*,'("FSDN FSUP FCDN FCUP: ",4E12.5)') ZFSDWN_i(1:klon,k+1),&
724!        ZFSUP_i(1:klon,k+1),ZFCDWN_i(1:klon,k+1),ZFCUP_i(1:klon,k+1)
725         ZLWFT(1:klon,k+1) =-ZFLUX_i(1:klon,2,k+1)-ZFLUX_i(1:klon,1,k+1)
726         ZLWFT0(1:klon,k+1)=-ZFLUC_i(1:klon,2,k+1)-ZFLUC_i(1:klon,1,k+1)
727!        print *,'FLUX2 FLUX1 FLUC2 FLUC1',ZFLUX_i(1:klon,2,k+1),
728!    s ZFLUX_i(1:klon,1,k+1),ZFLUC_i(1:klon,2,k+1),ZFLUC_i(1:klon,1,k+1)
729      ENDDO 
730      print*,'OK1'
731! ---------
732! On renseigne les champs LMDz, pour avoir la meme chose qu'en sortie de
733! LW_LMDAR4 et SW_LMDAR4
734      DO i = 1, kdlon
735         zsolsw(i)    = ZSWFT(i,1)
736         zsolsw0(i)   = ZSWFT0(i,1)
737!        zsolsw0(i)   = ZFSDN0(i,1)     -ZFSUP0(i,1)
738         ztopsw(i)    = ZSWFT(i,klev+1)
739         ztopsw0(i)   = ZSWFT0(i,klev+1)
740!        ztopsw0(i)   = ZFSDN0(i,klev+1)-ZFSUP0(i,klev+1)
741!         
742!        zsollw(i)    = ZFLDN(i,1)      -ZFLUP(i,1)
743!        zsollw0(i)   = ZFLDN0(i,1)     -ZFLUP0(i,1)
744!        ztoplw(i)    = ZFLDN(i,klev+1) -ZFLUP(i,klev+1)
745!        ztoplw0(i)   = ZFLDN0(i,klev+1)-ZFLUP0(i,klev+1)
746         zsollw(i)    = ZLWFT(i,1)
747         zsollw0(i)   = ZLWFT0(i,1)
748         ztoplw(i)    = ZLWFT(i,klev+1)
749         ztoplw0(i)   = ZLWFT0(i,klev+1)
750!         
751         zalbpla(i)   = ZFSUP(i,klev+1)/ZFSDN(i,klev+1)
752         zsollwdown(i)= ZFLDN(i,1)
753      ENDDO
754      print*,'OK2'
755
756! extrait de SW_AR4
757!     DO k = 1, KFLEV
758!        kpl1 = k+1
759!        DO i = 1, KDLON
760!           PHEAT(i,k) = -(ZFSUP(i,kpl1)-ZFSUP(i,k)) -(ZFSDN(i,k)-ZFSDN(i,kpl1))
761!           PHEAT(i,k) = PHEAT(i,k) * RDAY*RG/RCPD / PDP(i,k)
762! ZLWFT(klon,k),ZSWFT
763
764      DO k=1,kflev
765         DO i=1,kdlon
766           zheat(i,k)=(ZSWFT(i,k+1)-ZSWFT(i,k))*RDAY*RG/RCPD/PDP(i,k)
767           zheat0(i,k)=(ZSWFT0(i,k+1)-ZSWFT0(i,k))*RDAY*RG/RCPD/PDP(i,k)
768           zcool(i,k)=(ZLWFT(i,k)-ZLWFT(i,k+1))*RDAY*RG/RCPD/PDP(i,k)
769           zcool0(i,k)=(ZLWFT0(i,k)-ZLWFT0(i,k+1))*RDAY*RG/RCPD/PDP(i,k)
770!          print *,'heat cool heat0 coOl0 '&
771!         ,zheat(i,k),zcool(i,k),zheat0(i,k),zcool0(i,k)
772         ENDDO
773      ENDDO
774#else
775
776            abort_message='You should compile with -rrtm if running with iflag_rrtm=1'
777            call abort_gcm(modname,abort_message,1)
778#endif
779      ENDIF   ! if(iflag_rrtm=0)
780         print*,'OK3'
781!======================================================================
782!  PSOLSW(i)   = ZFSDN(i,1) - ZFSUP(i,1)
783!  PSOLSW0(i)  = ZFSDN0(i,1) - ZFSUP0(i,1)
784!  PSOLSWAD(i) = ZFSDNAD(i,1) - ZFSUPAD(i,1)
785!  PSOLSWAI(i) = ZFSDNAI(i,1) - ZFSUPAI(i,1)
786!  PTOPSW(i)   = ZFSDN(i,KFLEV+1) - ZFSUP(i,KFLEV+1)
787!  PTOPSW0(i)  = ZFSDN0(i,KFLEV+1) - ZFSUP0(i,KFLEV+1)
788!  PTOPSWAD(i) = ZFSDNAD(i,KFLEV+1) - ZFSUPAD(i,KFLEV+1)
789!  PTOPSWAI(i) = ZFSDNAI(i,KFLEV+1) - ZFSUPAI(i,KFLEV+1)     
790!======================================================================
791    DO i = 1, kdlon
792      radsol(iof+i) = zsolsw(i) + zsollw(i)
793      topsw(iof+i) = ztopsw(i)
794      toplw(iof+i) = ztoplw(i)
795      solsw(iof+i) = zsolsw(i)
796      sollw(iof+i) = zsollw(i)
797      sollwdown(iof+i) = zsollwdown(i)
798      DO k = 1, kflev+1
799        lwdn0 ( iof+i,k)   = ZFLDN0 ( i,k)
800        lwdn  ( iof+i,k)   = ZFLDN  ( i,k)
801        lwup0 ( iof+i,k)   = ZFLUP0 ( i,k)
802        lwup  ( iof+i,k)   = ZFLUP  ( i,k)
803      ENDDO
804      topsw0(iof+i) = ztopsw0(i)
805      toplw0(iof+i) = ztoplw0(i)
806      solsw0(iof+i) = zsolsw0(i)
807      sollw0(iof+i) = zsollw0(i)
808      albpla(iof+i) = zalbpla(i)
809
810      DO k = 1, kflev+1
811        swdn0 ( iof+i,k)   = ZFSDN0 ( i,k)
812        swdn  ( iof+i,k)   = ZFSDN  ( i,k)
813        swup0 ( iof+i,k)   = ZFSUP0 ( i,k)
814        swup  ( iof+i,k)   = ZFSUP  ( i,k)
815      ENDDO
816    ENDDO
817         print*,'OK4'
818    !-transform the aerosol forcings, if they have
819    ! to be calculated
820    IF (ok_ade) THEN
821        DO i = 1, kdlon
822          topswad_aero(iof+i) = ztopswadaero(i)
823          topswad0_aero(iof+i) = ztopswad0aero(i)
824          solswad_aero(iof+i) = zsolswadaero(i)
825          solswad0_aero(iof+i) = zsolswad0aero(i)
826          topsw_aero(iof+i,:) = ztopsw_aero(iof+i,:)
827          topsw0_aero(iof+i,:) = ztopsw0_aero(iof+i,:)
828          solsw_aero(iof+i,:) = zsolsw_aero(iof+i,:)
829          solsw0_aero(iof+i,:) = zsolsw0_aero(iof+i,:)
830         
831        ENDDO
832    ELSE
833        DO i = 1, kdlon
834          topswad_aero(iof+i) = 0.0
835          solswad_aero(iof+i) = 0.0
836          topswad0_aero(iof+i) = 0.0
837          solswad0_aero(iof+i) = 0.0
838          topsw_aero(iof+i,:) = 0.
839          topsw0_aero(iof+i,:) =0.
840          solsw_aero(iof+i,:) = 0.
841          solsw0_aero(iof+i,:) = 0.
842        ENDDO
843    ENDIF
844    IF (ok_aie) THEN
845        DO i = 1, kdlon
846          topswai_aero(iof+i) = ztopswaiaero(i)
847          solswai_aero(iof+i) = zsolswaiaero(i)
848        ENDDO
849    ELSE
850        DO i = 1, kdlon
851          topswai_aero(iof+i) = 0.0
852          solswai_aero(iof+i) = 0.0
853        ENDDO
854    ENDIF
855         print*,'OK5'
856    DO k = 1, kflev
857      DO i = 1, kdlon
858        !        scale factor to take into account the difference between
859        !        dry air and watter vapour scpecifi! heat capacity
860        zznormcp=1.0+RVTMP2*PWV(i,k)
861        heat(iof+i,k) = zheat(i,k)/zznormcp
862        cool(iof+i,k) = zcool(i,k)/zznormcp
863        heat0(iof+i,k) = zheat0(i,k)/zznormcp
864        cool0(iof+i,k) = zcool0(i,k)/zznormcp
865      ENDDO
866    ENDDO
867         print*,'OK6'
868    !
869  ENDDO
870       print*,'OK7'
871
872
873ENDSUBROUTINE radlwsw_aero
874
875
Note: See TracBrowser for help on using the repository browser.