source: LMDZ6/branches/LMDZ-INCA-Dyn/libf/phymar/rrtm_rrtm_140gp.F90 @ 5428

Last change on this file since 5428 was 2089, checked in by Laurent Fairhead, 10 years ago

Inclusion de la physique de MAR


Integration of MAR physics

File size: 9.9 KB
Line 
1!***************************************************************************
2!                                                                          *
3!                RRTM :  RAPID RADIATIVE TRANSFER MODEL                    *
4!                                                                          *
5!             ATMOSPHERIC AND ENVIRONMENTAL RESEARCH, INC.                 *
6!                        840 MEMORIAL DRIVE                                *
7!                        CAMBRIDGE, MA 02139                               *
8!                                                                          *
9!                           ELI J. MLAWER                                  *
10!                         STEVEN J. TAUBMAN~                               *
11!                         SHEPARD A. CLOUGH                                *
12!                                                                          *
13!                        ~currently at GFDL                                *
14!                                                                          *
15!                       email:  mlawer@aer.com                             *
16!                                                                          *
17!        The authors wish to acknowledge the contributions of the          *
18!        following people:  Patrick D. Brown, Michael J. Iacono,           *
19!        Ronald E. Farren, Luke Chen, Robert Bergstrom.                    *
20!                                                                          *
21!***************************************************************************
22!     Reformatted for F90 by JJMorcrette, ECMWF, 980714                    *
23!                                                                          *
24!***************************************************************************
25! *** mji ***
26! *** This version of RRTM has been altered to interface with either
27!     the ECMWF numerical weather prediction model or the ECMWF column
28!     radiation model (ECRT) package.
29
30!     Revised, April, 1997;  Michael J. Iacono, AER, Inc.
31!          - initial implementation of RRTM in ECRT code
32!     Revised, June, 1999;  Michael J. Iacono and Eli J. Mlawer, AER, Inc.
33!          - to implement generalized maximum/random cloud overlap
34
35SUBROUTINE RRTM_RRTM_140GP &
36 &( KIDIA , KFDIA , KLON , KLEV &
37 &, PAER  , PAPH  , PAP &
38 &, PTS   , PTH   , PT &
39 &, ZEMIS , ZEMIW &
40 &, PQ    , PCCO2 , POZN &
41 &, PCLDF , PTAUCLD &
42 &, PEMIT , PFLUX , PFLUC, PTCLEAR &
43 &)
44
45! *** This program is the driver for RRTM, the AER rapid model. 
46!     For each atmosphere the user wishes to analyze, this routine
47!     a) calls ECRTATM to read in the atmospheric profile
48!     b) calls SETCOEF to calculate various quantities needed for
49!        the radiative transfer algorithm
50!     c) calls RTRN to do the radiative transfer calculation for
51!        clear or cloudy sky
52!     d) writes out the upward, downward, and net flux for each
53!        level and the heating rate for each layer
54
55
56#include "tsmbkind.h"
57
58USE PARRRTM  , ONLY : JPBAND   ,JPG      ,JPXSEC   ,JPGPT    ,JPLAY    ,&
59            &JPINPX
60USE YOERRTWN , ONLY : WAVENUM1 ,WAVENUM2 ,DELWAVE  ,NG       ,NSPA     ,NSPB
61
62!------------------------------Arguments--------------------------------
63
64! Input arguments
65
66
67IMPLICIT NONE
68INTEGER_M :: kidia                 ! First atmosphere index
69INTEGER_M :: kfdia                 ! Last atmosphere index
70INTEGER_M :: klon                  ! Number of atmospheres (longitudes)
71INTEGER_M :: klev                  ! Number of atmospheric layers
72REAL_B :: paer(klon,6,klev)        ! Aerosol optical thickness
73REAL_B :: pap(klon,klev)           ! Layer pressures (Pa)
74REAL_B :: paph(klon,klev+1)        ! Interface pressures (Pa)
75REAL_B :: pts(klon)                ! Surface temperature (K)
76REAL_B :: pt(klon,klev)            ! Layer temperature (K)
77REAL_B :: zemis(klon)              ! Non-window surface emissivity
78REAL_B :: zemiw(klon)              ! Window surface emissivity
79REAL_B :: pth(klon,klev+1)         ! Interface temperatures (K)
80REAL_B :: pq(klon,klev)            ! H2O specific humidity (mmr)
81REAL_B :: pozn(klon,klev)          ! O3 mass mixing ratio
82REAL_B :: pcco2                    ! CO2 mass mixing ratio
83REAL_B :: rch4                     ! CH4 mass mixing ratio
84REAL_B :: rn2o                     ! N2O mass mixing ratio
85REAL_B :: rcfc11                   ! CFC11 mass mixing ratio
86REAL_B :: rcfc12                   ! CFC12 mass mixing ratio
87REAL_B :: pcldf(klon,klev)         ! Cloud fraction
88REAL_B :: ptaucld(klon,klev,JPBAND)! Cloud optical depth
89REAL_B :: PFLUX(klon,2,klev+1)     ! LW total sky flux (1=up, 2=down)
90REAL_B :: PFLUC(klon,2,klev+1)     ! LW clear sky flux (1=up, 2=down)
91REAL_B :: PEMIT(klon)              ! Surface LW emissivity
92REAL_B :: PTCLEAR(klon)            ! clear-sky fraction of column
93
94INTEGER_M :: ICLDLYR(JPLAY)        ! Cloud indicator
95REAL_B :: CLDFRAC(JPLAY)           ! Cloud fraction
96REAL_B :: TAUCLD(JPLAY,JPBAND)     ! Spectral optical thickness
97
98REAL_B :: ABSS1 (JPGPT*JPLAY)
99REAL_B :: ATR1  (JPGPT,JPLAY)
100EQUIVALENCE (ABSS1(1),ATR1(1,1))
101
102REAL_B :: OD    (JPGPT,JPLAY)
103
104REAL_B :: TAUSF1(JPGPT*JPLAY)
105REAL_B :: TF1   (JPGPT,JPLAY)
106EQUIVALENCE (TAUSF1(1),TF1(1,1))
107
108REAL_B :: COLDRY(JPLAY)
109REAL_B :: WKL(JPINPX,JPLAY)
110
111REAL_B :: WX(JPXSEC,JPLAY)         ! Amount of trace gases
112
113REAL_B :: CLFNET  (0:JPLAY)
114REAL_B :: CLHTR   (0:JPLAY)
115REAL_B :: FNET    (0:JPLAY)
116REAL_B :: HTR     (0:JPLAY)
117REAL_B :: TOTDFLUC(0:JPLAY)
118REAL_B :: TOTDFLUX(0:JPLAY)
119REAL_B :: TOTUFLUC(0:JPLAY)
120REAL_B :: TOTUFLUX(0:JPLAY)
121
122!     LOCAL INTEGER SCALARS
123INTEGER_M :: i, icld, iplon, K
124INTEGER_M :: ISTART
125INTEGER_M :: IEND
126
127!     LOCAL REAL SCALARS
128REAL_B :: FLUXFAC, HEATFAC, PI, ZEPSEC, ZTCLEAR
129
130!- from AER
131REAL_B :: TAUAERL(JPLAY,JPBAND)
132
133!- from INTFAC     
134REAL_B :: FAC00(JPLAY)
135REAL_B :: FAC01(JPLAY)
136REAL_B :: FAC10(JPLAY)
137REAL_B :: FAC11(JPLAY)
138REAL_B :: FORFAC(JPLAY)
139
140!- from INTIND
141INTEGER_M :: JP(JPLAY)
142INTEGER_M :: JT(JPLAY)
143INTEGER_M :: JT1(JPLAY)
144
145!- from PRECISE             
146REAL_B :: ONEMINUS
147
148!- from PROFDATA             
149REAL_B :: COLH2O(JPLAY)
150REAL_B :: COLCO2(JPLAY)
151REAL_B :: COLO3 (JPLAY)
152REAL_B :: COLN2O(JPLAY)
153REAL_B :: COLCH4(JPLAY)
154REAL_B :: COLO2 (JPLAY)
155REAL_B :: CO2MULT(JPLAY)
156INTEGER_M :: LAYTROP
157INTEGER_M :: LAYSWTCH
158INTEGER_M :: LAYLOW
159
160!- from PROFILE             
161REAL_B :: PAVEL(JPLAY)
162REAL_B :: TAVEL(JPLAY)
163REAL_B :: PZ(0:JPLAY)
164REAL_B :: TZ(0:JPLAY)
165REAL_B :: TBOUND
166INTEGER_M :: NLAYERS
167
168!- from SELF             
169REAL_B :: SELFFAC(JPLAY)
170REAL_B :: SELFFRAC(JPLAY)
171INTEGER_M :: INDSELF(JPLAY)
172
173!- from SP             
174REAL_B :: PFRAC(JPGPT,JPLAY)
175
176!- from SURFACE             
177REAL_B :: SEMISS(JPBAND)
178REAL_B :: SEMISLW
179INTEGER_M :: IREFLECT
180
181
182!     HEATFAC is the factor by which one must multiply delta-flux/
183!     delta-pressure, with flux in w/m-2 and pressure in mbar, to get
184!     the heating rate in units of degrees/day.  It is equal to
185!           (g)x(#sec/day)x(1e-5)/(specific heat of air at const. p)
186!        =  (9.8066)(86400)(1e-5)/(1.004)
187
188ZEPSEC = 1.E-06_JPRB
189ONEMINUS = _ONE_ - ZEPSEC
190PI = _TWO_*ASIN(_ONE_)
191FLUXFAC = PI * 2.D4
192HEATFAC = 8.4391_JPRB
193
194! *** mji ***
195! For use with ECRT, this loop is over atmospheres (or longitudes)
196DO iplon = kidia,kfdia
197
198! *** mji ***
199!- Prepare atmospheric profile from ECRT for use in RRTM, and define
200!  other RRTM input parameters.  Arrays are passed back through the
201!  existing RRTM commons and arrays.
202  ZTCLEAR=_ONE_
203
204!  print *,'before RRTM_ECRT_140GP'
205
206  CALL RRTM_ECRT_140GP &
207   &( iplon, klon , klev, icld &
208   &, paer , paph , pap &
209   &, pts  , pth  , pt &
210   &, zemis, zemiw &
211   &, pq   , pcco2, pozn, pcldf, ptaucld, ztclear &
212   &, CLDFRAC,TAUCLD,COLDRY,WKL,WX &
213   &, TAUAERL,PAVEL,TAVEL,PZ,TZ,TBOUND,NLAYERS,SEMISS,IREFLECT)
214
215  PTCLEAR(iplon)=ztclear
216
217  ISTART = 1
218  IEND   = 16
219
220!  Calculate information needed by the radiative transfer routine
221!  that is specific to this atmosphere, especially some of the
222!  coefficients and indices needed to compute the optical depths
223!  by interpolating data from stored reference atmospheres.
224
225!  print *,'before RRTM_SETCOEF_140GP'
226 
227    CALL RRTM_SETCOEF_140GP (KLEV,COLDRY,WKL &
228   &, FAC00,FAC01,FAC10,FAC11,FORFAC,JP,JT,JT1 &
229   &, COLH2O,COLCO2,COLO3,COLN2O,COLCH4,COLO2,CO2MULT &
230   &, LAYTROP,LAYSWTCH,LAYLOW,PAVEL,TAVEL,SELFFAC,SELFFRAC,INDSELF)
231
232! print *,'before RRTM_GASABS1A_140GP'
233
234  CALL RRTM_GASABS1A_140GP (KLEV,ATR1,OD,TF1,COLDRY,WX &
235   &, TAUAERL,FAC00,FAC01,FAC10,FAC11,FORFAC,JP,JT,JT1,ONEMINUS &
236   &, COLH2O,COLCO2,COLO3,COLN2O,COLCH4,COLO2,CO2MULT &
237   &, LAYTROP,LAYSWTCH,LAYLOW,SELFFAC,SELFFRAC,INDSELF,PFRAC)
238
239!- Call the radiative transfer routine.
240
241! *** mji ***
242!  Check for cloud in column.  Use ECRT threshold set as flag icld in
243!  routine ECRTATM.  If icld=1 then column is cloudy, otherwise it is
244!  clear.  Also, set up flag array, icldlyr, for use in radiative
245!  transfer.  Set icldlyr to one for each layer with non-zero cloud
246!  fraction.
247
248  DO K = 1, KLEV
249    IF (ICLD == 1.AND.CLDFRAC(K) > ZEPSEC) THEN
250      ICLDLYR(K) = 1
251    ELSE
252      ICLDLYR(K) = 0
253    ENDIF
254  ENDDO
255
256!  Clear and cloudy parts of column are treated together in RTRN.
257!  Clear radiative transfer is done for clear layers and cloudy radiative
258!  transfer is done for cloudy layers as identified by icldlyr.
259 
260! print *,'before RRTM_RTRN1A_140GP'
261 
262  CALL RRTM_RTRN1A_140GP (KLEV,ISTART,IEND,ICLDLYR,CLDFRAC,TAUCLD,ABSS1 &
263   &, OD,TAUSF1,CLFNET,CLHTR,FNET,HTR,TOTDFLUC,TOTDFLUX,TOTUFLUC,TOTUFLUX &
264   &, TAVEL,PZ,TZ,TBOUND,PFRAC,SEMISS,SEMISLW,IREFLECT)
265
266! ***   Pass clear sky and total sky up and down flux profiles to ECRT
267!       output arrays (zflux, zfluc). Array indexing from bottom to top
268!       is preserved for ECRT.
269!       Invert down flux arrays for consistency with ECRT sign conventions.
270
271  pemit(iplon) = SEMISLW
272  DO i = 0, KLEV
273    pfluc(iplon,1,i+1) =  TOTUFLUC(i)*FLUXFAC
274    pfluc(iplon,2,i+1) = -TOTDFLUC(i)*FLUXFAC
275    pflux(iplon,1,i+1) =  TOTUFLUX(i)*FLUXFAC
276    pflux(iplon,2,i+1) = -TOTDFLUX(i)*FLUXFAC
277  ENDDO
278ENDDO
279
280RETURN
281END SUBROUTINE RRTM_RRTM_140GP
Note: See TracBrowser for help on using the repository browser.