source: LMDZ6/branches/LMDZ_ECRad/libf/phylmd/ecrad/radiation_scheme.F90 @ 4388

Last change on this file since 4388 was 4388, checked in by idelkadi, 20 months ago

Code cleanup:

  • Removal of unnecessary arguments in the LMDZ-ECRAD interface routine.
  • deletion of useless comments
File size: 27.9 KB
Line 
1! AI mars 2021
2! ====================== Interface between ECRAD and LMDZ ====================
3! radiation_scheme.F90 appelee dans radlwsw_m.F90 si iflag_rttm = 2
4! revoir toutes les parties avec "AI ATTENTION"
5! Mars 2021 :
6!             - Revoir toutes les parties commentees AI ATTENTION
7!             1. Traitement des aerosols
8!             2. Verifier les parametres times issus de LMDZ (calcul issed)
9!             3. Configuration a partir de namelist
10!             4. frac_std = 0.75     
11! ============================================================================
12
13SUBROUTINE RADIATION_SCHEME &
14! Inputs
15     & (KIDIA, KFDIA, KLON, KLEV, KAEROSOL, NSW, &
16     &  IDAY, TIME, &
17     &  PSOLAR_IRRADIANCE, &
18     &  PMU0, PTEMPERATURE_SKIN, &
19     &  PALBEDO_DIF, PALBEDO_DIR, &
20     &  PEMIS, PEMIS_WINDOW, &
21     &  PGELAM, PGEMU, &
22     &  PPRESSURE_H, PTEMPERATURE_H, PQ, PQSAT, &
23     &  PCO2, PCH4, PN2O, PNO2, PCFC11, PCFC12, PHCFC22, &
24     &  PCCL4, PO3, PO2, &
25     &  PCLOUD_FRAC, PQ_LIQUID, PQ_ICE, PQ_SNOW, &
26     &  ZRE_LIQUID_UM, ZRE_ICE_UM, &
27     &  PAEROSOL_OLD, PAEROSOL, &
28! Outputs
29     &  PFLUX_SW, PFLUX_LW, PFLUX_SW_CLEAR, PFLUX_LW_CLEAR, &
30     &  PFLUX_SW_DN, PFLUX_LW_DN, PFLUX_SW_DN_CLEAR, PFLUX_LW_DN_CLEAR, &
31     &  PFLUX_SW_UP, PFLUX_LW_UP, PFLUX_SW_UP_CLEAR, PFLUX_LW_UP_CLEAR, &
32     &  PFLUX_DIR, PFLUX_DIR_CLEAR, PFLUX_DIR_INTO_SUN, &
33     &  PFLUX_UV, PFLUX_PAR, PFLUX_PAR_CLEAR, &
34     &  PEMIS_OUT, PLWDERIVATIVE, &
35     &  PSWDIFFUSEBAND, PSWDIRECTBAND)
36
37
38! RADIATION_SCHEME - Interface to modular radiation scheme
39!
40! (C) Copyright 2015- ECMWF.
41!
42! This software is licensed under the terms of the Apache Licence Version 2.0
43! which can be obtained at http://www.apache.org/licenses/LICENSE-2.0.
44!
45! In applying this licence, ECMWF does not waive the privileges and immunities
46! granted to it by virtue of its status as an intergovernmental organisation
47! nor does it submit to any jurisdiction.
48!
49! PURPOSE
50! -------
51!   The modular radiation scheme is contained in a separate
52!   library. This routine puts the the IFS arrays into appropriate
53!   objects, computing the additional data that is required, and sends
54!   it to the radiation scheme.  It returns net fluxes and surface
55!   flux components needed by the rest of the model.
56!
57!   Lower case is used for variables and types taken from the
58!   radiation library
59!
60! INTERFACE
61! ---------
62!    RADIATION_SCHEME is called from RADLSWR. The
63!    SETUP_RADIATION_SCHEME routine (in the RADIATION_SETUP module)
64!    should have been run first.
65!
66! AUTHOR
67! ------
68!   Robin Hogan, ECMWF
69!   Original: 2015-09-16
70!
71! MODIFICATIONS
72! -------------
73!
74! TO DO
75! -----
76!
77!-----------------------------------------------------------------------
78
79! Modules from ifs or ifsaux libraries
80USE PARKIND1 , ONLY : JPIM, JPRB
81USE YOMHOOK  , ONLY : LHOOK, DR_HOOK
82! AI ATTENTION module propre a ifs commentes
83!USE YOERAD   , ONLY : YRERAD
84USE RADIATION_SETUP, ONLY : SETUP_RADIATION_SCHEME, rad_config, &
85!USE RADIATION_SETUP, ONLY : &
86     &  NWEIGHT_UV,  IBAND_UV,  WEIGHT_UV, &
87     &  NWEIGHT_PAR, IBAND_PAR, WEIGHT_PAR, &
88     &  ITYPE_TROP_BG_AER,  TROP_BG_AER_MASS_EXT, &
89     &  ITYPE_STRAT_BG_AER, STRAT_BG_AER_MASS_EXT
90! Commentes : jour, date de la simulation
91!USE YOMRIP0  , ONLY : NINDAT
92!USE YOMCT3   , ONLY : NSTEP
93!USE YOMRIP   , ONLY : YRRIP
94USE YOMCST   , ONLY : RSIGMA ! Stefan-Boltzmann constant
95
96! Modules from radiation library
97USE radiation_single_level,   ONLY : single_level_type
98USE radiation_thermodynamics, ONLY : thermodynamics_type
99USE radiation_gas
100USE radiation_cloud,          ONLY : cloud_type
101USE radiation_aerosol,        ONLY : aerosol_type
102USE radiation_flux,           ONLY : flux_type
103USE radiation_interface,      ONLY : radiation, set_gas_units
104USE radiation_save,           ONLY : save_inputs
105
106USE mod_phys_lmdz_para
107
108IMPLICIT NONE
109
110! INPUT ARGUMENTS
111! *** Array dimensions and ranges
112INTEGER(KIND=JPIM),INTENT(IN) :: KIDIA    ! Start column to process
113INTEGER(KIND=JPIM),INTENT(IN) :: KFDIA    ! End column to process
114!INTEGER, INTENT(IN) :: KIDIA, KFDIA
115INTEGER(KIND=JPIM),INTENT(IN) :: KLON     ! Number of columns
116INTEGER(KIND=JPIM),INTENT(IN) :: KLEV     ! Number of levels
117!INTEGER, INTENT(IN) :: KLON, KLEV
118!INTEGER(KIND=JPIM),INTENT(IN) :: KAEROLMDZ ! Number of aerosol types
119INTEGER(KIND=JPIM),INTENT(IN) :: KAEROSOL
120INTEGER(KIND=JPIM),INTENT(IN) :: NSW ! Numbe of bands
121
122! AI ATTENTION
123!INTEGER, PARAMETER :: KAEROSOL = 12
124
125! *** Single-level fields
126REAL(KIND=JPRB),   INTENT(IN) :: PSOLAR_IRRADIANCE ! (W m-2)
127REAL(KIND=JPRB),   INTENT(IN) :: PMU0(KLON) ! Cosine of solar zenith ang
128REAL(KIND=JPRB),   INTENT(IN) :: PTEMPERATURE_SKIN(KLON) ! (K)
129! Diffuse and direct components of surface shortwave albedo
130!REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIF(KLON,YRERAD%NSW)
131!REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIR(KLON,YRERAD%NSW)
132REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIF(KLON,NSW)
133REAL(KIND=JPRB),   INTENT(IN) :: PALBEDO_DIR(KLON,NSW)
134! Longwave emissivity outside and inside the window region
135REAL(KIND=JPRB),   INTENT(IN) :: PEMIS(KLON)
136REAL(KIND=JPRB),   INTENT(IN) :: PEMIS_WINDOW(KLON)
137! Longitude (radians), sine of latitude
138REAL(KIND=JPRB),   INTENT(IN) :: PGELAM(KLON)
139REAL(KIND=JPRB),   INTENT(IN) :: PGEMU(KLON)
140! Land-sea mask
141!REAL(KIND=JPRB),   INTENT(IN) :: PLAND_SEA_MASK(KLON)
142
143! *** Variables on full levels
144!REAL(KIND=JPRB),   INTENT(IN) :: PPRESSURE(KLON,KLEV)    ! (Pa)
145!REAL(KIND=JPRB),   INTENT(IN) :: PTEMPERATURE(KLON,KLEV) ! (K)
146! *** Variables on half levels
147REAL(KIND=JPRB),   INTENT(IN) :: PPRESSURE_H(KLON,KLEV+1)    ! (Pa)
148REAL(KIND=JPRB),   INTENT(IN) :: PTEMPERATURE_H(KLON,KLEV+1) ! (K)
149
150! *** Gas mass mixing ratios on full levels
151REAL(KIND=JPRB),   INTENT(IN) :: PQ(KLON,KLEV)
152! AI
153REAL(KIND=JPRB),   INTENT(IN) :: PQSAT(KLON,KLEV)
154REAL(KIND=JPRB),   INTENT(IN) :: PCO2
155REAL(KIND=JPRB),   INTENT(IN) :: PCH4
156REAL(KIND=JPRB),   INTENT(IN) :: PN2O
157REAL(KIND=JPRB),   INTENT(IN) :: PNO2
158REAL(KIND=JPRB),   INTENT(IN) :: PCFC11
159REAL(KIND=JPRB),   INTENT(IN) :: PCFC12
160REAL(KIND=JPRB),   INTENT(IN) :: PHCFC22
161REAL(KIND=JPRB),   INTENT(IN) :: PCCL4
162REAL(KIND=JPRB),   INTENT(IN) :: PO3(KLON,KLEV) ! AI (kg/kg) ATTENTION (Pa*kg/kg)
163REAL(KIND=JPRB),   INTENT(IN) :: PO2
164
165! *** Cloud fraction and hydrometeor mass mixing ratios
166REAL(KIND=JPRB),   INTENT(IN) :: PCLOUD_FRAC(KLON,KLEV)
167REAL(KIND=JPRB),   INTENT(IN) :: PQ_LIQUID(KLON,KLEV)
168REAL(KIND=JPRB),   INTENT(IN) :: PQ_ICE(KLON,KLEV)
169!REAL(KIND=JPRB),   INTENT(IN) :: PQ_RAIN(KLON,KLEV)
170REAL(KIND=JPRB),   INTENT(IN) :: PQ_SNOW(KLON,KLEV)
171
172! *** Aerosol mass mixing ratios
173REAL(KIND=JPRB),   INTENT(IN) :: PAEROSOL_OLD(KLON,6,KLEV)
174REAL(KIND=JPRB),   INTENT(IN) :: PAEROSOL(KLON,KLEV,KAEROSOL)
175
176!REAL(KIND=JPRB),   INTENT(IN) :: PCCN_LAND(KLON)
177!REAL(KIND=JPRB),   INTENT(IN) :: PCCN_SEA(KLON)
178
179!AI mars 2021
180INTEGER(KIND=JPIM), INTENT(IN) :: IDAY
181REAL(KIND=JPRB), INTENT(IN)    :: TIME
182
183
184! OUTPUT ARGUMENTS
185
186! *** Net fluxes on half-levels (W m-2)
187REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW(KLON,KLEV+1)
188REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW(KLON,KLEV+1)
189REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_CLEAR(KLON,KLEV+1)
190REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_CLEAR(KLON,KLEV+1)
191
192!*** DN and UP flux on half-levels (W m-2)
193REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN(KLON,KLEV+1)
194REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN(KLON,KLEV+1)
195REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN_CLEAR(KLON,KLEV+1)
196REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN_CLEAR(KLON,KLEV+1)
197REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP(KLON,KLEV+1)
198REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP(KLON,KLEV+1)
199REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP_CLEAR(KLON,KLEV+1)
200REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP_CLEAR(KLON,KLEV+1)
201
202! *** Surface flux components (W m-2)
203! AI ATTENTION
204!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN_SURF(KLON)
205!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN_SURF(KLON)
206!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN_CLEAR_SURF(KLON)
207!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN_CLEAR_SURF(KLON)
208!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP_SURF(KLON)
209!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP_SURF(KLON)
210!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP_CLEAR_SURF(KLON)
211!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP_CLEAR_SURF(KLON)
212
213! Direct component of surface flux into horizontal plane
214REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR(KLON)
215REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_CLEAR(KLON)
216! As PFLUX_DIR but into a plane perpendicular to the sun
217REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_DIR_INTO_SUN(KLON)
218
219! *** Ultraviolet and photosynthetically active radiation (W m-2)
220REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_UV(KLON)
221REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_PAR(KLON)
222REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_PAR_CLEAR(KLON)
223
224! *** Other single-level diagnostics
225! Top-of-atmosphere incident solar flux (W m-2)
226! AI ATTENTION
227!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_DN_TOA(KLON)
228!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_SW_UP_TOA(KLON)
229!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_DN_TOA(KLON)
230!REAL(KIND=JPRB),  INTENT(OUT) :: PFLUX_LW_UP_TOA(KLON)
231
232! Diagnosed longwave surface emissivity across the whole spectrum
233REAL(KIND=JPRB),  INTENT(OUT) :: PEMIS_OUT(KLON)   
234
235! Partial derivative of total-sky longwave upward flux at each level
236! with respect to upward flux at surface, used to correct heating
237! rates at gridpoints/timesteps between calls to the full radiation
238! scheme.  Note that this version uses the convention of level index
239! increasing downwards, unlike the local variable ZLwDerivative that
240! is returned from the LW radiation scheme.
241REAL(KIND=JPRB),  INTENT(OUT) :: PLWDERIVATIVE(KLON,KLEV+1)
242
243! Surface diffuse and direct downwelling shortwave flux in each
244! shortwave albedo band, used in RADINTG to update the surface fluxes
245! accounting for high-resolution albedo information
246REAL(KIND=JPRB),  INTENT(OUT) :: PSWDIFFUSEBAND(KLON,NSW)
247REAL(KIND=JPRB),  INTENT(OUT) :: PSWDIRECTBAND (KLON,NSW)
248
249! LOCAL VARIABLES
250! AI ATTENTION
251!type(config_type)         :: rad_config
252TYPE(single_level_type)   :: single_level
253TYPE(thermodynamics_type) :: thermodynamics
254TYPE(gas_type)            :: gas
255TYPE(cloud_type)          :: cloud
256TYPE(aerosol_type)        :: aerosol
257TYPE(flux_type)           :: flux
258
259! Mass mixing ratio of ozone (kg/kg)
260REAL(KIND=JPRB)           :: ZO3(KLON,KLEV)
261
262! Cloud effective radii in microns
263REAL(KIND=JPRB)           :: ZRE_LIQUID_UM(KLON,KLEV)
264REAL(KIND=JPRB)           :: ZRE_ICE_UM(KLON,KLEV)
265
266! Cloud overlap decorrelation length for cloud boundaries in km
267REAL(KIND=JPRB)           :: ZDECORR_LEN_KM(KLON)
268
269! Ratio of cloud overlap decorrelation length for cloud water
270! inhomogeneities to that for cloud boundaries (typically 0.5)
271REAL(KIND=JPRB)           :: ZDECORR_LEN_RATIO
272
273! The surface net longwave flux if the surface was a black body, used
274! to compute the effective broadband surface emissivity
275REAL(KIND=JPRB)           :: ZBLACK_BODY_NET_LW(KIDIA:KFDIA)
276
277! Layer mass in kg m-2
278REAL(KIND=JPRB)           :: ZLAYER_MASS(KIDIA:KFDIA,KLEV)
279
280! Time integers
281INTEGER :: ITIM
282
283! Loop indices
284INTEGER :: JLON, JLEV, JBAND, JB_ALBEDO, JAER
285
286REAL(KIND=JPRB) :: ZHOOK_HANDLE
287
288! AI ATTENTION traitement aerosols
289INTEGER, PARAMETER :: NAERMACC = 1
290
291! AI ATTENTION
292real(jprb), parameter    :: frac_std = 0.75
293
294! Name of file names specified on command line
295character(len=512) :: file_name
296
297logical :: loutput=.true.
298logical :: lprint_input=.false.
299logical :: lprint_config=.false.
300logical, save :: debut_ecrad=.true.
301!$OMP THREADPRIVATE(debut_ecrad)
302
303! Import time functions for iseed calculation
304! AI ATTENTION propre a ifs
305!#include "fcttim.func.h"
306!#include "liquid_effective_radius.intfb.h"
307!#include "ice_effective_radius.intfb.h"
308!#include "cloud_overlap_decorr_len.intfb.h"
309!#include "satur.intfb.h"
310
311! Verifier les inputs
312print*,'=============== dans radiation_scheme : ==================='
313if (lprint_input) then
314  print*,'********** Verification des entrees *************'
315  print*,'KIDIA, KFDIA, KLON, KLEV, KAEROSOL, NSW =', &
316        KIDIA, KFDIA, KLON, KLEV, KAEROSOL, NSW
317  print*,'IDAY, TIME =', IDAY, TIME
318  print*,'PSOLAR_IRRADIANCE =', PSOLAR_IRRADIANCE
319  print*,'PMU0 =', PMU0
320  print*,'PTEMPERATURE_SKIN =',PTEMPERATURE_SKIN
321  print*,'PEMIS, PEMIS_WINDOW =', PEMIS, PEMIS_WINDOW
322  print*,'PGELAM, PGEMU =', PGELAM, PGEMU
323  print*,'PPRESSURE_H =', PPRESSURE_H
324  print*,'PTEMPERATURE_H =', PTEMPERATURE_H
325  print*,'PQ =', PQ
326  print*,'PQSAT=',PQSAT
327  print*,'PCO2, PCH4, PN2O, PNO2, PCFC11, PCFC12, PHCFC22, PCCL4 =', &
328        PCO2, PCH4, PN2O, PNO2, PCFC11, PCFC12, PHCFC22, PCCL4
329  print*,'PO3 =',PO3
330  print*,'PCLOUD_FRAC, PQ_LIQUID, PQ_ICE, PQ_SNOW =', &
331        PCLOUD_FRAC, PQ_LIQUID, PQ_ICE, PQ_SNOW
332  print*,'ZRE_LIQUID_UM, ZRE_ICE_UM =', &
333        ZRE_LIQUID_UM, ZRE_ICE_UM
334  print*,'PAEROSOL_OLD, PAEROSOL =', PAEROSOL_OLD, PAEROSOL
335endif
336
337IF (LHOOK) CALL DR_HOOK('RADIATION_SCHEME',0,ZHOOK_HANDLE)
338print*,'Entree dans radiation_scheme'
339
340!$OMP MASTER
341if (debut_ecrad) then
342! AI appel radiation_setup
343call SETUP_RADIATION_SCHEME(loutput)
344
345 if (lprint_config) then
346  print*,'************* Parametres de configuration  ********************'
347  print*,'rad_config%iverbosesetup = ',rad_config%iverbosesetup
348  print*,'rad_config%iverbose = ',rad_config%iverbose
349  print*,'rad_config%directory_name =', rad_config%directory_name
350  print*,'rad_config%do_lw_derivatives =',rad_config%do_lw_derivatives
351  print*,'rad_config%do_surface_sw_spectral_flux =', &
352        rad_config%do_surface_sw_spectral_flux
353  print*,'rad_config%do_setup_ifsrrtm =', rad_config%do_setup_ifsrrtm
354  print*,'rad_config%i_liq_model =',rad_config%i_liq_model
355  print*,'rad_config%i_ice_model =',rad_config%i_ice_model
356  print*,'rad_config%i_overlap_scheme =', rad_config%i_overlap_scheme
357  print*,'rad_config%use_aerosols = ', rad_config%use_aerosols
358  print*,'rad_config%n_aerosol_types = ', rad_config%n_aerosol_types
359  print*,'rad_config%i_solver_lw =',rad_config%i_solver_lw
360  print*,'rad_config%i_solver_sw =',rad_config%i_solver_sw
361  print*,'rad_config%do_3d_effects =', rad_config%do_3d_effects
362  print*,'rad_config%do_sw_delta_scaling_with_gases =', &
363       rad_config%do_sw_delta_scaling_with_gases
364  print*,'rad_config%do_lw_aerosol_scattering =', &
365       rad_config%do_lw_aerosol_scattering
366  print*,'rad_config%i_albedo_from_band_sw = ', &
367       rad_config%i_albedo_from_band_sw
368  print*,'n_bands_lw =', rad_config%n_bands_lw
369  print*,'rad_config%i_emiss_from_band_lw =', rad_config%i_emiss_from_band_lw
370 endif
371 debut_ecrad=.false.
372endif
373!$OMP END MASTER
374!$OMP BARRIER
375! Fin partie initialisation et configuration
376
377! AI : allocation des tableaux pour chaque partie (thermo, ...)
378!      passage des champs LMDZ aux structures Ecrad
379!      calculs Ecrad
380! AI ATTENTION
381! Allocate memory in radiation objects
382! emissivite avec une seule bande
383CALL single_level%allocate(KLON, NSW, 1, &
384     &                     use_sw_albedo_direct=.TRUE.)
385
386print*,'************* THERMO (allocate + input) ************************************'
387! Set thermodynamic profiles: simply copy over the half-level
388! pressure and temperature
389!print*,'Appel allocate thermo'
390CALL thermodynamics%allocate(KLON, KLEV, use_h2o_sat=.true.)
391!print*,'Definir les champs thermo'
392! AI
393! pressure_hl > paprs
394! temperature_hl calculee dans radlsw de la meme facon que pour RRTM
395thermodynamics%pressure_hl   (KIDIA:KFDIA,:) = PPRESSURE_H   (KIDIA:KFDIA,:)
396thermodynamics%temperature_hl(KIDIA:KFDIA,:) = PTEMPERATURE_H(KIDIA:KFDIA,:)
397
398!print*,'Compute saturation specific humidity'
399! Compute saturation specific humidity, used to hydrate aerosols. The
400! "2" for the last argument indicates that the routine is not being
401! called from within the convection scheme.
402!CALL SATUR(KIDIA, KFDIA, KLON, 1, KLEV, &
403!     &  PPRESSURE, PTEMPERATURE, thermodynamics%h2o_sat_liq, 2)
404! Alternative approximate version using temperature and pressure from
405! the thermodynamics structure
406!CALL thermodynamics%calc_saturation_wrt_liquid(KIDIA, KFDIA)
407!AI ATTENTION
408thermodynamics%h2o_sat_liq = PQSAT
409
410print*,'********** SINGLE LEVEL VARS **********************************'
411!AI ATTENTION
412! Set single-level fileds
413single_level%solar_irradiance              = PSOLAR_IRRADIANCE
414single_level%cos_sza(KIDIA:KFDIA)          = PMU0(KIDIA:KFDIA)
415single_level%skin_temperature(KIDIA:KFDIA) = PTEMPERATURE_SKIN(KIDIA:KFDIA)
416single_level%sw_albedo(KIDIA:KFDIA,:)      = PALBEDO_DIF(KIDIA:KFDIA,:)
417single_level%sw_albedo_direct(KIDIA:KFDIA,:)=PALBEDO_DIR(KIDIA:KFDIA,:)
418single_level%lw_emissivity(KIDIA:KFDIA,1)  = PEMIS(KIDIA:KFDIA)
419!single_level%lw_emissivity(KIDIA:KFDIA,2)  = PEMIS_WINDOW(KIDIA:KFDIA)
420
421! Create the relevant seed from date and time get the starting day
422! and number of minutes since start
423!IDAY = NDD(NINDAT)
424!cur_day
425!ITIM = NINT(NSTEP * YRRIP%TSTEP / 60.0_JPRB)
426ITIM = NINT(TIME / 60.0_JPRB)
427!current_time
428!allocate(single_level%iseed(KIDIA:KFDIA))
429DO JLON = KIDIA, KFDIA
430  ! This method gives a unique value for roughly every 1-km square
431  ! on the globe and every minute.  ASIN(PGEMU)*60 gives rough
432  ! latitude in degrees, which we multiply by 100 to give a unique
433  ! value for roughly every km. PGELAM*60*100 gives a unique number
434  ! for roughly every km of longitude around the equator, which we
435  ! multiply by 180*100 so there is no overlap with the latitude
436  ! values.  The result can be contained in a 32-byte integer (but
437  ! since random numbers are generated with the help of integer
438  ! overflow, it should not matter if the number did overflow).
439  single_level%iseed(JLON) = ITIM + IDAY &
440       &  +  NINT(PGELAM(JLON)*108000000.0_JPRB &
441       &          + ASIN(PGEMU(JLON))*6000.0_JPRB)
442ENDDO
443
444print*,'********** CLOUDS (allocate + input) *******************************************'
445!print*,'Appel Allocate clouds'
446CALL cloud%allocate(KLON, KLEV)
447! Set cloud fields
448cloud%q_liq(KIDIA:KFDIA,:)    = PQ_LIQUID(KIDIA:KFDIA,:)
449cloud%q_ice(KIDIA:KFDIA,:)    = PQ_ICE(KIDIA:KFDIA,:) + PQ_SNOW(KIDIA:KFDIA,:)
450cloud%fraction(KIDIA:KFDIA,:) = PCLOUD_FRAC(KIDIA:KFDIA,:)
451
452!AI ATTENTION a voir avec JL
453! Compute effective radii and convert to metres
454!CALL LIQUID_EFFECTIVE_RADIUS(KIDIA, KFDIA, KLON, KLEV, &
455!     &  PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_LIQUID, PQ_RAIN, &
456!     &  PLAND_SEA_MASK, PCCN_LAND, PCCN_SEA, &
457!     &  ZRE_LIQUID_UM)
458cloud%re_liq(KIDIA:KFDIA,:) = ZRE_LIQUID_UM(KIDIA:KFDIA,:)
459
460!CALL ICE_EFFECTIVE_RADIUS(KIDIA, KFDIA, KLON, KLEV, &
461!     &  PPRESSURE, PTEMPERATURE, PCLOUD_FRAC, PQ_ICE, PQ_SNOW, PGEMU, &
462!     &  ZRE_ICE_UM)
463cloud%re_ice(KIDIA:KFDIA,:) = ZRE_ICE_UM(KIDIA:KFDIA,:)
464
465! Get the cloud overlap decorrelation length (for cloud boundaries),
466! in km, according to the parameterization specified by NDECOLAT,
467! and insert into the "cloud" object. Also get the ratio of
468! decorrelation lengths for cloud water content inhomogeneities and
469! cloud boundaries, and set it in the "rad_config" object.
470!CALL CLOUD_OVERLAP_DECORR_LEN(KIDIA, KFDIA, KLON, PGEMU, YRERAD%NDECOLAT, &
471!     &    ZDECORR_LEN_KM, PDECORR_LEN_RATIO=ZDECORR_LEN_RATIO)
472
473! AI ATTENTION (valeur lue dans namelist)
474!ZDECORR_LEN_RATIO = 0.5_JPRB
475!rad_config%cloud_inhom_decorr_scaling = ZDECORR_LEN_RATIO
476!AI ATTENTION meme valeur que dans offline
477ZDECORR_LEN_KM = 2000.0_JPRB
478DO JLON = KIDIA,KFDIA
479  CALL cloud%set_overlap_param(thermodynamics, &
480       &                       ZDECORR_LEN_KM(JLON), &
481       &                       istartcol=JLON, iendcol=JLON)
482ENDDO
483
484! Cloud water content fractional standard deviation is configurable
485! from namelist NAERAD but must be globally constant. Before it was
486! hard coded at 1.0.
487!CALL cloud%create_fractional_std(KLON, KLEV, YRERAD%RCLOUD_FRAC_STD)
488! AI ATTENTION frac_std=0.75 meme valeur que dans la version offline
489CALL cloud%create_fractional_std(KLON, KLEV, frac_std)
490
491! AI ! Read cloud properties needed by SPARTACUS
492! By default mid and high cloud effective size is 10 km
493!CALL cloud%create_inv_cloud_effective_size(KLON,KLEV,1.0_JPRB/10000.0_JPRB)
494! But for boundary clouds (eta > 0.8) we set it to 1 km
495!DO JLEV = 1,KLEV
496!  DO JLON = KIDIA,KFDIA
497!    IF (PPRESSURE(JLON,JLEV) > 0.8_JPRB * PPRESSURE_H(JLON,KLEV+1)) THEN
498!      cloud%inv_cloud_effective_size(JLON,JLEV) = 1.0e-3_JPRB
499!    ENDIF
500!  ENDDO
501!ENDDO
502!AI ATTENTION meme traitement dans le version offline
503call cloud%create_inv_cloud_effective_size_eta(KLON, KLEV, &
504               &  thermodynamics%pressure_hl, &
505               &  0.005_JPRB, &
506               &  0.0001_JPRB, &
507               &  0.0001, 0.8_jprb, 0.45_jprb)     
508
509print*,'******** AEROSOLS (allocate + input) **************************************'
510!IF (NAERMACC > 0) THEN
511  CALL aerosol%allocate(KLON, 1, KLEV, KAEROSOL) ! MACC climatology
512!ELSE
513!  CALL aerosol%allocate(KLON, 1, KLEV, 6) ! Tegen climatology
514!ENDIF
515! Compute the dry mass of each layer neglecting humidity effects, in
516! kg m-2, needed to scale some of the aerosol inputs
517! AI commente ATTENTION
518!CALL thermodynamics%get_layer_mass(ZLAYER_MASS)
519
520! Copy over aerosol mass mixing ratio
521!IF (NAERMACC > 0) THEN
522
523  ! MACC aerosol climatology - this is already in mass mixing ratio
524  ! units with the required array orientation so we can copy it over
525  ! directly
526  aerosol%mixing_ratio(KIDIA:KFDIA,:,:) = PAEROSOL(KIDIA:KFDIA,:,:)
527
528  ! Add the tropospheric and stratospheric backgrounds contained in the
529  ! old Tegen arrays - this is very ugly!
530! AI ATTENTION
531!  IF (TROP_BG_AER_MASS_EXT > 0.0_JPRB) THEN
532!    aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_TROP_BG_AER) &
533!         &  = aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_TROP_BG_AER) &
534!         &  + PAEROSOL_OLD(KIDIA:KFDIA,1,:) &
535!         &  / (ZLAYER_MASS * TROP_BG_AER_MASS_EXT)
536!  ENDIF
537!  IF (STRAT_BG_AER_MASS_EXT > 0.0_JPRB) THEN
538!    aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_STRAT_BG_AER) &
539!         &  = aerosol%mixing_ratio(KIDIA:KFDIA,:,ITYPE_STRAT_BG_AER) &
540!         &  + PAEROSOL_OLD(KIDIA:KFDIA,6,:) &
541!         &  / (ZLAYER_MASS * STRAT_BG_AER_MASS_EXT)
542!  ENDIF
543
544!ELSE
545
546  ! Tegen aerosol climatology - the array PAEROSOL_OLD contains the
547  ! 550-nm optical depth in each layer. The optics data file
548  ! aerosol_ifs_rrtm_tegen.nc does not contain mass extinction
549  ! coefficient, but a scaling factor that the 550-nm optical depth
550  ! should be multiplied by to obtain the optical depth in each
551  ! spectral band.  Therefore, in order for the units to work out, we
552  ! need to divide by the layer mass (in kg m-2) to obtain the 550-nm
553  ! cross-section per unit mass of dry air (so in m2 kg-1).  We also
554  ! need to permute the array.
555!  DO JLEV = 1,KLEV
556!    DO JAER = 1,6
557!      aerosol%mixing_ratio(KIDIA:KFDIA,JLEV,JAER) &
558!         &  = PAEROSOL_OLD(KIDIA:KFDIA,JAER,JLEV) &
559!         &  / ZLAYER_MASS(KIDIA:KFDIA,JLEV)
560!    ENDDO
561!  ENDDO
562
563!ENDIF
564
565print*,'********** GAS (allocate + input) ************************************************'
566!print*,'Appel Allocate gas'
567CALL gas%allocate(KLON, KLEV)
568
569! Convert ozone Pa*kg/kg to kg/kg
570! AI ATTENTION
571!DO JLEV = 1,KLEV
572!  DO JLON = KIDIA,KFDIA
573!    ZO3(JLON,JLEV) = PO3_DP(JLON,JLEV) &
574!         &         / (PPRESSURE_H(JLON,JLEV+1)-PPRESSURE_H(JLON,JLEV))
575!  ENDDO
576!ENDDO
577
578!  Insert gas mixing ratios
579!print*,'Insert gas mixing ratios'
580CALL gas%put(IH2O,    IMassMixingRatio, PQ)
581CALL gas%put(IO3,     IMassMixingRatio, PO3)
582CALL gas%put_well_mixed(ICO2,    IMAssMixingRatio, PCO2)
583CALL gas%put_well_mixed(ICH4,    IMassMixingRatio, PCH4)
584CALL gas%put_well_mixed(IN2O,    IMassMixingRatio, PN2O)
585CALL gas%put_well_mixed(ICFC11,  IMassMixingRatio, PCFC11)
586CALL gas%put_well_mixed(ICFC12,  IMassMixingRatio, PCFC12)
587CALL gas%put_well_mixed(IHCFC22, IMassMixingRatio, PHCFC22)
588CALL gas%put_well_mixed(ICCL4,   IMassMixingRatio, PCCL4)
589CALL gas%put_well_mixed(IO2,     IMassMixingRatio, PO2)
590! Ensure the units of the gas mixing ratios are what is required by
591! the gas absorption model
592call set_gas_units(rad_config, gas)
593
594print*,'************** FLUX (allocate) ***********************'
595CALL flux%allocate(rad_config, 1, KLON, KLEV)
596
597! Call radiation scheme
598print*,'******** Appel radiation scheme **************************'
599CALL radiation(KLON, KLEV, KIDIA, KFDIA, rad_config, &
600     &  single_level, thermodynamics, gas, cloud, aerosol, flux)
601
602! Compute required output fluxes
603! DN and UP flux
604PFLUX_SW_DN(KIDIA:KFDIA,:) = flux%sw_dn(KIDIA:KFDIA,:)
605PFLUX_SW_UP(KIDIA:KFDIA,:) = flux%sw_up(KIDIA:KFDIA,:)
606PFLUX_LW_DN(KIDIA:KFDIA,:) = flux%lw_dn(KIDIA:KFDIA,:)
607PFLUX_LW_UP(KIDIA:KFDIA,:) = flux%lw_up(KIDIA:KFDIA,:)
608PFLUX_SW_DN_CLEAR(KIDIA:KFDIA,:) = flux%sw_dn_clear(KIDIA:KFDIA,:)
609PFLUX_SW_UP_CLEAR(KIDIA:KFDIA,:) = flux%sw_up_clear(KIDIA:KFDIA,:)
610PFLUX_LW_DN_CLEAR(KIDIA:KFDIA,:) = flux%lw_dn_clear(KIDIA:KFDIA,:)
611PFLUX_LW_UP_CLEAR(KIDIA:KFDIA,:) = flux%lw_up_clear(KIDIA:KFDIA,:)
612
613! First the net fluxes
614PFLUX_SW(KIDIA:KFDIA,:) = flux%sw_dn(KIDIA:KFDIA,:) - flux%sw_up(KIDIA:KFDIA,:)
615PFLUX_LW(KIDIA:KFDIA,:) = flux%lw_dn(KIDIA:KFDIA,:) - flux%lw_up(KIDIA:KFDIA,:)
616PFLUX_SW_CLEAR(KIDIA:KFDIA,:) &
617     &  = flux%sw_dn_clear(KIDIA:KFDIA,:) - flux%sw_up_clear(KIDIA:KFDIA,:)
618PFLUX_LW_CLEAR(KIDIA:KFDIA,:) &
619     &  = flux%lw_dn_clear(KIDIA:KFDIA,:) - flux%lw_up_clear(KIDIA:KFDIA,:)
620
621! Now the surface fluxes
622!PFLUX_SW_DN_SURF(KIDIA:KFDIA) = flux%sw_dn(KIDIA:KFDIA,KLEV+1)
623!PFLUX_LW_DN_SURF(KIDIA:KFDIA) = flux%lw_dn(KIDIA:KFDIA,KLEV+1)
624!PFLUX_SW_UP_SURF(KIDIA:KFDIA) = flux%sw_up(KIDIA:KFDIA,KLEV+1)
625!PFLUX_LW_UP_SURF(KIDIA:KFDIA) = flux%lw_up(KIDIA:KFDIA,KLEV+1)
626!PFLUX_SW_DN_CLEAR_SURF(KIDIA:KFDIA) = flux%sw_dn_clear(KIDIA:KFDIA,KLEV+1)
627!PFLUX_LW_DN_CLEAR_SURF(KIDIA:KFDIA) = flux%lw_dn_clear(KIDIA:KFDIA,KLEV+1)
628!PFLUX_SW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%sw_up_clear(KIDIA:KFDIA,KLEV+1)
629!PFLUX_LW_UP_CLEAR_SURF(KIDIA:KFDIA) = flux%lw_up_clear(KIDIA:KFDIA,KLEV+1)
630PFLUX_DIR(KIDIA:KFDIA) = flux%sw_dn_direct(KIDIA:KFDIA,KLEV+1)
631PFLUX_DIR_CLEAR(KIDIA:KFDIA) = flux%sw_dn_direct_clear(KIDIA:KFDIA,KLEV+1)
632PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = 0.0_JPRB
633WHERE (PMU0(KIDIA:KFDIA) > EPSILON(1.0_JPRB))
634  PFLUX_DIR_INTO_SUN(KIDIA:KFDIA) = PFLUX_DIR(KIDIA:KFDIA) / PMU0(KIDIA:KFDIA)
635END WHERE
636
637! Top-of-atmosphere downwelling flux
638!PFLUX_SW_DN_TOA(KIDIA:KFDIA) = flux%sw_dn(KIDIA:KFDIA,1)
639!PFLUX_SW_UP_TOA(KIDIA:KFDIA) = flux%sw_up(KIDIA:KFDIA,1)
640!PFLUX_LW_DN_TOA(KIDIA:KFDIA) = flux%lw_dn(KIDIA:KFDIA,1)
641!PFLUX_LW_UP_TOA(KIDIA:KFDIA) = flux%lw_up(KIDIA:KFDIA,1)
642
643! Compute UV fluxes as weighted sum of appropriate shortwave bands
644PFLUX_UV       (KIDIA:KFDIA) = 0.0_JPRB
645DO JBAND = 1,NWEIGHT_UV
646  PFLUX_UV(KIDIA:KFDIA) = PFLUX_UV(KIDIA:KFDIA) + WEIGHT_UV(JBAND) &
647       &  * flux%sw_dn_surf_band(IBAND_UV(JBAND),KIDIA:KFDIA)
648ENDDO
649
650! Compute photosynthetically active radiation similarly
651PFLUX_PAR      (KIDIA:KFDIA) = 0.0_JPRB
652PFLUX_PAR_CLEAR(KIDIA:KFDIA) = 0.0_JPRB
653DO JBAND = 1,NWEIGHT_PAR
654  PFLUX_PAR(KIDIA:KFDIA) = PFLUX_PAR(KIDIA:KFDIA) + WEIGHT_PAR(JBAND) &
655       &  * flux%sw_dn_surf_band(IBAND_PAR(JBAND),KIDIA:KFDIA)
656  PFLUX_PAR_CLEAR(KIDIA:KFDIA) = PFLUX_PAR_CLEAR(KIDIA:KFDIA) &
657       &  + WEIGHT_PAR(JBAND) &
658       &  * flux%sw_dn_surf_clear_band(IBAND_PAR(JBAND),KIDIA:KFDIA)
659ENDDO
660
661! Compute effective broadband emissivity
662ZBLACK_BODY_NET_LW = flux%lw_dn(KIDIA:KFDIA,KLEV+1) &
663     &  - RSIGMA*PTEMPERATURE_SKIN(KIDIA:KFDIA)**4
664PEMIS_OUT(KIDIA:KFDIA) = PEMIS(KIDIA:KFDIA)
665WHERE (ABS(ZBLACK_BODY_NET_LW) > 1.0E-5)
666  PEMIS_OUT(KIDIA:KFDIA) = PFLUX_LW(KIDIA:KFDIA,KLEV+1) / ZBLACK_BODY_NET_LW
667END WHERE
668
669! Copy longwave derivatives
670! AI ATTENTION
671!IF (YRERAD%LAPPROXLWUPDATE) THEN
672IF (rad_config%do_lw_derivatives) THEN
673  PLWDERIVATIVE(KIDIA:KFDIA,:) = flux%lw_derivatives(KIDIA:KFDIA,:)
674END IF
675
676! Store the shortwave downwelling fluxes in each albedo band
677!AI ATTENTION
678!IF (YRERAD%LAPPROXSWUPDATE) THEN
679IF (rad_config%do_surface_sw_spectral_flux) THEN
680  PSWDIFFUSEBAND(KIDIA:KFDIA,:) = 0.0_JPRB
681  PSWDIRECTBAND (KIDIA:KFDIA,:) = 0.0_JPRB
682  DO JBAND = 1,rad_config%n_bands_sw
683    JB_ALBEDO = rad_config%i_albedo_from_band_sw(JBAND)
684    DO JLON = KIDIA,KFDIA
685      PSWDIFFUSEBAND(JLON,JB_ALBEDO) = PSWDIFFUSEBAND(JLON,JB_ALBEDO) &
686           &  + flux%sw_dn_surf_band(JBAND,JLON) &
687           &  - flux%sw_dn_direct_surf_band(JBAND,JLON)
688      PSWDIRECTBAND(JLON,JB_ALBEDO)  = PSWDIRECTBAND(JLON,JB_ALBEDO) &
689           &  + flux%sw_dn_direct_surf_band(JBAND,JLON)
690    ENDDO
691  ENDDO
692ENDIF
693
694CALL single_level%deallocate
695CALL thermodynamics%deallocate
696CALL gas%deallocate
697CALL cloud%deallocate
698CALL aerosol%deallocate
699CALL flux%deallocate
700
701IF (LHOOK) CALL DR_HOOK('RADIATION_SCHEME',1,ZHOOK_HANDLE)
702
703END SUBROUTINE RADIATION_SCHEME
Note: See TracBrowser for help on using the repository browser.