source: branches/publications/ORCHIDEE_CAMEO_gmd_2022/src_sechiba/sechiba.f90 @ 7693

Last change on this file since 7693 was 6608, checked in by josefine.ghattas, 4 years ago

Followed coding guidelines: remove mcs and mcfc as public variables in hydrol module as introduced in [6606].

  • Property svn:keywords set to Revision Date HeadURL Date Author Revision
File size: 117.9 KB
Line 
1!  ==============================================================================================================================\n
2!  MODULE       : sechiba
3!
4!  CONTACT      : orchidee-help _at_ listes.ipsl.fr
5!
6!  LICENCE      : IPSL (2006)
7!  This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
8!
9!>\BRIEF        Structures the calculation of atmospheric and hydrological
10!! variables by calling diffuco_main, enerbil_main, hydrol_main,
11!! condveg_main and thermosoil_main. Note that sechiba_main
12!! calls slowproc_main and thus indirectly calculates the biogeochemical
13!! processes as well.
14!!
15!!\n DESCRIPTION  : :: shumdiag, :: litterhumdiag and :: stempdiag are not
16!! saved in the restart file because at the first time step because they
17!! are recalculated. However, they must be saved as they are in slowproc
18!! which is called before the modules which calculate them.
19!!
20!! RECENT CHANGE(S): None
21!!
22!! REFERENCE(S) : None
23!!   
24!! SVN     :
25!! $HeadURL$
26!! $Date$
27!! $Revision$
28!! \n
29!_ ================================================================================================================================
30 
31MODULE sechiba
32 
33  USE ioipsl
34  USE xios_orchidee
35  USE constantes
36  USE time, ONLY : one_day, dt_sechiba
37  USE constantes_soil
38  USE pft_parameters
39  USE grid
40  USE diffuco
41  USE condveg
42  USE enerbil
43  USE hydrol
44  USE thermosoil
45  USE sechiba_io_p
46  USE slowproc
47  USE routing
48  USE ioipsl_para
49  USE chemistry
50
51  IMPLICIT NONE
52
53  PRIVATE
54  PUBLIC sechiba_main, sechiba_initialize, sechiba_clear, sechiba_interface_orchidee_inca, sechiba_xios_initialize
55
56  INTEGER(i_std), SAVE                             :: printlev_loc   !! local printlev for this module
57!$OMP THREADPRIVATE(printlev_loc)
58  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexveg       !! indexing array for the 3D fields of vegetation
59!$OMP THREADPRIVATE(indexveg)
60  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexlai       !! indexing array for the 3D fields of vegetation
61!$OMP THREADPRIVATE(indexlai)
62  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexnobio     !! indexing array for the 3D fields of other surfaces (ice,
63                                                                     !! lakes, ...)
64!$OMP THREADPRIVATE(indexnobio)
65  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexsoil      !! indexing array for the 3D fields of soil types (kjpindex*nstm)
66!$OMP THREADPRIVATE(indexsoil)
67  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexgrnd      !! indexing array for the 3D ground heat profiles (kjpindex*ngrnd)
68!$OMP THREADPRIVATE(indexgrnd)
69  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexlayer     !! indexing array for the 3D fields of soil layers in CWRR (kjpindex*nslm)
70!$OMP THREADPRIVATE(indexlayer)
71  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexnslm      !! indexing array for the 3D fields of diagnostic soil layers (kjpindex*nslm)
72!$OMP THREADPRIVATE(indexnslm)
73  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexalb       !! indexing array for the 2 fields of albedo
74!$OMP THREADPRIVATE(indexalb)
75  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: indexsnow      !! indexing array for the 3D fields snow layers
76!$OMP THREADPRIVATE(indexsnow)
77  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: veget          !! Fraction of vegetation type (unitless, 0-1)       
78!$OMP THREADPRIVATE(veget)
79  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: veget_max      !! Max. fraction of vegetation type (LAI -> infty, unitless)
80!$OMP THREADPRIVATE(veget_max)
81  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: tot_bare_soil  !! Total evaporating bare soil fraction
82!$OMP THREADPRIVATE(tot_bare_soil)
83  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: height         !! Vegetation Height (m)
84!$OMP THREADPRIVATE(height)
85  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: totfrac_nobio  !! Total fraction of continental ice+lakes+cities+...
86                                                                     !! (unitless, 0-1)
87!$OMP THREADPRIVATE(totfrac_nobio)
88  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: floodout       !! Flow out of floodplains from hydrol
89!$OMP THREADPRIVATE(floodout)
90  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: runoff         !! Surface runoff calculated by hydrol
91                                                                     !! @tex $(kg m^{-2})$ @endtex
92!$OMP THREADPRIVATE(runoff)
93  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: drainage       !! Deep drainage calculatedd by hydrol
94                                                                     !! @tex $(kg m^{-2})$ @endtex
95!$OMP THREADPRIVATE(drainage)
96  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: returnflow     !! Water flow from lakes and swamps which returns to
97                                                                     !! the grid box @tex $(kg m^{-2})$ @endtex
98!$OMP THREADPRIVATE(returnflow)
99  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: reinfiltration !! Routed water which returns into the soil
100!$OMP THREADPRIVATE(reinfiltration)
101  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: irrigation     !! Irrigation flux taken from the routing reservoirs and
102                                                                     !! being put into the upper layers of the soil
103                                                                     !! @tex $(kg m^{-2})$ @endtex
104!$OMP THREADPRIVATE(irrigation)
105  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: emis           !! Surface emissivity (unitless)
106!$OMP THREADPRIVATE(emis)
107  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: z0h           !! Surface roughness for heat (m)
108!$OMP THREADPRIVATE(z0h)
109  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: z0m           !! Surface roughness for momentum (m)
110!$OMP THREADPRIVATE(z0m)
111  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: roughheight    !! Effective height for roughness (m)
112!$OMP THREADPRIVATE(roughheight)
113  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: reinf_slope    !! slope coefficient (reinfiltration)
114!$OMP THREADPRIVATE(reinf_slope)
115  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: shumdiag       !! Mean relative soil moisture in the different levels used
116                                                                     !! by thermosoil.f90 (unitless, 0-1)
117!$OMP THREADPRIVATE(shumdiag)
118  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: shumdiag_perma !! Saturation degree of the soil
119!$OMP THREADPRIVATE(shumdiag_perma)
120  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: k_litt         !! litter cond.
121!$OMP THREADPRIVATE(k_litt)
122  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: litterhumdiag  !! Litter dryness factor (unitless, 0-1)
123!$OMP THREADPRIVATE(litterhumdiag)
124  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: stempdiag      !! Temperature which controls canopy evolution (K)
125!$OMP THREADPRIVATE(stempdiag)
126  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: qsintveg       !! Water on vegetation due to interception
127                                                                     !! @tex $(kg m^{-2})$ @endtex
128!$OMP THREADPRIVATE(qsintveg)
129  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: vbeta2         !! Interception resistance (unitless,0-1)
130!$OMP THREADPRIVATE(vbeta2)
131  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: vbeta3         !! Vegetation resistance (unitless,0-1)
132!$OMP THREADPRIVATE(vbeta3)
133  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: vbeta3pot      !! Potential vegetation resistance
134!$OMP THREADPRIVATE(vbeta3pot)
135  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: gsmean         !! Mean stomatal conductance for CO2 (mol m-2 s-1)
136!$OMP THREADPRIVATE(gsmean)
137  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: cimean         !! STOMATE: mean intercellular CO2 concentration (ppm)
138!$OMP THREADPRIVATE(cimean)
139  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: vevapwet       !! Interception loss over each PFT
140                                                                     !! @tex $(kg m^{-2} days^{-1})$ @endtex
141!$OMP THREADPRIVATE(vevapwet)
142  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: transpir       !! Transpiration @tex $(kg m^{-2} days^{-1})$ @endtex
143!$OMP THREADPRIVATE(transpir)
144  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: transpot       !! Potential Transpiration (needed for irrigation)
145!$OMP THREADPRIVATE(transpot)
146  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: qsintmax       !! Maximum amount of water in the canopy interception
147                                                                     !! reservoir @tex $(kg m^{-2})$ @endtex
148!$OMP THREADPRIVATE(qsintmax)
149  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: rveget         !! Surface resistance for the vegetation
150                                                                     !! @tex $(s m^{-1})$ @endtex
151!$OMP THREADPRIVATE(rveget)
152  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: rstruct        !! Vegetation structural resistance
153!$OMP THREADPRIVATE(rstruct)
154  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: snow_nobio     !! Snow mass of non-vegetative surfaces
155                                                                     !! @tex $(kg m^{-2})$ @endtex
156!$OMP THREADPRIVATE(snow_nobio)
157  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: snow_nobio_age !! Snow age on non-vegetative surfaces (days)
158!$OMP THREADPRIVATE(snow_nobio_age)
159  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: frac_nobio     !! Fraction of non-vegetative surfaces (continental ice,
160                                                                     !! lakes, ...) (unitless, 0-1)
161!$OMP THREADPRIVATE(frac_nobio)
162  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:,:):: assim_param    !! vcmax, nue, and leaf N for photosynthesis  for photosynthesis
163!$OMP THREADPRIVATE(assim_param)
164  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: lai            !! Surface foliaire
165!$OMP THREADPRIVATE(lai)
166  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: gpp            !! STOMATE: GPP. gC/m**2 of total area
167!$OMP THREADPRIVATE(gpp)
168  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)     :: temp_growth    !! Growth temperature (C) - Is equal to t2m_month
169!$OMP THREADPRIVATE(temp_growth)
170  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: humrel         !! Relative humidity
171!$OMP THREADPRIVATE(humrel)
172  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: vegstress      !! Vegetation moisture stress (only for vegetation growth)
173!$OMP THREADPRIVATE(vegstress)
174  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:,:):: frac_age       !! Age efficacity from STOMATE for isoprene
175!$OMP THREADPRIVATE(frac_age)
176  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: soiltile       !! Fraction of each soil tile (0-1, unitless)
177!$OMP THREADPRIVATE(soiltile)
178  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: fraclut        !! Fraction of each landuse tile (0-1, unitless)
179!$OMP THREADPRIVATE(fraclut)
180  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: nwdFraclut     !! Fraction of non-woody vegetation in each landuse tile (0-1, unitless)
181!$OMP THREADPRIVATE(nwdFraclut)
182  INTEGER(i_std), ALLOCATABLE, SAVE, DIMENSION (:) :: njsc           !! Index of the dominant soil textural class in the grid cell (1-nscm, unitless)
183!$OMP THREADPRIVATE(njsc)
184  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vbeta1         !! Snow resistance
185!$OMP THREADPRIVATE(vbeta1)
186  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vbeta4         !! Bare soil resistance
187!$OMP THREADPRIVATE(vbeta4)
188  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vbeta5         !! Floodplains resistance
189!$OMP THREADPRIVATE(vbeta5)
190  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: soilcap        !!
191!$OMP THREADPRIVATE(soilcap)
192  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: soilflx        !!
193!$OMP THREADPRIVATE(soilflx)
194  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: temp_sol       !! Soil temperature
195!$OMP THREADPRIVATE(temp_sol)
196  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: qsurf          !! near soil air moisture
197!$OMP THREADPRIVATE(qsurf)
198  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: flood_res      !! flood reservoir estimate
199!$OMP THREADPRIVATE(flood_res)
200  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: flood_frac     !! flooded fraction
201!$OMP THREADPRIVATE(flood_frac)
202  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: snow           !! Snow mass [Kg/m^2]
203!$OMP THREADPRIVATE(snow)
204  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: snow_age       !! Snow age @tex ($d$) @endtex
205!$OMP THREADPRIVATE(snow_age)
206  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: drysoil_frac   !! Fraction of visibly (albedo) Dry soil (Between 0 and 1)
207!$OMP THREADPRIVATE(drysoil_frac)
208  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: evap_bare_lim  !! Bare soil stress
209!$OMP THREADPRIVATE(evap_bare_lim)
210  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: evap_bare_lim_ns !! Bare soil stress
211!$OMP THREADPRIVATE(evap_bare_lim_ns)
212  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: co2_flux       !! CO2 flux (gC/m**2 of average ground/one_day)
213!$OMP THREADPRIVATE(co2_flux)
214  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: co2_to_bm      !! virtual CO2 flux (gC/m**2 of average ground/s)
215!$OMP THREADPRIVATE(co2_to_bm)
216  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: evapot         !! Soil Potential Evaporation
217!$OMP THREADPRIVATE(evapot)
218  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: evapot_corr    !! Soil Potential Evaporation Correction (Milly 1992)
219!$OMP THREADPRIVATE(evapot_corr)
220  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vevapflo       !! Floodplains evaporation
221!$OMP THREADPRIVATE(vevapflo)
222  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vevapsno       !! Snow evaporation
223!$OMP THREADPRIVATE(vevapsno)
224  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vevapnu        !! Bare soil evaporation
225!$OMP THREADPRIVATE(vevapnu)
226  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: tot_melt       !! Total melt
227!$OMP THREADPRIVATE(tot_melt)
228  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: vbeta          !! Resistance coefficient
229!$OMP THREADPRIVATE(vbeta)
230  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: rau            !! Density
231!$OMP THREADPRIVATE(rau)
232  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: deadleaf_cover !! Fraction of soil covered by dead leaves
233!$OMP THREADPRIVATE(deadleaf_cover)
234  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)    :: ptnlev1        !! 1st level Different levels soil temperature
235!$OMP THREADPRIVATE(ptnlev1)
236  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: mc_layh        !! Volumetric soil moisture for each layer in hydrol(liquid + ice) (m3/m3)
237!$OMP THREADPRIVATE(mc_layh)
238  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: mcl_layh       !! Volumetric soil moisture for each layer in hydrol(liquid) (m3/m3)
239!$OMP THREADPRIVATE(mcl_layh)
240  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)  :: soilmoist      !! Total soil moisture content for each layer in hydrol(liquid + ice) (mm)
241!$OMP THREADPRIVATE(soilmoist)
242  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:,:):: mc_layh_s      !! Volumetric soil moisture for each layer in hydrol per soiltile (liquid + ice) (m3/m3)
243!$OMP THREADPRIVATE(mc_layh_s)
244  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:,:):: mcl_layh_s     !! Volumetric soil moisture for each layer in hydrol per soiltile (liquid) (m3/m3)
245!$OMP THREADPRIVATE(mcl_layh_s)
246  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:,:):: soilmoist_s    !! Total soil moisture content for each layer in hydrol per soiltile (liquid + ice) (mm)
247!$OMP THREADPRIVATE(soilmoist_s)
248
249  LOGICAL, SAVE                                    :: l_first_sechiba = .TRUE. !! Flag controlling the intialisation (true/false)
250!$OMP THREADPRIVATE(l_first_sechiba)
251
252  ! Variables related to snow processes calculations 
253
254  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)    :: frac_snow_veg   !! Snow cover fraction on vegetation (unitless)
255!$OMP THREADPRIVATE(frac_snow_veg)
256  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: frac_snow_nobio !! Snow cover fraction on continental ice, lakes, etc (unitless)
257!$OMP THREADPRIVATE(frac_snow_nobio)
258  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: snowrho      !! snow density for each layer (Kg/m^3)
259!$OMP THREADPRIVATE(snowrho)
260  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: snowheat     !! snow heat content for each layer (J/m2)
261!$OMP THREADPRIVATE(snowheat)
262  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: snowgrain    !! snow grain size (m)
263!$OMP THREADPRIVATE(snowgrain)
264  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: snowtemp     !! snow temperature profile (K)
265!$OMP THREADPRIVATE(snowtemp)
266  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: snowdz       !! snow layer thickness (m)
267!$OMP THREADPRIVATE(snowdz)
268  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)    :: gtemp        !! soil surface temperature
269!$OMP THREADPRIVATE(gtemp)
270  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:)   :: pgflux       !! net energy into snow pack
271!$OMP THREADPRIVATE(pgflux)
272  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: cgrnd_snow   !! Integration coefficient for snow numerical scheme
273!$OMP THREADPRIVATE(cgrnd_snow)
274  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:)  :: dgrnd_snow   !! Integration coefficient for snow numerical scheme
275!$OMP THREADPRIVATE(dgrnd_snow)
276  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)    :: lambda_snow  !! Coefficient of the linear extrapolation of surface temperature
277                                                                  !! from the first and second snow layers
278!$OMP THREADPRIVATE(lambda_snow)
279  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)    :: temp_sol_add !! Additional energy to melt snow for snow ablation case (K)
280!$OMP THREADPRIVATE(temp_sol_add)
281
282! Variables related to carbon soil discretization
283  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:,:)     :: tdeep            !! Deep temperature profile (K)
284!$OMP THREADPRIVATE(tdeep)
285  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:,:)     :: hsdeep           !! Deep soil humidity profile (unitless)
286!$OMP THREADPRIVATE(hsdeep)
287  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:,:)     :: heat_Zimov       !! Heating associated with decomposition [W/m**3 soil]
288!$OMP THREADPRIVATE(heat_Zimov)
289  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)         :: sfluxCH4_deep    !! Surface flux of CH4 to atmosphere from permafrost
290!$OMP THREADPRIVATE(sfluxCH4_deep)
291  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)         :: sfluxCO2_deep    !! Surface flux of CO2 to atmosphere from permafrost
292!$OMP THREADPRIVATE(sfluxCO2_deep)
293  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:)   :: som_total        !! total soil organic matter for use in thermal calcs (g/m**3)
294!$OMP THREADPRIVATE(som_total)
295  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION (:,:)      :: altmax           !! Maximul active layer thickness (m). Be careful, here active means non frozen.
296                                                                           !! Not related with the active soil carbon pool.
297!$OMP THREADPRIVATE(altmax)
298  REAL(r_std), ALLOCATABLE, SAVE, DIMENSION(:)         :: depth_organic_soil !! Depth at which there is still organic matter (m)
299!$OMP THREADPRIVATE(depth_organic_soil)
300
301CONTAINS
302
303
304!!  =============================================================================================================================
305!! SUBROUTINE:    sechiba_xios_initialize
306!!
307!>\BRIEF          Initialize xios dependant defintion before closing context defintion
308!!
309!! DESCRIPTION:   Initialize xios dependant defintion before closing context defintion
310!!
311!! \n
312!_ ==============================================================================================================================
313
314  SUBROUTINE sechiba_xios_initialize
315    LOGICAL :: lerr
316   
317    IF (xios_orchidee_ok) THEN 
318      lerr=xios_orchidee_setvar('min_sechiba',min_sechiba)
319      CALL slowproc_xios_initialize
320      CALL condveg_xios_initialize
321      CALL chemistry_xios_initialize
322      CALL thermosoil_xios_initialize
323    END IF
324    IF (printlev_loc>=3) WRITE(numout,*) 'End sechiba_xios_initialize'
325
326  END SUBROUTINE sechiba_xios_initialize
327
328
329
330
331!!  =============================================================================================================================
332!! SUBROUTINE:    sechiba_initialize
333!!
334!>\BRIEF          Initialize all prinicipal modules by calling their "_initialize" subroutines
335!!
336!! DESCRIPTION:   Initialize all prinicipal modules by calling their "_initialize" subroutines
337!!
338!! \n
339!_ ==============================================================================================================================
340
341  SUBROUTINE sechiba_initialize( &
342       kjit,         kjpij,        kjpindex,     index,                   &
343       lalo,         contfrac,     neighbours,   resolution, zlev,        &
344       u,            v,            qair,         temp_air,    &
345       petAcoef,     peqAcoef,     petBcoef,     peqBcoef,                &
346       precip_rain,  precip_snow,  lwdown,       swnet,      swdown,      &
347       pb,           rest_id,      hist_id,      hist2_id,                &
348       rest_id_stom, hist_id_stom, hist_id_stom_IPCC,                     &
349       coastalflow,  riverflow,    tsol_rad,     vevapp,       qsurf_out, &
350       z0m_out,      z0h_out,      albedo,       fluxsens,     fluxlat,      emis_out,  &
351       temp_sol_new, tq_cdrag)
352
353!! 0.1 Input variables
354    INTEGER(i_std), INTENT(in)                               :: kjit              !! Time step number (unitless)
355    INTEGER(i_std), INTENT(in)                               :: kjpij             !! Total size of the un-compressed grid
356                                                                                  !! (unitless)
357    INTEGER(i_std), INTENT(in)                               :: kjpindex          !! Domain size - terrestrial pixels only
358                                                                                  !! (unitless)
359    INTEGER(i_std),INTENT (in)                               :: rest_id           !! _Restart_ file identifier (unitless)
360    INTEGER(i_std),INTENT (in)                               :: hist_id           !! _History_ file identifier (unitless)
361    INTEGER(i_std),INTENT (in)                               :: hist2_id          !! _History_ file 2 identifier (unitless)
362    INTEGER(i_std),INTENT (in)                               :: rest_id_stom      !! STOMATE's _Restart_ file identifier
363                                                                                  !! (unitless)
364    INTEGER(i_std),INTENT (in)                               :: hist_id_stom      !! STOMATE's _History_ file identifier
365                                                                                  !! (unitless)
366    INTEGER(i_std),INTENT(in)                                :: hist_id_stom_IPCC !! STOMATE's IPCC _history_ file file
367                                                                                  !! identifier (unitless)
368    REAL(r_std),DIMENSION (:,:), INTENT (in)          :: lalo              !! Geographic coordinates (latitude,longitude)
369                                                                                  !! for grid cells (degrees)
370    REAL(r_std),DIMENSION (:), INTENT (in)            :: contfrac          !! Fraction of continent in the grid
371                                                                                  !! (unitless, 0-1)
372    INTEGER(i_std),DIMENSION (:), INTENT (in)         :: index             !! Indices of the pixels on the map.
373                                                                                  !! Sechiba uses a reduced grid excluding oceans
374                                                                                  !! ::index contains the indices of the
375                                                                                  !! terrestrial pixels only! (unitless)
376    INTEGER(i_std), DIMENSION (:,:), INTENT(in):: neighbours        !! Neighboring grid points if land!(unitless)
377    REAL(r_std), DIMENSION (:,:), INTENT(in)          :: resolution        !! Size in x and y of the grid (m)
378   
379    REAL(r_std),DIMENSION (:), INTENT (in)            :: u                 !! Lowest level wind speed in direction u
380                                                                                  !! @tex $(m.s^{-1})$ @endtex
381    REAL(r_std),DIMENSION (:), INTENT (in)            :: v                 !! Lowest level wind speed in direction v
382                                                                                  !! @tex $(m.s^{-1})$ @endtex
383    REAL(r_std),DIMENSION (:), INTENT (in)            :: zlev              !! Height of first layer (m)
384    REAL(r_std),DIMENSION (:), INTENT (in)            :: qair              !! Lowest level specific humidity
385                                                                                  !! @tex $(kg kg^{-1})$ @endtex
386    REAL(r_std),DIMENSION (:), INTENT (in)            :: precip_rain       !! Rain precipitation
387                                                                                  !! @tex $(kg m^{-2})$ @endtex
388    REAL(r_std),DIMENSION (:), INTENT (in)            :: precip_snow       !! Snow precipitation
389                                                                                  !! @tex $(kg m^{-2})$ @endtex
390    REAL(r_std),DIMENSION (:), INTENT (in)            :: lwdown            !! Down-welling long-wave flux
391                                                                                  !! @tex $(W m^{-2})$ @endtex
392    REAL(r_std),DIMENSION (:), INTENT (in)            :: swnet             !! Net surface short-wave flux
393                                                                                  !! @tex $(W m^{-2})$ @endtex
394    REAL(r_std),DIMENSION (:), INTENT (in)            :: swdown            !! Down-welling surface short-wave flux
395                                                                                  !! @tex $(W m^{-2})$ @endtex
396    REAL(r_std),DIMENSION (:), INTENT (in)            :: temp_air          !! Air temperature (K)
397    REAL(r_std),DIMENSION (:), INTENT (in)            :: petAcoef          !! Coefficients A for T from the Planetary
398                                                                                  !! Boundary Layer
399    REAL(r_std),DIMENSION (:), INTENT (in)            :: peqAcoef          !! Coefficients A for q from the Planetary
400                                                                                  !! Boundary Layer
401    REAL(r_std),DIMENSION (:), INTENT (in)            :: petBcoef          !! Coefficients B for T from the Planetary
402                                                                                  !! Boundary Layer
403    REAL(r_std),DIMENSION (:), INTENT (in)            :: peqBcoef          !! Coefficients B for q from the Planetary
404                                                                                  !! Boundary Layer
405    REAL(r_std),DIMENSION (:), INTENT (in)            :: pb                !! Surface pressure (hPa)
406
407
408!! 0.2 Output variables
409    REAL(r_std),DIMENSION (:), INTENT (out)           :: coastalflow       !! Outflow on coastal points by small basins.
410                                                                                  !! This is the water which flows in a disperse
411                                                                                  !! way into the ocean
412                                                                                  !! @tex $(kg dt_routing^{-1})$ @endtex
413    REAL(r_std),DIMENSION (:), INTENT (out)           :: riverflow         !! Outflow of the major rivers.
414                                                                                  !! The flux will be located on the continental
415                                                                                  !! grid but this should be a coastal point 
416                                                                                  !! @tex $(kg dt_routing^{-1})$ @endtex
417    REAL(r_std),DIMENSION (:), INTENT (out)           :: tsol_rad          !! Radiative surface temperature
418                                                                                  !! @tex $(W m^{-2})$ @endtex
419    REAL(r_std),DIMENSION (:), INTENT (out)           :: vevapp            !! Total of evaporation
420                                                                                  !! @tex $(kg m^{-2} days^{-1})$ @endtex
421   
422    REAL(r_std),DIMENSION (:), INTENT (out)           :: qsurf_out         !! Surface specific humidity
423                                                                                  !! @tex $(kg kg^{-1})$ @endtex
424    REAL(r_std),DIMENSION (:), INTENT (out)           :: z0m_out           !! Surface roughness momentum (output diagnostic, m)
425    REAL(r_std),DIMENSION (:), INTENT (out)           :: z0h_out           !! Surface roughness heat (output diagnostic, m)
426    REAL(r_std),DIMENSION (:,:), INTENT (out)         :: albedo            !! Surface albedo for visible and near-infrared
427                                                                                  !! (unitless, 0-1)
428    REAL(r_std),DIMENSION (:), INTENT (out)           :: fluxsens          !! Sensible heat flux
429                                                                                  !! @tex $(W m^{-2})$ @endtex
430    REAL(r_std),DIMENSION (:), INTENT (out)           :: fluxlat           !! Latent heat flux
431                                                                                  !! @tex $(W m^{-2})$ @endtex
432    REAL(r_std),DIMENSION (:), INTENT (out)           :: emis_out          !! Emissivity (unitless)
433    REAL(r_std),DIMENSION (:), INTENT (out)           :: temp_sol_new      !! New ground temperature (K)
434
435!! 0.3 Modified
436    REAL(r_std),DIMENSION (:), INTENT (inout)         :: tq_cdrag          !! Surface drag coefficient (-)
437
438!! 0.4 Local variables
439    INTEGER(i_std)                                           :: ji, jv            !! Index (unitless)
440    REAL(r_std), DIMENSION(kjpindex)                         :: zmaxh_glo         !! 2D field of constant soil depth (zmaxh) (m)
441    CHARACTER(LEN=80)                                        :: var_name          !! To store variables names for I/O (unitless)
442    REAL(r_std), DIMENSION(kjpindex,nslm,nvm)                :: mc_layh_pft       !! As mc_layh per pft
443    REAL(r_std), DIMENSION(kjpindex,nslm,nvm)                :: mcl_layh_pft      !! As mcl_layh per pft
444    REAL(r_std), DIMENSION(kjpindex,nslm,nvm)                :: soilmoist_pft     !! As soilmoist per pft
445
446!_ ================================================================================================================================
447
448    IF (printlev>=3) WRITE(numout,*) 'Start sechiba_initialize'
449
450    !! 1. Initialize variables on first call
451
452    !! 1.1 Initialize most of sechiba's variables
453    CALL sechiba_init (kjit, kjpij, kjpindex, index, rest_id, lalo)
454   
455    !! 1.3 Initialize stomate's variables
456
457    CALL slowproc_initialize (kjit,          kjpij,        kjpindex,                          &
458                              rest_id,       rest_id_stom, hist_id_stom,   hist_id_stom_IPCC, &
459                              index,         indexveg,     lalo,           neighbours,        &
460                              resolution,    contfrac,     temp_air,                          &
461                              soiltile,      reinf_slope,  deadleaf_cover, assim_param,       &
462                              lai,           frac_age,     height,         veget,             &
463                              frac_nobio,    njsc,         veget_max,      fraclut,           &
464                              nwdfraclut,    tot_bare_soil,totfrac_nobio,  qsintmax,          &
465                              co2_to_bm,     temp_growth,                                     &
466                              som_total,     heat_Zimov,   altmax,         depth_organic_soil)
467   
468    !! 1.4 Initialize diffusion coefficients
469    CALL diffuco_initialize (kjit,    kjpindex, index,                  &
470                             rest_id, lalo,     neighbours, resolution, &
471                             rstruct, tq_cdrag)
472   
473    !! 1.5 Initialize remaining variables of energy budget
474    CALL enerbil_initialize (kjit,     kjpindex,     index,    rest_id,  &
475                             qair,                                       &
476                             temp_sol, temp_sol_new, tsol_rad,           &
477                             evapot,   evapot_corr,  qsurf,    fluxsens, &
478                             fluxlat,  vevapp )
479   
480   
481    !! 1.7 Initialize remaining hydrological variables
482    CALL hydrol_initialize ( kjit,           kjpindex,  index,         rest_id,          &
483         njsc,           soiltile,  veget,         veget_max,        &
484         altmax,                                                     &
485         humrel,         vegstress, drysoil_frac,                    &
486         shumdiag_perma,    qsintveg,                                &
487         evap_bare_lim, evap_bare_lim_ns,  snow,      snow_age,      snow_nobio,       &
488         snow_nobio_age, snowrho,   snowtemp,      snowgrain,        &
489         snowdz,         snowheat,  &
490         mc_layh,    mcl_layh,  soilmoist, mc_layh_s,  mcl_layh_s, soilmoist_s)
491
492   
493    !! 1.9 Initialize surface parameters (emissivity, albedo and roughness)
494    CALL condveg_initialize (kjit, kjpindex, index, rest_id, &
495         lalo, neighbours, resolution, contfrac, veget, veget_max, frac_nobio, totfrac_nobio, &
496         zlev, snow, snow_age, snow_nobio, snow_nobio_age, &
497         drysoil_frac, height, snowdz,snowrho, tot_bare_soil, &
498         temp_air, pb, u, v, lai, &
499         emis, albedo, z0m, z0h, roughheight, &
500         frac_snow_veg,frac_snow_nobio)
501   
502
503    !! 1.10 Initialization of soil thermodynamics
504    DO jv = 1,nvm
505       mc_layh_pft(:,:,jv) = mc_layh_s(:,:,pref_soil_veg(jv))
506       mcl_layh_pft(:,:,jv) = mcl_layh_s(:,:,pref_soil_veg(jv))
507       soilmoist_pft(:,:,jv) = soilmoist_s(:,:,pref_soil_veg(jv))
508    END DO
509    CALL thermosoil_initialize (kjit, kjpindex, rest_id,  &
510         temp_sol_new, snow,       shumdiag_perma,        &
511         soilcap,      soilflx,    depth_organic_soil,    & 
512         stempdiag,    gtemp,                             &
513         mc_layh,  mcl_layh,   soilmoist,       njsc ,    &
514         frac_snow_veg,frac_snow_nobio,totfrac_nobio,     &
515         snowdz, snowrho, snowtemp, lambda_snow, cgrnd_snow, dgrnd_snow, pb, &
516         som_total, &
517         veget_max, mc_layh_pft, mcl_layh_pft, soilmoist_pft)
518
519
520    !! 1.12 Initialize river routing
521    IF ( river_routing .AND. nbp_glo .GT. 1) THEN
522       !! 1.12.1 Initialize river routing
523       CALL routing_initialize( kjit,        kjpindex,       index,                 &
524            rest_id,     hist_id,        hist2_id,   lalo,      &
525            neighbours,  resolution,     contfrac,   stempdiag, &
526            returnflow,  reinfiltration, irrigation, riverflow, &
527            coastalflow, flood_frac,     flood_res )
528    ELSE
529       !! 1.12.2 No routing, set variables to zero
530       riverflow(:) = zero
531       coastalflow(:) = zero
532       returnflow(:) = zero
533       reinfiltration(:) = zero
534       irrigation(:) = zero
535       flood_frac(:) = zero
536       flood_res(:) = zero
537    ENDIF
538   
539    !! 1.13 Write internal variables to output fields
540    z0m_out(:) = z0m(:)
541    z0h_out(:) = z0h(:)
542    emis_out(:) = emis(:) 
543    qsurf_out(:) = qsurf(:)
544
545    !! 2. Output variables only once
546    zmaxh_glo(:) = zmaxh
547    CALL xios_orchidee_send_field("zmaxh",zmaxh_glo)
548
549    IF (printlev_loc>=3) WRITE(numout,*) 'sechiba_initialize done'
550
551  END SUBROUTINE sechiba_initialize
552
553!! ==============================================================================================================================\n
554!! SUBROUTINE   : sechiba_main
555!!
556!>\BRIEF        Main routine for the sechiba module performing three functions:
557!! calculating temporal evolution of all variables and preparation of output and
558!! restart files (during the last call only)
559!!
560!!\n DESCRIPTION : Main routine for the sechiba module.
561!! One time step evolution consists of:
562!! - call sechiba_var_init to do some initialization,
563!! - call slowproc_main to do some daily calculations
564!! - call diffuco_main for diffusion coefficient calculation,
565!! - call enerbil_main for energy budget calculation,
566!! - call hydrol_main for hydrologic processes calculation,
567!! - call condveg_main for surface conditions such as roughness, albedo, and emmisivity,
568!! - call thermosoil_main for soil thermodynamic calculation,
569!! - call sechiba_end to swap previous to new fields.
570!!
571!! RECENT CHANGE(S): None
572!!
573!! MAIN OUTPUT VARIABLE(S): Hydrological variables (:: coastalflow and :: riverflow),
574!! components of the energy budget (:: tsol_rad, :: vevapp, :: fluxsens,
575!! :: temp_sol_new and :: fluxlat), surface characteristics (:: z0_out, :: emis_out,
576!! :: tq_cdrag and :: albedo) and land use related CO2 fluxes (:: netco2flux and
577!! :: fco2_lu, :: fco2_wh, ::fco2_ha)           
578!!
579!! REFERENCE(S) :
580!!
581!! FLOWCHART    :
582!! \latexonly
583!! \includegraphics[scale = 0.5]{sechibamainflow.png}
584!! \endlatexonly
585!! \n
586!_ ================================================================================================================================
587
588  SUBROUTINE sechiba_main (kjit, kjpij, kjpindex, index, &
589       & ldrestart_read, ldrestart_write, &
590       & lalo, contfrac, neighbours, resolution,&
591       & zlev, u, v, qair, temp_air, epot_air, ccanopy, &
592       & tq_cdrag, petAcoef, peqAcoef, petBcoef, peqBcoef, &
593       & precip_rain, precip_snow, lwdown, swnet, swdown, coszang, pb, &
594       & vevapp, fluxsens, fluxlat, coastalflow, riverflow, &
595       & netco2flux, fco2_lu, fco2_wh, fco2_ha, &
596       & tsol_rad, temp_sol_new, qsurf_out, albedo, emis_out, z0m_out, z0h_out,&
597       & veget_out, lai_out, height_out, &
598       & rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC)
599
600!! 0.1 Input variables
601   
602    INTEGER(i_std), INTENT(in)                               :: kjit              !! Time step number (unitless)
603    INTEGER(i_std), INTENT(in)                               :: kjpij             !! Total size of the un-compressed grid
604                                                                                  !! (unitless)
605    INTEGER(i_std), INTENT(in)                               :: kjpindex          !! Domain size - terrestrial pixels only
606                                                                                  !! (unitless)
607    INTEGER(i_std),INTENT (in)                               :: rest_id           !! _Restart_ file identifier (unitless)
608    INTEGER(i_std),INTENT (in)                               :: hist_id           !! _History_ file identifier (unitless)
609    INTEGER(i_std),INTENT (in)                               :: hist2_id          !! _History_ file 2 identifier (unitless)
610    INTEGER(i_std),INTENT (in)                               :: rest_id_stom      !! STOMATE's _Restart_ file identifier
611                                                                                  !! (unitless)
612    INTEGER(i_std),INTENT (in)                               :: hist_id_stom      !! STOMATE's _History_ file identifier
613                                                                                  !! (unitless)
614    INTEGER(i_std),INTENT(in)                                :: hist_id_stom_IPCC !! STOMATE's IPCC _history_ file file
615                                                                                  !! identifier (unitless)
616    LOGICAL, INTENT(in)                                      :: ldrestart_read    !! Logical for _restart_ file to read
617                                                                                  !! (true/false)
618    LOGICAL, INTENT(in)                                      :: ldrestart_write   !! Logical for _restart_ file to write
619                                                                                  !! (true/false)
620    REAL(r_std),DIMENSION (:,:), INTENT (in)          :: lalo              !! Geographic coordinates (latitude,longitude)
621                                                                                  !! for grid cells (degrees)
622    REAL(r_std),DIMENSION (:), INTENT (in)            :: contfrac          !! Fraction of continent in the grid
623                                                                                  !! (unitless, 0-1)
624    INTEGER(i_std),DIMENSION (:), INTENT (in)         :: index             !! Indices of the pixels on the map.
625                                                                                  !! Sechiba uses a reduced grid excluding oceans
626                                                                                  !! ::index contains the indices of the
627                                                                                  !! terrestrial pixels only! (unitless)
628    INTEGER(i_std), DIMENSION(:,:), INTENT(in) :: neighbours        !! Neighboring grid points if land!(unitless)
629    REAL(r_std), DIMENSION (:,:), INTENT(in)          :: resolution        !! Size in x and y of the grid (m)
630   
631    REAL(r_std),DIMENSION (:), INTENT (in)            :: u                 !! Lowest level wind speed in direction u
632                                                                                  !! @tex $(m.s^{-1})$ @endtex
633    REAL(r_std),DIMENSION (:), INTENT (in)            :: v                 !! Lowest level wind speed in direction v
634                                                                                  !! @tex $(m.s^{-1})$ @endtex
635    REAL(r_std),DIMENSION (:), INTENT (in)            :: zlev              !! Height of first layer (m)
636    REAL(r_std),DIMENSION (:), INTENT (in)            :: qair              !! Lowest level specific humidity
637                                                                                  !! @tex $(kg kg^{-1})$ @endtex
638    REAL(r_std),DIMENSION (:), INTENT (in)            :: precip_rain       !! Rain precipitation
639                                                                                  !! @tex $(kg m^{-2})$ @endtex
640    REAL(r_std),DIMENSION (:), INTENT (in)            :: precip_snow       !! Snow precipitation
641                                                                                  !! @tex $(kg m^{-2})$ @endtex
642    REAL(r_std),DIMENSION (:), INTENT (in)            :: lwdown            !! Down-welling long-wave flux
643                                                                                  !! @tex $(W m^{-2})$ @endtex
644    REAL(r_std),DIMENSION (:), INTENT (in)            :: coszang           !! Cosine of the solar zenith angle (unitless)
645    REAL(r_std),DIMENSION (:), INTENT (in)            :: swnet             !! Net surface short-wave flux
646                                                                                  !! @tex $(W m^{-2})$ @endtex
647    REAL(r_std),DIMENSION (:), INTENT (in)            :: swdown            !! Down-welling surface short-wave flux
648                                                                                  !! @tex $(W m^{-2})$ @endtex
649    REAL(r_std),DIMENSION (:), INTENT (in)            :: temp_air          !! Air temperature (K)
650    REAL(r_std),DIMENSION (:), INTENT (in)            :: epot_air          !! Air potential energy (??J)
651    REAL(r_std),DIMENSION (:), INTENT (in)            :: ccanopy           !! CO2 concentration in the canopy (ppm)
652    REAL(r_std),DIMENSION (:), INTENT (in)            :: petAcoef          !! Coefficients A for T from the Planetary
653                                                                                  !! Boundary Layer
654    REAL(r_std),DIMENSION (:), INTENT (in)            :: peqAcoef          !! Coefficients A for q from the Planetary
655                                                                                  !! Boundary Layer
656    REAL(r_std),DIMENSION (:), INTENT (in)            :: petBcoef          !! Coefficients B for T from the Planetary
657                                                                                  !! Boundary Layer
658    REAL(r_std),DIMENSION (:), INTENT (in)            :: peqBcoef          !! Coefficients B for q from the Planetary
659                                                                                  !! Boundary Layer
660    REAL(r_std),DIMENSION (:), INTENT (in)            :: pb                !! Surface pressure (hPa)
661
662
663!! 0.2 Output variables
664
665    REAL(r_std),DIMENSION (:), INTENT (out)           :: coastalflow       !! Outflow on coastal points by small basins.
666                                                                                  !! This is the water which flows in a disperse
667                                                                                  !! way into the ocean
668                                                                                  !! @tex $(kg dt_routing^{-1})$ @endtex
669    REAL(r_std),DIMENSION (:), INTENT (out)           :: riverflow         !! Outflow of the major rivers.
670                                                                                  !! The flux will be located on the continental
671                                                                                  !! grid but this should be a coastal point 
672                                                                                  !! @tex $(kg dt_routing^{-1})$ @endtex
673    REAL(r_std),DIMENSION (:), INTENT (out)           :: tsol_rad          !! Radiative surface temperature
674                                                                                  !! @tex $(W m^{-2})$ @endtex
675    REAL(r_std),DIMENSION (:), INTENT (out)           :: vevapp            !! Total of evaporation
676                                                                                  !! @tex $(kg m^{-2} days^{-1})$ @endtex
677   
678    REAL(r_std),DIMENSION (:), INTENT (out)           :: qsurf_out         !! Surface specific humidity
679                                                                                  !! @tex $(kg kg^{-1})$ @endtex
680    REAL(r_std),DIMENSION (:), INTENT (out)           :: z0m_out           !! Surface roughness momentum (output diagnostic, m)
681    REAL(r_std),DIMENSION (:), INTENT (out)           :: z0h_out           !! Surface roughness heat (output diagnostic, m)
682    REAL(r_std),DIMENSION (:,:), INTENT (out)         :: albedo            !! Surface albedo for visible and near-infrared
683                                                                                  !! (unitless, 0-1)
684    REAL(r_std),DIMENSION (:), INTENT (out)           :: fluxsens          !! Sensible heat flux
685                                                                                  !! @tex $(W m^{-2})$ @endtex
686    REAL(r_std),DIMENSION (:), INTENT (out)           :: fluxlat           !! Latent heat flux
687                                                                                  !! @tex $(W m^{-2})$ @endtex
688    REAL(r_std),DIMENSION (:), INTENT (out)           :: emis_out          !! Emissivity (unitless)
689    REAL(r_std),DIMENSION (:), INTENT (out)           :: netco2flux        !! Sum CO2 flux over PFTs
690                                                                                  !! (gC/m2/dt_sechiba)
691    REAL(r_std),DIMENSION (:), INTENT (out)           :: fco2_lu           !! Land Cover Change CO2 flux (gC/m2/one_day)
692    REAL(r_std),DIMENSION (:), INTENT (out)           :: fco2_wh           !! Wood harvest CO2 flux (gC/m2/one_day)
693    REAL(r_std),DIMENSION (:), INTENT (out)           :: fco2_ha           !! Crop harvest CO2 flux (gC/m2/one_day)
694    REAL(r_std),DIMENSION (:,:), INTENT (out)       :: veget_out         !! Fraction of vegetation type (unitless, 0-1)
695    REAL(r_std),DIMENSION (:,:), INTENT (out)       :: lai_out           !! Leaf area index (m^2 m^{-2})
696    REAL(r_std),DIMENSION (:,:), INTENT (out)       :: height_out        !! Vegetation Height (m)
697    REAL(r_std),DIMENSION (:), INTENT (out)           :: temp_sol_new      !! New ground temperature (K)
698
699!! 0.3 Modified
700
701    REAL(r_std),DIMENSION (:), INTENT (inout)         :: tq_cdrag          !! Surface drag coefficient  (-)
702
703!! 0.4 local variables
704
705    INTEGER(i_std)                                           :: ji, jv            !! Index (unitless)
706    REAL(r_std), DIMENSION(kjpindex)                         :: histvar           !! Computations for history files (unitless)
707    REAL(r_std), DIMENSION(kjpindex,nlut)                    :: histvar2          !! Computations for history files (unitless)
708    CHARACTER(LEN=80)                                        :: var_name          !! To store variables names for I/O (unitless)
709    REAL(r_std), DIMENSION(kjpindex)                         :: sum_treefrac      !! Total fraction occupied by trees (0-1, uniless)
710    REAL(r_std), DIMENSION(kjpindex)                         :: sum_grassfracC3   !! Total fraction occupied by C3 grasses (0-1, unitless)
711    REAL(r_std), DIMENSION(kjpindex)                         :: sum_grassfracC4   !! Total fraction occupied by C4 grasses (0-1, unitless)
712    REAL(r_std), DIMENSION(kjpindex)                         :: sum_cropfracC3    !! Total fraction occupied by C3 crops (0-1, unitess)
713    REAL(r_std), DIMENSION(kjpindex)                         :: sum_cropfracC4    !! Total fraction occupied by C4 crops (0-1, unitess)
714    REAL(r_std), DIMENSION(kjpindex)                         :: sum_treeFracNdlEvg!! Total fraction occupied by treeFracNdlEvg (0-1, unitess)
715    REAL(r_std), DIMENSION(kjpindex)                         :: sum_treeFracBdlEvg!! Total fraction occupied by treeFracBdlEvg (0-1, unitess)
716    REAL(r_std), DIMENSION(kjpindex)                         :: sum_treeFracNdlDcd!! Total fraction occupied by treeFracNdlDcd (0-1, unitess)
717    REAL(r_std), DIMENSION(kjpindex)                         :: sum_treeFracBdlDcd!! Total fraction occupied by treeFracBdlDcd (0-1, unitess)
718    REAL(r_std),DIMENSION (kjpindex,nvm)                     :: tmc_pft           !! Total soil water per PFT (mm/m2)
719    REAL(r_std),DIMENSION (kjpindex,nvm)                     :: drainage_pft     !! Drainage per PFT (mm/m2)   
720    REAL(r_std),DIMENSION (kjpindex,nvm)                     :: runoff_pft       !! Runoff per PFT (mm/m2)   
721    REAL(r_std),DIMENSION (kjpindex,nvm)                     :: swc_pft          !! Relative Soil water content [tmcr:tmcs] per pft (-)   
722    REAL(r_std), DIMENSION(kjpindex)                         :: grndflux          !! Net energy into soil (W/m2)
723    REAL(r_std), DIMENSION(kjpindex,nsnow)                   :: snowliq           !! Liquid water content (m)
724    REAL(r_std), DIMENSION(kjpindex)                         :: snow_age_diag     !! Only for diag, contains xios_default_val
725    REAL(r_std), DIMENSION(kjpindex,nnobio)                  :: snow_nobio_age_diag !! Only for diag, contains xios_default_val
726    REAL(r_std), DIMENSION(kjpindex)                         :: snowage_glob      !! Snow age on total area including snow on vegetated and bare soil and nobio area @tex ($d$) @endtex
727    REAL(r_std), DIMENSION(kjpindex,nlut)                    :: gpplut            !! GPP on landuse tile, only for diagnostics
728    REAL(r_std), DIMENSION(kjpindex,nslm,nvm)                :: mc_layh_pft       !! As mc_layh per pft
729    REAL(r_std), DIMENSION(kjpindex,nslm,nvm)                :: mcl_layh_pft      !! As mcl_layh per pft
730    REAL(r_std), DIMENSION(kjpindex,nslm,nvm)                :: soilmoist_pft     !! As soilmoist per pft
731    REAL(r_std), DIMENSION(nscm)                             :: mcs_hydrol        !! Saturated volumetric water content output to be used in stomate_soilcarbon
732    REAL(r_std), DIMENSION(nscm)                             :: mcfc_hydrol       !! Volumetric water content at field capacity output to be used in stomate_soilcarbon
733
734
735!_ ================================================================================================================================
736
737    IF (printlev_loc>=3) WRITE(numout,*) 'Start sechiba_main kjpindex =',kjpindex
738
739    !! 1. Initialize variables at each time step
740    CALL sechiba_var_init (kjpindex, rau, pb, temp_air) 
741
742    !! 2. Compute diffusion coefficients
743    CALL diffuco_main (kjit, kjpindex, index, indexveg, indexlai, u, v, &
744         & zlev, z0m, z0h, roughheight, temp_sol, temp_air, temp_growth, rau, tq_cdrag, qsurf, qair, pb ,  &
745         & evap_bare_lim, evap_bare_lim_ns, evapot, evapot_corr, snow, flood_frac, flood_res, &
746         & frac_nobio, snow_nobio, totfrac_nobio, &
747         & swnet, swdown, coszang, ccanopy, humrel, veget, veget_max, lai, qsintveg, qsintmax, assim_param, &
748         & vbeta, vbeta1, vbeta2, vbeta3, vbeta3pot, vbeta4, vbeta5, gsmean, rveget, rstruct, cimean, gpp, co2_to_bm, &
749         & lalo, neighbours, resolution, ptnlev1, precip_rain, frac_age, tot_bare_soil, frac_snow_veg, frac_snow_nobio, &
750         & hist_id, hist2_id)
751   
752    !! 3. Compute energy balance
753    CALL enerbil_main (kjit, kjpindex, &
754         & index, indexveg, zlev, lwdown, swnet, epot_air, temp_air, u, v, petAcoef, petBcoef, &
755         & qair, peqAcoef, peqBcoef, pb, rau, vbeta, vbeta1, vbeta2, vbeta3, vbeta3pot, vbeta4, vbeta5, &
756         & emis, soilflx, soilcap, tq_cdrag, humrel, fluxsens, fluxlat, &
757         & vevapp, transpir, transpot, vevapnu, vevapwet, vevapsno, vevapflo, temp_sol, tsol_rad, &
758         & temp_sol_new, qsurf, evapot, evapot_corr, rest_id, hist_id, hist2_id, &
759         & precip_rain,  pgflux, snowdz, temp_sol_add)
760
761 
762    !! 4. Compute hydrology
763    !! 4.1 Water balance from CWRR module (11 soil layers)
764    CALL hydrol_main (kjit, kjpindex, &
765         & index, indexveg, indexsoil, indexlayer, indexnslm, &
766         & temp_sol_new, floodout, runoff, drainage, frac_nobio, totfrac_nobio, vevapwet, veget, veget_max, njsc, &
767         & qsintmax, qsintveg, vevapnu, vevapsno, vevapflo, snow, snow_age, snow_nobio, snow_nobio_age,  &
768         & tot_melt, transpir, precip_rain, precip_snow, returnflow, reinfiltration, irrigation, &
769         & humrel, vegstress, drysoil_frac, evapot, evapot_corr, evap_bare_lim, evap_bare_lim_ns, flood_frac, flood_res, &
770         & shumdiag,shumdiag_perma, k_litt, litterhumdiag, soilcap, soiltile, fraclut, reinf_slope,&
771         & rest_id, hist_id, hist2_id,&
772         & contfrac, stempdiag, &
773         & temp_air, pb, u, v, tq_cdrag, swnet, pgflux, &
774         & snowrho,snowtemp,snowgrain,snowdz,snowheat,snowliq, &
775         & grndflux,gtemp,tot_bare_soil, &
776         & lambda_snow,cgrnd_snow,dgrnd_snow,frac_snow_veg,temp_sol_add, &
777         & mc_layh, mcl_layh, tmc_pft, drainage_pft, runoff_pft, swc_pft, soilmoist, mc_layh_s, mcl_layh_s, soilmoist_s, &
778         & mcs_hydrol, mcfc_hydrol) 
779
780         
781    !! 6. Compute surface variables (emissivity, albedo and roughness)
782    CALL condveg_main (kjit, kjpindex, index, rest_id, hist_id, hist2_id, &
783         lalo, neighbours, resolution, contfrac, veget, veget_max, frac_nobio, totfrac_nobio, &
784         zlev, snow, snow_age, snow_nobio, snow_nobio_age, &
785         drysoil_frac, height, snowdz, snowrho, tot_bare_soil, &
786         temp_air, pb, u, v, lai, &
787         emis, albedo, z0m, z0h, roughheight, &
788         frac_snow_veg, frac_snow_nobio)
789
790
791    !! 7. Compute soil thermodynamics
792    DO jv = 1,nvm
793      mc_layh_pft(:,:,jv) = mc_layh_s(:,:,pref_soil_veg(jv))
794      mcl_layh_pft(:,:,jv) = mcl_layh_s(:,:,pref_soil_veg(jv))
795      soilmoist_pft(:,:,jv) = soilmoist_s(:,:,pref_soil_veg(jv))
796    END DO
797    CALL thermosoil_main (kjit, kjpindex, &
798         index, indexgrnd, &
799         temp_sol_new, snow, soilcap, soilflx, &
800         shumdiag_perma, stempdiag, ptnlev1, rest_id, hist_id, hist2_id, &
801         snowdz,snowrho,snowtemp,gtemp,pb,&
802         mc_layh, mcl_layh, soilmoist, &
803         mc_layh_pft, mcl_layh_pft,soilmoist_pft, njsc, &
804         depth_organic_soil, heat_Zimov, tdeep, hsdeep,&
805         som_total, veget_max, &
806         frac_snow_veg,frac_snow_nobio,totfrac_nobio,temp_sol_add, &
807         lambda_snow, cgrnd_snow, dgrnd_snow)
808
809
810    !! 8. Compute river routing
811    IF ( river_routing .AND. nbp_glo .GT. 1) THEN
812       !! 8.1 River routing
813       CALL routing_main (kjit, kjpindex, index, &
814            & lalo, neighbours, resolution, contfrac, totfrac_nobio, veget_max, floodout, runoff, &
815            & drainage, transpot, precip_rain, humrel, k_litt, flood_frac, flood_res, &
816            & stempdiag, reinf_slope, returnflow, reinfiltration, irrigation, riverflow, coastalflow, rest_id, hist_id, hist2_id)
817    ELSE
818       !! 8.2 No routing, set variables to zero
819       riverflow(:) = zero
820       coastalflow(:) = zero
821       returnflow(:) = zero
822       reinfiltration(:) = zero
823       irrigation(:) = zero
824       flood_frac(:) = zero
825       flood_res(:) = zero
826
827       CALL xios_orchidee_send_field("coastalflow",coastalflow/dt_sechiba)
828       CALL xios_orchidee_send_field("riverflow",riverflow/dt_sechiba)
829    ENDIF
830
831    !! 9. Compute slow processes (i.e. 'daily' and annual time step)
832    CALL slowproc_main (kjit, kjpij, kjpindex, njsc, &
833         index, indexveg, lalo, neighbours, resolution, contfrac, soiltile, fraclut, nwdFraclut, &
834         temp_air, temp_sol, stempdiag, &
835         vegstress, shumdiag, litterhumdiag, precip_rain, precip_snow, pb, gpp, &
836         tmc_pft, drainage_pft, runoff_pft, swc_pft, deadleaf_cover, &
837         assim_param, &
838         lai, frac_age, height, veget, frac_nobio, veget_max, totfrac_nobio, qsintmax, &
839         rest_id, hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, &
840         co2_flux, fco2_lu, fco2_wh, fco2_ha, co2_to_bm, temp_growth, tot_bare_soil, &
841         tdeep, hsdeep, snow, heat_Zimov, &
842         sfluxCH4_deep, sfluxCO2_deep, &
843         som_total,snowdz,snowrho, altmax, depth_organic_soil, mcs_hydrol, mcfc_hydrol)
844
845
846    !! 9.2 Compute global CO2 flux
847    netco2flux(:) = zero
848    DO jv = 2,nvm
849      netco2flux(:) = netco2flux(:) + co2_flux(:,jv)*(1-totfrac_nobio)
850    ENDDO
851 
852    !! 10. Update the temperature (temp_sol) with newly computed values
853    CALL sechiba_end (kjpindex, temp_sol_new, temp_sol)
854
855   
856    !! 11. Write internal variables to output fields
857    z0m_out(:) = z0m(:)
858    z0h_out(:) = z0h(:)
859    emis_out(:) = emis(:)
860    qsurf_out(:) = qsurf(:)
861    veget_out(:,:)  = veget(:,:)
862    lai_out(:,:)    = lai(:,:)
863    height_out(:,:) = height(:,:)
864 
865    !! 12. Write global variables to history files
866    sum_treefrac(:) = zero
867    sum_grassfracC3(:) = zero
868    sum_grassfracC4(:) = zero
869    sum_cropfracC3(:) = zero
870    sum_cropfracC4(:) = zero
871    sum_treeFracNdlEvg(:) = zero
872    sum_treeFracBdlEvg(:) = zero
873    sum_treeFracNdlDcd(:) = zero
874    sum_treeFracBdlDcd(:) = zero
875
876    DO jv = 2, nvm 
877       IF (is_tree(jv) .AND. natural(jv)) THEN
878          sum_treefrac(:) = sum_treefrac(:) + veget_max(:,jv)
879       ELSE IF ((.NOT. is_tree(jv))  .AND. natural(jv)) THEN
880          ! Grass
881          IF (is_c4(jv)) THEN
882             sum_grassfracC4(:) = sum_grassfracC4(:) + veget_max(:,jv)
883          ELSE
884             sum_grassfracC3(:) = sum_grassfracC3(:) + veget_max(:,jv)
885          END IF
886       ELSE 
887          ! Crop and trees not natural
888          IF (is_c4(jv)) THEN
889             sum_cropfracC4(:) = sum_cropfracC4(:) + veget_max(:,jv)
890          ELSE
891             sum_cropfracC3(:) = sum_cropfracC3(:) + veget_max(:,jv)
892          END IF
893       ENDIF
894
895       IF (is_tree(jv)) THEN
896          IF (is_evergreen(jv)) THEN
897             IF (is_needleleaf(jv)) THEN
898                ! Fraction for needleleaf evergreen trees (treeFracNdlEvg)
899                sum_treeFracNdlEvg(:) = sum_treeFracNdlEvg(:) + veget_max(:,jv)
900             ELSE
901                ! Fraction for broadleaf evergreen trees (treeFracBdlEvg)
902                sum_treeFracBdlEvg(:) = sum_treeFracBdlEvg(:) + veget_max(:,jv)
903             END IF
904          ELSE IF (is_deciduous(jv)) THEN
905             IF (is_needleleaf(jv)) THEN
906                ! Fraction for needleleaf deciduous trees (treeFracNdlDcd)
907                sum_treeFracNdlDcd(:) = sum_treeFracNdlDcd(:) + veget_max(:,jv)
908             ELSE 
909                ! Fraction for broadleafs deciduous trees (treeFracBdlDcd)
910                sum_treeFracBdlDcd(:) = sum_treeFracBdlDcd(:) + veget_max(:,jv)
911             END IF
912          END IF
913       END IF
914    ENDDO         
915
916    histvar(:)=zero
917    DO jv = 2, nvm
918       IF (is_deciduous(jv)) THEN
919          histvar(:) = histvar(:) + veget_max(:,jv)*100*contfrac
920       ENDIF
921    ENDDO
922    CALL xios_orchidee_send_field("treeFracPrimDec",histvar)
923
924    histvar(:)=zero
925    DO jv = 2, nvm
926       IF (is_evergreen(jv)) THEN
927          histvar(:) = histvar(:) + veget_max(:,jv)*100*contfrac
928       ENDIF
929    ENDDO
930    CALL xios_orchidee_send_field("treeFracPrimEver",histvar)
931
932    histvar(:)=zero
933    DO jv = 2, nvm
934       IF ( .NOT.(is_c4(jv)) ) THEN
935          histvar(:) = histvar(:) + veget_max(:,jv)*100*contfrac
936       ENDIF
937    ENDDO
938    CALL xios_orchidee_send_field("c3PftFrac",histvar)
939
940    histvar(:)=zero
941    DO jv = 2, nvm
942       IF ( is_c4(jv) ) THEN
943          histvar(:) = histvar(:) + veget_max(:,jv)*100*contfrac
944       ENDIF
945    ENDDO
946    CALL xios_orchidee_send_field("c4PftFrac",histvar)
947
948    CALL xios_orchidee_send_field("temp_sol_new",temp_sol_new)
949    CALL xios_orchidee_send_field("fluxsens",fluxsens)
950    CALL xios_orchidee_send_field("fluxlat",fluxlat)
951
952
953    ! Add XIOS default value where no snow
954    DO ji=1,kjpindex 
955       IF (snow(ji) .GT. zero) THEN
956          snow_age_diag(ji) = snow_age(ji)
957          snow_nobio_age_diag(ji,:) = snow_nobio_age(ji,:)
958       
959          snowage_glob(ji) = snow_age(ji)*frac_snow_veg(ji)*(1-totfrac_nobio(ji)) + &
960               SUM(snow_nobio_age(ji,:)*frac_snow_nobio(ji,:)*frac_nobio(ji,:))
961          IF (snowage_glob(ji) .NE. 0) snowage_glob(ji) = snowage_glob(ji) / &
962               (frac_snow_veg(ji)*(1-totfrac_nobio(ji)) + SUM(frac_snow_nobio(ji,:)*frac_nobio(ji,:)))
963       ELSE
964          snow_age_diag(ji) = xios_default_val
965          snow_nobio_age_diag(ji,:) = xios_default_val
966          snowage_glob(ji) = xios_default_val
967       END IF
968    END DO
969   
970    CALL xios_orchidee_send_field("snow",snow)
971    CALL xios_orchidee_send_field("snowage",snow_age_diag)
972    CALL xios_orchidee_send_field("snownobio",snow_nobio)
973    CALL xios_orchidee_send_field("snownobioage",snow_nobio_age_diag)
974    CALL xios_orchidee_send_field("snowage_glob",snowage_glob)
975
976    CALL xios_orchidee_send_field("frac_snow", SUM(frac_snow_nobio,2)*totfrac_nobio+frac_snow_veg*(1-totfrac_nobio))
977    CALL xios_orchidee_send_field("frac_snow_veg", frac_snow_veg)
978    CALL xios_orchidee_send_field("frac_snow_nobio", frac_snow_nobio)
979    CALL xios_orchidee_send_field("reinf_slope",reinf_slope)
980    CALL xios_orchidee_send_field("njsc",REAL(njsc, r_std))
981    CALL xios_orchidee_send_field("vegetfrac",veget)
982    CALL xios_orchidee_send_field("maxvegetfrac",veget_max)
983    CALL xios_orchidee_send_field("nobiofrac",frac_nobio)
984    CALL xios_orchidee_send_field("soiltile",soiltile)
985    CALL xios_orchidee_send_field("rstruct",rstruct)
986    CALL xios_orchidee_send_field("gpp",gpp/dt_sechiba)
987    CALL xios_orchidee_send_field("gpp_ipcc2",SUM(gpp,dim=2)/dt_sechiba)
988
989    histvar(:)=zero
990    DO jv = 2, nvm
991       IF ( .NOT. is_tree(jv) .AND. natural(jv) ) THEN
992          histvar(:) = histvar(:) + gpp(:,jv)
993       ENDIF
994    ENDDO
995    CALL xios_orchidee_send_field("gppgrass",histvar/dt_sechiba)
996
997    histvar(:)=zero
998    DO jv = 2, nvm
999       IF ( (.NOT. is_tree(jv)) .AND. (.NOT. natural(jv)) ) THEN
1000          histvar(:) = histvar(:) + gpp(:,jv)
1001       ENDIF
1002    ENDDO
1003    CALL xios_orchidee_send_field("gppcrop",histvar/dt_sechiba)
1004
1005    histvar(:)=zero
1006    DO jv = 2, nvm
1007       IF ( is_tree(jv) ) THEN
1008          histvar(:) = histvar(:) + gpp(:,jv)
1009       ENDIF
1010    ENDDO
1011    CALL xios_orchidee_send_field("gpptree",histvar/dt_sechiba)
1012    CALL xios_orchidee_send_field("nee",co2_flux/1.e3/one_day)
1013    CALL xios_orchidee_send_field("drysoil_frac",drysoil_frac)
1014    CALL xios_orchidee_send_field("vevapflo",vevapflo/dt_sechiba)
1015    CALL xios_orchidee_send_field("k_litt",k_litt)
1016    CALL xios_orchidee_send_field("beta",vbeta)
1017    CALL xios_orchidee_send_field("vbeta1",vbeta1)
1018    CALL xios_orchidee_send_field("vbeta2",vbeta2)
1019    CALL xios_orchidee_send_field("vbeta3",vbeta3)
1020    CALL xios_orchidee_send_field("vbeta4",vbeta4)
1021    CALL xios_orchidee_send_field("vbeta5",vbeta5)
1022    CALL xios_orchidee_send_field("gsmean",gsmean)
1023    CALL xios_orchidee_send_field("cimean",cimean)
1024    CALL xios_orchidee_send_field("rveget",rveget)
1025 
1026    histvar(:)=SUM(vevapwet(:,:),dim=2)
1027    CALL xios_orchidee_send_field("evspsblveg",histvar/dt_sechiba)
1028    histvar(:)= vevapnu(:)+vevapsno(:)
1029    CALL xios_orchidee_send_field("evspsblsoi",histvar/dt_sechiba)
1030    histvar(:)=SUM(transpir(:,:),dim=2)
1031    CALL xios_orchidee_send_field("tran",histvar/dt_sechiba)
1032
1033    ! For CMIP6 data request: the following fractions are fractions of the total grid-cell,
1034    ! which explains the multiplication by contfrac
1035    CALL xios_orchidee_send_field("treeFrac",sum_treefrac(:)*100*contfrac(:))
1036    CALL xios_orchidee_send_field("grassFracC3",sum_grassfracC3(:)*100*contfrac(:))
1037    CALL xios_orchidee_send_field("grassFracC4",sum_grassfracC4(:)*100*contfrac(:))
1038    CALL xios_orchidee_send_field("cropFracC3",sum_cropfracC3(:)*100*contfrac(:))
1039    CALL xios_orchidee_send_field("cropFracC4",sum_cropfracC4(:)*100*contfrac(:))
1040    CALL xios_orchidee_send_field("treeFracNdlEvg",sum_treeFracNdlEvg(:)*100*contfrac(:))
1041    CALL xios_orchidee_send_field("treeFracBdlEvg",sum_treeFracBdlEvg(:)*100*contfrac(:))
1042    CALL xios_orchidee_send_field("treeFracNdlDcd",sum_treeFracNdlDcd(:)*100*contfrac(:))
1043    CALL xios_orchidee_send_field("treeFracBdlDcd",sum_treeFracBdlDcd(:)*100*contfrac(:))
1044
1045    histvar(:)=veget_max(:,1)*100*contfrac(:)
1046    CALL xios_orchidee_send_field("baresoilFrac",histvar)
1047    histvar(:)=SUM(frac_nobio(:,1:nnobio),dim=2)*100*contfrac(:)
1048    CALL xios_orchidee_send_field("residualFrac",histvar)
1049
1050    ! For CMIP6 data request: cnc = canopy cover fraction over land area
1051    histvar(:)=zero
1052    DO jv=2,nvm
1053       histvar(:) = histvar(:) + veget_max(:,jv)*100
1054    END DO
1055    CALL xios_orchidee_send_field("cnc",histvar)
1056   
1057    CALL xios_orchidee_send_field("tsol_rad",tsol_rad-273.15)
1058    CALL xios_orchidee_send_field("qsurf",qsurf)
1059    CALL xios_orchidee_send_field("emis",emis)
1060    CALL xios_orchidee_send_field("z0m",z0m)
1061    CALL xios_orchidee_send_field("z0h",z0h)
1062    CALL xios_orchidee_send_field("roughheight",roughheight)
1063    CALL xios_orchidee_send_field("lai",lai)
1064    histvar(:)=zero   
1065    DO ji = 1, kjpindex
1066       IF (SUM(veget_max(ji,:)) > zero) THEN
1067         DO jv=2,nvm
1068            histvar(ji) = histvar(ji) + veget_max(ji,jv)*lai(ji,jv)/SUM(veget_max(ji,:))
1069         END DO
1070       END IF
1071    END DO
1072    CALL xios_orchidee_send_field("LAImean",histvar)
1073    CALL xios_orchidee_send_field("vevapsno",vevapsno/dt_sechiba)
1074    CALL xios_orchidee_send_field("vevapp",vevapp/dt_sechiba)
1075    CALL xios_orchidee_send_field("vevapnu",vevapnu/dt_sechiba)
1076    CALL xios_orchidee_send_field("transpir",transpir*one_day/dt_sechiba)
1077    CALL xios_orchidee_send_field("inter",vevapwet*one_day/dt_sechiba)
1078    histvar(:)=zero
1079    DO jv=1,nvm
1080      histvar(:) = histvar(:) + vevapwet(:,jv)
1081    ENDDO
1082    CALL xios_orchidee_send_field("ECanop",histvar/dt_sechiba)
1083    histvar(:)=zero
1084    DO jv=1,nvm
1085      histvar(:) = histvar(:) + transpir(:,jv)
1086    ENDDO
1087    CALL xios_orchidee_send_field("TVeg",histvar/dt_sechiba)
1088    CALL xios_orchidee_send_field("ACond",tq_cdrag)
1089
1090
1091    !! Calculate diagnostic variables on Landuse tiles for LUMIP/CMIP6
1092
1093    ! Calculate fraction of landuse tiles related to the whole grid cell
1094    DO jv=1,nlut
1095       histvar2(:,jv) = fraclut(:,jv) * contfrac(:)
1096    END DO
1097    CALL xios_orchidee_send_field("fraclut",histvar2)
1098
1099    CALL xios_orchidee_send_field("nwdFraclut",nwdFraclut(:,:))
1100   
1101    ! Calculate GPP on landuse tiles
1102    ! val_exp is used as missing value where the values are not known i.e. where the tile is not represented
1103    ! or for pasture (id_pst) or urban land (id_urb).
1104    gpplut(:,:)=0
1105    DO jv=1,nvm
1106       IF (natural(jv)) THEN
1107          gpplut(:,id_psl) = gpplut(:,id_psl) + gpp(:,jv)
1108       ELSE
1109          gpplut(:,id_crp) = gpplut(:,id_crp) + gpp(:,jv)
1110       ENDIF
1111    END DO
1112
1113    ! Transform from gC/m2/s into kgC/m2/s
1114    WHERE (fraclut(:,id_psl)>min_sechiba)
1115       gpplut(:,id_psl) = gpplut(:,id_psl)/fraclut(:,id_psl)/1000
1116    ELSEWHERE
1117       gpplut(:,id_psl) = xios_default_val
1118    END WHERE
1119    WHERE (fraclut(:,id_crp)>min_sechiba)
1120       gpplut(:,id_crp) = gpplut(:,id_crp)/fraclut(:,id_crp)/1000
1121    ELSEWHERE
1122       gpplut(:,id_crp) = xios_default_val
1123    END WHERE
1124    gpplut(:,id_pst) = xios_default_val
1125    gpplut(:,id_urb) = xios_default_val
1126
1127    CALL xios_orchidee_send_field("gpplut",gpplut)
1128
1129
1130    IF ( .NOT. almaoutput ) THEN
1131       ! Write history file in IPSL-format
1132       CALL histwrite_p(hist_id, 'beta', kjit, vbeta, kjpindex, index)
1133       CALL histwrite_p(hist_id, 'z0m', kjit, z0m, kjpindex, index)
1134       CALL histwrite_p(hist_id, 'z0h', kjit, z0h, kjpindex, index)
1135       CALL histwrite_p(hist_id, 'roughheight', kjit, roughheight, kjpindex, index)
1136       CALL histwrite_p(hist_id, 'vegetfrac', kjit, veget, kjpindex*nvm, indexveg)
1137       CALL histwrite_p(hist_id, 'maxvegetfrac', kjit, veget_max, kjpindex*nvm, indexveg)
1138       CALL histwrite_p(hist_id, 'nobiofrac', kjit, frac_nobio, kjpindex*nnobio, indexnobio)
1139       CALL histwrite_p(hist_id, 'lai', kjit, lai, kjpindex*nvm, indexveg)
1140       histvar(:)=zero   
1141       DO ji = 1, kjpindex
1142          IF (SUM(veget_max(ji,:)) > zero) THEN
1143            DO jv=2,nvm
1144               histvar(ji) = histvar(ji) + veget_max(ji,jv)*lai(ji,jv)/SUM(veget_max(ji,:))
1145            END DO
1146          END IF
1147       END DO
1148
1149       CALL histwrite_p(hist_id, 'LAImean', kjit, histvar, kjpindex, index)
1150
1151       CALL histwrite_p(hist_id, 'subli', kjit, vevapsno, kjpindex, index)
1152       CALL histwrite_p(hist_id, 'evapnu', kjit, vevapnu, kjpindex, index)
1153       CALL histwrite_p(hist_id, 'transpir', kjit, transpir, kjpindex*nvm, indexveg)
1154       CALL histwrite_p(hist_id, 'inter', kjit, vevapwet, kjpindex*nvm, indexveg)
1155       CALL histwrite_p(hist_id, 'vbeta1', kjit, vbeta1, kjpindex, index)
1156       CALL histwrite_p(hist_id, 'vbeta2', kjit, vbeta2, kjpindex*nvm, indexveg)
1157       CALL histwrite_p(hist_id, 'vbeta3', kjit, vbeta3, kjpindex*nvm, indexveg)
1158       CALL histwrite_p(hist_id, 'vbeta4', kjit, vbeta4, kjpindex, index)   
1159       CALL histwrite_p(hist_id, 'vbeta5', kjit, vbeta5, kjpindex, index)   
1160       CALL histwrite_p(hist_id, 'drysoil_frac', kjit, drysoil_frac, kjpindex, index)
1161       CALL histwrite_p(hist_id, 'rveget', kjit, rveget, kjpindex*nvm, indexveg)
1162       CALL histwrite_p(hist_id, 'rstruct', kjit, rstruct, kjpindex*nvm, indexveg)
1163       CALL histwrite_p(hist_id, 'snow', kjit, snow, kjpindex, index)
1164       CALL histwrite_p(hist_id, 'snowage', kjit, snow_age, kjpindex, index)
1165       CALL histwrite_p(hist_id, 'snownobio', kjit, snow_nobio, kjpindex*nnobio, indexnobio)
1166       CALL histwrite_p(hist_id, 'snownobioage', kjit, snow_nobio_age, kjpindex*nnobio, indexnobio)
1167
1168       CALL histwrite_p(hist_id, 'grndflux', kjit, grndflux, kjpindex,index)
1169       CALL histwrite_p(hist_id, 'snowtemp',kjit,snowtemp,kjpindex*nsnow,indexsnow)
1170       CALL histwrite_p(hist_id, 'snowliq', kjit,snowliq,kjpindex*nsnow,indexsnow)
1171       CALL histwrite_p(hist_id, 'snowdz', kjit,snowdz,kjpindex*nsnow,indexsnow)
1172       CALL histwrite_p(hist_id, 'snowrho', kjit,snowrho,kjpindex*nsnow,indexsnow)
1173       CALL histwrite_p(hist_id, 'snowgrain',kjit,snowgrain,kjpindex*nsnow,indexsnow)
1174       CALL histwrite_p(hist_id, 'snowheat',kjit,snowheat,kjpindex*nsnow,indexsnow)
1175       
1176       CALL histwrite_p(hist_id,'frac_snow_veg',kjit,frac_snow_veg,kjpindex,index)
1177       CALL histwrite_p(hist_id, 'frac_snow_nobio', kjit,frac_snow_nobio,kjpindex*nnobio, indexnobio)
1178       CALL histwrite_p(hist_id, 'pgflux',kjit,pgflux,kjpindex,index)
1179       CALL histwrite_p(hist_id, 'soiltile',  kjit, soiltile, kjpindex*nstm, indexsoil)
1180
1181       CALL histwrite_p(hist_id, 'soilindex',  kjit, REAL(njsc, r_std), kjpindex, index)
1182       CALL histwrite_p(hist_id, 'reinf_slope',  kjit, reinf_slope, kjpindex, index)
1183       CALL histwrite_p(hist_id, 'k_litt', kjit, k_litt, kjpindex, index)
1184       
1185       IF ( do_floodplains ) THEN
1186          CALL histwrite_p(hist_id, 'evapflo', kjit, vevapflo, kjpindex, index)
1187          CALL histwrite_p(hist_id, 'flood_frac', kjit, flood_frac, kjpindex, index)
1188       ENDIF
1189       
1190       CALL histwrite_p(hist_id, 'gsmean', kjit, gsmean, kjpindex*nvm, indexveg)   
1191       CALL histwrite_p(hist_id, 'gpp', kjit, gpp, kjpindex*nvm, indexveg)
1192       CALL histwrite_p(hist_id, 'cimean', kjit, cimean, kjpindex*nvm, indexveg)   
1193       
1194       IF ( ok_stomate ) THEN
1195          CALL histwrite_p(hist_id, 'nee', kjit, co2_flux/1.e3/one_day, kjpindex*nvm, indexveg)   
1196       ENDIF
1197
1198       histvar(:)=SUM(vevapwet(:,:),dim=2)
1199       CALL histwrite_p(hist_id, 'evspsblveg', kjit, histvar, kjpindex, index)
1200
1201       histvar(:)= vevapnu(:)+vevapsno(:)
1202       CALL histwrite_p(hist_id, 'evspsblsoi', kjit, histvar, kjpindex, index)
1203
1204       histvar(:)=SUM(transpir(:,:),dim=2)
1205       CALL histwrite_p(hist_id, 'tran', kjit, histvar, kjpindex, index)
1206
1207       histvar(:)= sum_treefrac(:)*100*contfrac(:)
1208       CALL histwrite_p(hist_id, 'treeFrac', kjit, histvar, kjpindex, index) 
1209
1210       histvar(:)= (sum_grassfracC3(:)+sum_grassfracC4(:))*100*contfrac(:)
1211       CALL histwrite_p(hist_id, 'grassFrac', kjit, histvar, kjpindex, index) 
1212
1213       histvar(:)= (sum_cropfracC3(:)+sum_cropfracC4(:))*100*contfrac(:)
1214       CALL histwrite_p(hist_id, 'cropFrac', kjit, histvar, kjpindex, index)
1215
1216       histvar(:)=veget_max(:,1)*100*contfrac(:)
1217       CALL histwrite_p(hist_id, 'baresoilFrac', kjit, histvar, kjpindex, index)
1218
1219       histvar(:)=SUM(frac_nobio(:,1:nnobio),dim=2)*100*contfrac(:)
1220       CALL histwrite_p(hist_id, 'residualFrac', kjit, histvar, kjpindex, index)
1221    ELSE
1222       ! Write history file in ALMA format
1223       CALL histwrite_p(hist_id, 'vegetfrac', kjit, veget, kjpindex*nvm, indexveg)
1224       CALL histwrite_p(hist_id, 'maxvegetfrac', kjit, veget_max, kjpindex*nvm, indexveg)
1225       CALL histwrite_p(hist_id, 'nobiofrac', kjit, frac_nobio, kjpindex*nnobio, indexnobio)
1226       CALL histwrite_p(hist_id, 'lai', kjit, lai, kjpindex*nvm, indexveg)
1227       CALL histwrite_p(hist_id, 'ESoil', kjit, vevapnu, kjpindex, index)
1228       CALL histwrite_p(hist_id, 'EWater', kjit, vevapflo, kjpindex, index)
1229       CALL histwrite_p(hist_id, 'SWE', kjit, snow, kjpindex, index)
1230       histvar(:)=zero
1231       DO jv=1,nvm
1232          histvar(:) = histvar(:) + transpir(:,jv)
1233       ENDDO
1234       CALL histwrite_p(hist_id, 'TVeg', kjit, histvar, kjpindex, index)
1235       histvar(:)=zero
1236       DO jv=1,nvm
1237          histvar(:) = histvar(:) + vevapwet(:,jv)
1238       ENDDO
1239       CALL histwrite_p(hist_id, 'ECanop', kjit, histvar, kjpindex, index)
1240       CALL histwrite_p(hist_id, 'ACond', kjit, tq_cdrag, kjpindex, index)
1241       CALL histwrite_p(hist_id, 'SnowFrac', kjit, vbeta1, kjpindex, index)
1242       !
1243       CALL histwrite_p(hist_id, 'Z0m', kjit, z0m, kjpindex, index)
1244       CALL histwrite_p(hist_id, 'Z0h', kjit, z0h, kjpindex, index)
1245       CALL histwrite_p(hist_id, 'EffectHeight', kjit, roughheight, kjpindex, index)
1246       !
1247       IF ( do_floodplains ) THEN
1248          CALL histwrite_p(hist_id, 'Qflood', kjit, vevapflo, kjpindex, index)
1249          CALL histwrite_p(hist_id, 'FloodFrac', kjit, flood_frac, kjpindex, index)
1250       ENDIF
1251
1252       CALL histwrite_p(hist_id, 'gsmean', kjit, gsmean, kjpindex*nvm, indexveg)   
1253       CALL histwrite_p(hist_id, 'cimean', kjit, cimean, kjpindex*nvm, indexveg)   
1254       CALL histwrite_p(hist_id, 'GPP', kjit, gpp, kjpindex*nvm, indexveg)
1255       
1256       IF ( ok_stomate ) THEN
1257             CALL histwrite_p(hist_id, 'NEE', kjit, co2_flux, kjpindex*nvm, indexveg)   
1258       ENDIF
1259    ENDIF ! almaoutput
1260   
1261    !! 13. Write additional output file with higher frequency
1262    IF ( hist2_id > 0 ) THEN
1263       IF ( .NOT. almaoutput ) THEN
1264          ! Write history file in IPSL-format
1265          CALL histwrite_p(hist2_id, 'tsol_rad', kjit, tsol_rad, kjpindex, index)
1266          CALL histwrite_p(hist2_id, 'qsurf', kjit, qsurf, kjpindex, index)
1267          CALL histwrite_p(hist2_id, 'albedo', kjit, albedo, kjpindex*2, indexalb)
1268          CALL histwrite_p(hist2_id, 'emis', kjit, emis, kjpindex, index)
1269          CALL histwrite_p(hist2_id, 'beta', kjit, vbeta, kjpindex, index)
1270          CALL histwrite_p(hist2_id, 'z0m', kjit, z0m, kjpindex, index)
1271          CALL histwrite_p(hist2_id, 'z0h', kjit, z0h, kjpindex, index)
1272          CALL histwrite_p(hist2_id, 'roughheight', kjit, roughheight, kjpindex, index)
1273          CALL histwrite_p(hist2_id, 'vegetfrac', kjit, veget, kjpindex*nvm, indexveg)
1274          CALL histwrite_p(hist2_id, 'maxvegetfrac', kjit, veget_max, kjpindex*nvm, indexveg)
1275          CALL histwrite_p(hist2_id, 'nobiofrac', kjit, frac_nobio, kjpindex*nnobio, indexnobio)
1276          CALL histwrite_p(hist2_id, 'lai', kjit, lai, kjpindex*nvm, indexveg)
1277          CALL histwrite_p(hist2_id, 'subli', kjit, vevapsno, kjpindex, index)
1278          IF ( do_floodplains ) THEN
1279             CALL histwrite_p(hist2_id, 'vevapflo', kjit, vevapflo, kjpindex, index)
1280             CALL histwrite_p(hist2_id, 'flood_frac', kjit, flood_frac, kjpindex, index)
1281          ENDIF
1282          CALL histwrite_p(hist2_id, 'vevapnu', kjit, vevapnu, kjpindex, index)
1283          CALL histwrite_p(hist2_id, 'transpir', kjit, transpir, kjpindex*nvm, indexveg)
1284          CALL histwrite_p(hist2_id, 'inter', kjit, vevapwet, kjpindex*nvm, indexveg)
1285          CALL histwrite_p(hist2_id, 'vbeta1', kjit, vbeta1, kjpindex, index)
1286          CALL histwrite_p(hist2_id, 'vbeta2', kjit, vbeta2, kjpindex*nvm, indexveg)
1287          CALL histwrite_p(hist2_id, 'vbeta3', kjit, vbeta3, kjpindex*nvm, indexveg)
1288          CALL histwrite_p(hist2_id, 'vbeta4', kjit, vbeta4, kjpindex, index)   
1289          CALL histwrite_p(hist2_id, 'vbeta5', kjit, vbeta5, kjpindex, index)   
1290          CALL histwrite_p(hist2_id, 'drysoil_frac', kjit, drysoil_frac, kjpindex, index)
1291          CALL histwrite_p(hist2_id, 'rveget', kjit, rveget, kjpindex*nvm, indexveg)
1292          CALL histwrite_p(hist2_id, 'rstruct', kjit, rstruct, kjpindex*nvm, indexveg)
1293          CALL histwrite_p(hist2_id, 'snow', kjit, snow, kjpindex, index)
1294          CALL histwrite_p(hist2_id, 'snowage', kjit, snow_age, kjpindex, index)
1295          CALL histwrite_p(hist2_id, 'snownobio', kjit, snow_nobio, kjpindex*nnobio, indexnobio)
1296          CALL histwrite_p(hist2_id, 'snownobioage', kjit, snow_nobio_age, kjpindex*nnobio, indexnobio)
1297          CALL histwrite_p(hist2_id, 'soilindex',  kjit, REAL(njsc, r_std), kjpindex, index)
1298          CALL histwrite_p(hist2_id, 'reinf_slope',  kjit, reinf_slope, kjpindex, index)
1299         
1300          CALL histwrite_p(hist2_id, 'gsmean', kjit, gsmean, kjpindex*nvm, indexveg)   
1301          CALL histwrite_p(hist2_id, 'gpp', kjit, gpp, kjpindex*nvm, indexveg)
1302          CALL histwrite_p(hist2_id, 'cimean', kjit, cimean, kjpindex*nvm, indexveg)   
1303         
1304          IF ( ok_stomate ) THEN
1305             CALL histwrite_p(hist2_id, 'nee', kjit, co2_flux/1.e3/one_day, kjpindex*nvm, indexveg)   
1306          ENDIF
1307       ELSE
1308          ! Write history file in ALMA format
1309          CALL histwrite_p(hist2_id, 'vegetfrac', kjit, veget, kjpindex*nvm, indexveg)
1310          CALL histwrite_p(hist2_id, 'maxvegetfrac', kjit, veget_max, kjpindex*nvm, indexveg)
1311          CALL histwrite_p(hist2_id, 'nobiofrac', kjit, frac_nobio, kjpindex*nnobio, indexnobio)
1312          CALL histwrite_p(hist2_id, 'ESoil', kjit, vevapnu, kjpindex, index)
1313          IF ( do_floodplains ) THEN
1314             CALL histwrite_p(hist2_id, 'EWater', kjit, vevapflo, kjpindex, index)
1315             CALL histwrite_p(hist2_id, 'FloodFrac', kjit, flood_frac, kjpindex, index)
1316          ENDIF
1317          CALL histwrite_p(hist2_id, 'SWE', kjit, snow, kjpindex, index)
1318          histvar(:)=zero
1319          DO jv=1,nvm
1320             histvar(:) = histvar(:) + transpir(:,jv)
1321          ENDDO
1322          CALL histwrite_p(hist2_id, 'TVeg', kjit, histvar, kjpindex, index)
1323          histvar(:)=zero
1324          DO jv=1,nvm
1325             histvar(:) = histvar(:) + vevapwet(:,jv)
1326          ENDDO
1327          CALL histwrite_p(hist2_id, 'ECanop', kjit, histvar, kjpindex, index)
1328          CALL histwrite_p(hist2_id, 'ACond', kjit, tq_cdrag, kjpindex, index)
1329          CALL histwrite_p(hist2_id, 'SnowFrac', kjit, vbeta1, kjpindex, index)
1330          CALL histwrite_p(hist2_id, 'GPP', kjit, gpp, kjpindex*nvm, indexveg)
1331         
1332          IF ( ok_stomate ) THEN
1333             CALL histwrite_p(hist2_id, 'NEE', kjit, co2_flux, kjpindex*nvm, indexveg)   
1334          ENDIF
1335       ENDIF ! almaoutput
1336    ENDIF ! hist2_id
1337
1338
1339    !! Change the vegetation fractions if a new map was read in slowproc. This is done
1340    !! after lcchange has been done in stomatelpj
1341    IF (done_stomate_lcchange) THEN
1342       CALL slowproc_change_frac(kjpindex, lai, &
1343                                 veget_max, veget, frac_nobio, totfrac_nobio, tot_bare_soil, soiltile, fraclut, nwdFraclut)
1344       done_stomate_lcchange=.FALSE.
1345    END IF
1346
1347    !! 14. If it is the last time step, write restart files
1348    IF (ldrestart_write) THEN
1349       CALL sechiba_finalize( &
1350            kjit,     kjpij,  kjpindex, index,   rest_id, &
1351            tq_cdrag, vevapp, fluxsens, fluxlat, tsol_rad)
1352    END IF
1353
1354  END SUBROUTINE sechiba_main
1355
1356
1357!!  =============================================================================================================================
1358!! SUBROUTINE:    sechiba_finalize
1359!!
1360!>\BRIEF          Finalize all modules by calling their "_finalize" subroutines.
1361!!
1362!! DESCRIPTION:   Finalize all modules by calling their "_finalize" subroutines. These subroutines will write variables to
1363!!                restart file.
1364!!
1365!! \n
1366!_ ==============================================================================================================================
1367
1368  SUBROUTINE sechiba_finalize( &
1369       kjit,     kjpij,  kjpindex, index,   rest_id, &
1370       tq_cdrag, vevapp, fluxsens, fluxlat, tsol_rad)
1371
1372!! 0.1 Input variables   
1373    INTEGER(i_std), INTENT(in)                               :: kjit              !! Time step number (unitless)
1374    INTEGER(i_std), INTENT(in)                               :: kjpij             !! Total size of the un-compressed grid
1375                                                                                  !! (unitless)
1376    INTEGER(i_std), INTENT(in)                               :: kjpindex          !! Domain size - terrestrial pixels only
1377                                                                                  !! (unitless)
1378    INTEGER(i_std),INTENT (in)                               :: rest_id           !! _Restart_ file identifier (unitless)
1379    INTEGER(i_std),DIMENSION (kjpindex), INTENT (in)         :: index             !! Indices of the pixels on the map.
1380                                                                                  !! Sechiba uses a reduced grid excluding oceans
1381                                                                                  !! ::index contains the indices of the
1382                                                                                  !! terrestrial pixels only! (unitless)
1383    REAL(r_std),DIMENSION (kjpindex), INTENT (in)           :: tsol_rad           !! Radiative surface temperature
1384                                                                                  !! @tex $(W m^{-2})$ @endtex
1385    REAL(r_std),DIMENSION (kjpindex), INTENT (in)           :: vevapp             !! Total of evaporation
1386                                                                                  !! @tex $(kg m^{-2} days^{-1})$ @endtex
1387    REAL(r_std),DIMENSION (kjpindex), INTENT (in)           :: fluxsens           !! Sensible heat flux
1388                                                                                  !! @tex $(W m^{-2})$ @endtex
1389    REAL(r_std),DIMENSION (kjpindex), INTENT (in)           :: fluxlat            !! Latent heat flux
1390                                                                                  !! @tex $(W m^{-2})$ @endtex
1391    REAL(r_std),DIMENSION (kjpindex), INTENT (in)           :: tq_cdrag           !! Surface drag coefficient (-)
1392
1393!! 0.2 Local variables
1394    INTEGER(i_std)                                          :: ji, jv             !! Index (unitless)
1395    REAL(r_std), DIMENSION(kjpindex)                        :: histvar            !! Computations for history files (unitless)
1396    CHARACTER(LEN=80)                                       :: var_name           !! To store variables names for I/O (unitless)
1397
1398
1399    !! Write restart file for the next simulation from SECHIBA and other modules
1400
1401    IF (printlev_loc>=3) WRITE (numout,*) 'Start sechiba_finalize for writing restart files'
1402
1403    !! 1. Call diffuco_finalize to write restart files
1404    CALL diffuco_finalize (kjit, kjpindex, rest_id, rstruct )
1405   
1406    !! 2. Call energy budget to write restart files
1407    CALL enerbil_finalize (kjit,   kjpindex,    rest_id,            &
1408                           evapot, evapot_corr, temp_sol, tsol_rad, &
1409                           qsurf,  fluxsens,    fluxlat,  vevapp )
1410   
1411    !! 3. Call hydrology to write restart files
1412    CALL hydrol_finalize( kjit,           kjpindex, rest_id,  vegstress,  &
1413         qsintveg,       humrel,   snow,     snow_age, snow_nobio, &
1414         snow_nobio_age, snowrho,  snowtemp, snowdz,     &
1415         snowheat,       snowgrain,    &
1416         drysoil_frac,   evap_bare_lim, evap_bare_lim_ns)
1417   
1418    !! 4. Call condveg to write surface variables to restart files
1419    CALL condveg_finalize (kjit, kjpindex, rest_id, z0m, z0h, roughheight)
1420   
1421    !! 5. Call soil thermodynamic to write restart files
1422    CALL thermosoil_finalize (kjit,    kjpindex, rest_id,   gtemp, &
1423         soilcap, soilflx, lambda_snow, cgrnd_snow, dgrnd_snow)
1424
1425
1426    !! 6. Add river routing to restart files 
1427    IF ( river_routing .AND. nbp_glo .GT. 1) THEN
1428       !! 6.1 Call river routing to write restart files
1429       CALL routing_finalize( kjit, kjpindex, rest_id, flood_frac, flood_res )
1430    ELSE
1431       !! 6.2 No routing, set variables to zero
1432       reinfiltration(:) = zero
1433       returnflow(:) = zero
1434       irrigation(:) = zero
1435       flood_frac(:) = zero
1436       flood_res(:) = zero
1437    ENDIF
1438   
1439    !! 7. Call slowproc_main to add 'daily' and annual variables to restart file
1440    CALL slowproc_finalize (kjit,       kjpindex,  rest_id,  index,      &
1441                            njsc,       lai,       height,   veget,      &
1442                            frac_nobio, veget_max, reinf_slope,          &
1443                            co2_to_bm,  assim_param, frac_age,   &
1444                            heat_Zimov, altmax,    depth_organic_soil)
1445
1446   
1447    IF (printlev_loc>=3) WRITE (numout,*) 'sechiba_finalize done'
1448   
1449  END SUBROUTINE sechiba_finalize
1450
1451 
1452!! ==============================================================================================================================\n
1453!! SUBROUTINE   : sechiba_init
1454!!
1455!>\BRIEF        Dynamic allocation of the variables, the dimensions of the
1456!! variables are determined by user-specified settings.
1457!!
1458!! DESCRIPTION  : The domain size (:: kjpindex) is used to allocate the correct
1459!! dimensions to all variables in sechiba. Depending on the variable, its
1460!! dimensions are also determined by the number of PFT's (::nvm), number of
1461!! soil types (::nstm), number of non-vegetative surface types (::nnobio),
1462!! number of soil levels (::ngrnd), number of soil layers in the hydrological
1463!! model (i.e. cwrr) (::nslm). Values for these variables are set in
1464!! constantes_soil.f90 and constantes_veg.f90.\n
1465!!
1466!! Memory is allocated for all Sechiba variables and new indexing tables
1467!! are build making use of both (::kjpij) and (::kjpindex). New indexing tables
1468!! are needed because a single pixel can contain several PFTs, soil types, etc.
1469!! The new indexing tables have separate indices for the different
1470!! PFTs, soil types, etc.\n
1471!!
1472!! RECENT CHANGE(S): None
1473!!
1474!! MAIN OUTPUT VARIABLE(S): Strictly speaking the subroutine has no output
1475!! variables. However, the routine allocates memory and builds new indexing
1476!! variables for later use.
1477!!
1478!! REFERENCE(S) : None
1479!!
1480!! FLOWCHART    : None
1481!! \n
1482!_ ================================================================================================================================
1483
1484  SUBROUTINE sechiba_init (kjit, kjpij, kjpindex, index, rest_id, lalo)
1485
1486!! 0.1 Input variables
1487 
1488    INTEGER(i_std), INTENT (in)                         :: kjit               !! Time step number (unitless)
1489    INTEGER(i_std), INTENT (in)                         :: kjpij              !! Total size of the un-compressed grid (unitless)
1490    INTEGER(i_std), INTENT (in)                         :: kjpindex           !! Domain size - terrestrial pixels only (unitless)
1491    INTEGER(i_std), INTENT (in)                         :: rest_id            !! _Restart_ file identifier (unitless)
1492    INTEGER(i_std),DIMENSION (kjpindex), INTENT (in)    :: index              !! Indeces of the points on the map (unitless)
1493    REAL(r_std),DIMENSION (kjpindex,2), INTENT (in)     :: lalo               !! Geographical coordinates (latitude,longitude)
1494                                                                              !! for pixels (degrees)
1495!! 0.2 Output variables
1496
1497!! 0.3 Modified variables
1498
1499!! 0.4 Local variables
1500
1501    INTEGER(i_std)                                      :: ier                !! Check errors in memory allocation (unitless)
1502    INTEGER(i_std)                                      :: ji, jv             !! Indeces (unitless)
1503!_ ==============================================================================================================================
1504
1505!! 1. Initialize variables
1506   
1507    ! Dynamic allocation with user-specified dimensions on first call
1508    IF (l_first_sechiba) THEN
1509       l_first_sechiba=.FALSE.
1510    ELSE
1511       CALL ipslerr_p(3,'sechiba_init',' l_first_sechiba false . we stop ','','')
1512    ENDIF
1513
1514    !! Initialize local printlev
1515    printlev_loc=get_printlev('sechiba')
1516   
1517
1518    !! 1.1 Initialize 3D vegetation indexation table
1519    ALLOCATE (indexveg(kjpindex*nvm),stat=ier)
1520    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexveg','','')
1521
1522    ALLOCATE (indexlai(kjpindex*(nlai+1)),stat=ier)
1523    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexlai','','')
1524
1525    ALLOCATE (indexsoil(kjpindex*nstm),stat=ier)
1526    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexsoil','','')
1527
1528    ALLOCATE (indexnobio(kjpindex*nnobio),stat=ier)
1529    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexnobio','','')
1530
1531    ALLOCATE (indexgrnd(kjpindex*ngrnd),stat=ier)
1532    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexgrnd','','')
1533
1534    ALLOCATE (indexsnow(kjpindex*nsnow),stat=ier)
1535    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexsnow','','')
1536
1537    ALLOCATE (indexlayer(kjpindex*nslm),stat=ier)
1538    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexlayer','','')
1539
1540    ALLOCATE (indexnslm(kjpindex*nslm),stat=ier)
1541    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexnslm','','')
1542
1543    ALLOCATE (indexalb(kjpindex*2),stat=ier)
1544    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for indexalb','','')
1545
1546    !! 1.2  Initialize 1D array allocation with restartable value
1547    ALLOCATE (flood_res(kjpindex),stat=ier)
1548    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for flood_res','','')
1549    flood_res(:) = undef_sechiba
1550
1551    ALLOCATE (flood_frac(kjpindex),stat=ier)
1552    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for kjpindex','','')
1553    flood_frac(:) = undef_sechiba
1554
1555    ALLOCATE (snow(kjpindex),stat=ier)
1556    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snow','','')
1557    snow(:) = undef_sechiba
1558
1559    ALLOCATE (snow_age(kjpindex),stat=ier)
1560    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snow_age','','')
1561    snow_age(:) = undef_sechiba
1562
1563    ALLOCATE (drysoil_frac(kjpindex),stat=ier)
1564    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for drysoil_frac','','')
1565
1566    ALLOCATE (evap_bare_lim(kjpindex),stat=ier)
1567    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for evap_bare_lim','','')
1568
1569    ALLOCATE (evap_bare_lim_ns(kjpindex,nstm),stat=ier)
1570    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for evap_bare_lim_ns','','')
1571
1572    ALLOCATE (evapot(kjpindex),stat=ier)
1573    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for evapot','','')
1574    evapot(:) = undef_sechiba
1575
1576    ALLOCATE (evapot_corr(kjpindex),stat=ier)
1577    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for evapot_corr','','')
1578
1579    ALLOCATE (humrel(kjpindex,nvm),stat=ier)
1580    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for humrel','','')
1581    humrel(:,:) = undef_sechiba
1582
1583    ALLOCATE (vegstress(kjpindex,nvm),stat=ier)
1584    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vegstress','','')
1585    vegstress(:,:) = undef_sechiba
1586
1587    ALLOCATE (njsc(kjpindex),stat=ier)
1588    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for njsc','','')
1589    njsc(:)= undef_int
1590
1591    ALLOCATE (soiltile(kjpindex,nstm),stat=ier)
1592    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for soiltile','','')
1593
1594    ALLOCATE (fraclut(kjpindex,nlut),stat=ier)
1595    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for fraclut','','')
1596
1597    ALLOCATE (nwdFraclut(kjpindex,nlut),stat=ier)
1598    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for nwdFraclut','','')
1599
1600    ALLOCATE (reinf_slope(kjpindex),stat=ier)
1601    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for reinf_slope','','')
1602
1603    ALLOCATE (vbeta1(kjpindex),stat=ier)
1604    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta1','','')
1605
1606    ALLOCATE (vbeta4(kjpindex),stat=ier)
1607    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta4','','')
1608
1609    ALLOCATE (vbeta5(kjpindex),stat=ier)
1610    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta5','','')
1611
1612    ALLOCATE (soilcap(kjpindex),stat=ier)
1613    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for soilcap','','')
1614
1615    ALLOCATE (soilflx(kjpindex),stat=ier)
1616    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for soilflx','','')
1617
1618    ALLOCATE (temp_sol(kjpindex),stat=ier)
1619    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for temp_sol','','')
1620    temp_sol(:) = undef_sechiba
1621
1622    ALLOCATE (qsurf(kjpindex),stat=ier)
1623    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for qsurf','','')
1624    qsurf(:) = undef_sechiba
1625
1626    !! 1.3 Initialize 2D array allocation with restartable value
1627    ALLOCATE (qsintveg(kjpindex,nvm),stat=ier)
1628    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for qsintveg','','')
1629    qsintveg(:,:) = undef_sechiba
1630
1631    ALLOCATE (vbeta2(kjpindex,nvm),stat=ier)
1632    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta2','','')
1633
1634    ALLOCATE (vbeta3(kjpindex,nvm),stat=ier)
1635    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta3','','')
1636
1637    ALLOCATE (vbeta3pot(kjpindex,nvm),stat=ier)
1638    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta3pot','','')
1639
1640    ALLOCATE (gsmean(kjpindex,nvm),stat=ier)
1641    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for gsmean','','')
1642
1643    ALLOCATE (cimean(kjpindex,nvm),stat=ier)
1644    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for cimean','','')
1645
1646    ALLOCATE (gpp(kjpindex,nvm),stat=ier)
1647    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for gpp','','')
1648    gpp(:,:) = undef_sechiba
1649
1650 
1651    ALLOCATE (temp_growth(kjpindex),stat=ier) 
1652    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for temp_growth','','')
1653    temp_growth(:) = undef_sechiba 
1654
1655    ALLOCATE (veget(kjpindex,nvm),stat=ier)
1656    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for veget','','')
1657    veget(:,:)=undef_sechiba
1658
1659    ALLOCATE (veget_max(kjpindex,nvm),stat=ier)
1660    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for veget_max','','')
1661
1662    ALLOCATE (tot_bare_soil(kjpindex),stat=ier)
1663    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for tot_bare_soil','','')
1664
1665    ALLOCATE (lai(kjpindex,nvm),stat=ier)
1666    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for lai','','')
1667    lai(:,:)=undef_sechiba
1668
1669    ALLOCATE (frac_age(kjpindex,nvm,nleafages),stat=ier)
1670    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for frac_age','','')
1671    frac_age(:,:,:)=undef_sechiba
1672
1673    ALLOCATE (height(kjpindex,nvm),stat=ier)
1674    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for height','','')
1675    height(:,:)=undef_sechiba
1676
1677    ALLOCATE (frac_nobio(kjpindex,nnobio),stat=ier)
1678    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for frac_nobio','','')
1679    frac_nobio(:,:) = undef_sechiba
1680
1681    ALLOCATE (snow_nobio(kjpindex,nnobio),stat=ier)
1682    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snow_nobio','','')
1683    snow_nobio(:,:) = undef_sechiba
1684
1685    ALLOCATE (snow_nobio_age(kjpindex,nnobio),stat=ier)
1686    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snow_nobio_age','','')
1687    snow_nobio_age(:,:) = undef_sechiba
1688
1689    ALLOCATE (assim_param(kjpindex,nvm,npco2),stat=ier)
1690    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for assim_param','','')
1691
1692    !! 1.4 Initialize 1D array allocation
1693    ALLOCATE (vevapflo(kjpindex),stat=ier)
1694    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vevapflo','','')
1695    vevapflo(:)=zero
1696
1697    ALLOCATE (vevapsno(kjpindex),stat=ier)
1698    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vevapsno','','')
1699
1700    ALLOCATE (vevapnu(kjpindex),stat=ier)
1701    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vevapnu','','')
1702
1703    ALLOCATE (totfrac_nobio(kjpindex),stat=ier)
1704    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for totfrac_nobio','','')
1705
1706    ALLOCATE (floodout(kjpindex),stat=ier)
1707    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for floodout','','')
1708
1709    ALLOCATE (runoff(kjpindex),stat=ier)
1710    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for runoff','','')
1711
1712    ALLOCATE (drainage(kjpindex),stat=ier)
1713    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for drainage','','')
1714
1715    ALLOCATE (returnflow(kjpindex),stat=ier)
1716    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for returnflow','','')
1717    returnflow(:) = zero
1718
1719    ALLOCATE (reinfiltration(kjpindex),stat=ier)
1720    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for reinfiltration','','')
1721    reinfiltration(:) = zero
1722
1723    ALLOCATE (irrigation(kjpindex),stat=ier)
1724    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for irrigation','','')
1725    irrigation(:) = zero
1726
1727    ALLOCATE (z0h(kjpindex),stat=ier)
1728    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for z0h','','')
1729
1730    ALLOCATE (z0m(kjpindex),stat=ier)
1731    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for z0m','','')
1732
1733    ALLOCATE (roughheight(kjpindex),stat=ier)
1734    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for roughheight','','')
1735
1736    ALLOCATE (emis(kjpindex),stat=ier)
1737    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for emis','','')
1738
1739    ALLOCATE (tot_melt(kjpindex),stat=ier)
1740    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for tot_melt','','')
1741
1742    ALLOCATE (vbeta(kjpindex),stat=ier)
1743    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vbeta','','')
1744
1745    ALLOCATE (rau(kjpindex),stat=ier)
1746    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for rau','','')
1747
1748    ALLOCATE (deadleaf_cover(kjpindex),stat=ier)
1749    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for deadleaf_cover','','')
1750
1751    ALLOCATE (stempdiag(kjpindex, nslm),stat=ier)
1752    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for stempdiag','','')
1753
1754    ALLOCATE (co2_flux(kjpindex,nvm),stat=ier)
1755    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for co2_flux','','')
1756    co2_flux(:,:)=zero
1757
1758    ALLOCATE (co2_to_bm(kjpindex,nvm),stat=ier)
1759    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for co2_to_bm','','')
1760
1761    ALLOCATE (shumdiag(kjpindex,nslm),stat=ier)
1762    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for shumdiag','','')
1763   
1764    ALLOCATE (shumdiag_perma(kjpindex,nslm),stat=ier)
1765    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for shumdiag_perma','','')
1766
1767    ALLOCATE (litterhumdiag(kjpindex),stat=ier)
1768    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for litterhumdiag','','')
1769
1770    ALLOCATE (ptnlev1(kjpindex),stat=ier)
1771    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for ptnlev1','','')
1772
1773    ALLOCATE (k_litt(kjpindex),stat=ier)
1774    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for k_litt','','')
1775
1776    !! 1.5 Initialize 2D array allocation
1777    ALLOCATE (vevapwet(kjpindex,nvm),stat=ier)
1778    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for vevapwet','','')
1779    vevapwet(:,:)=undef_sechiba
1780
1781    ALLOCATE (transpir(kjpindex,nvm),stat=ier)
1782    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for transpir','','')
1783
1784    ALLOCATE (transpot(kjpindex,nvm),stat=ier)
1785    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for transpot','','')
1786
1787    ALLOCATE (qsintmax(kjpindex,nvm),stat=ier)
1788    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for qsintmax','','')
1789
1790    ALLOCATE (rveget(kjpindex,nvm),stat=ier)
1791    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for rveget','','')
1792
1793    ALLOCATE (rstruct(kjpindex,nvm),stat=ier)
1794    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for rstruct','','')
1795
1796    ALLOCATE (pgflux(kjpindex),stat=ier)
1797    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for pgflux','','')
1798    pgflux(:)= 0.0
1799
1800    ALLOCATE (cgrnd_snow(kjpindex,nsnow),stat=ier)
1801    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for cgrnd_snow','','')
1802    cgrnd_snow(:,:) = 0
1803
1804    ALLOCATE (dgrnd_snow(kjpindex,nsnow),stat=ier)
1805    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for dgrnd_snow','','')
1806    dgrnd_snow(:,:) = 0
1807
1808    ALLOCATE (lambda_snow(kjpindex),stat=ier)
1809    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for lambda_snow','','')
1810    lambda_snow(:) = 0
1811
1812    ALLOCATE (temp_sol_add(kjpindex),stat=ier)
1813    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for temp_sol_add','','')
1814
1815    ALLOCATE (gtemp(kjpindex),stat=ier)
1816    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for gtemp','','')
1817
1818    ALLOCATE (frac_snow_veg(kjpindex),stat=ier)
1819    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for frac_snow_veg','','')
1820
1821    ALLOCATE (frac_snow_nobio(kjpindex,nnobio),stat=ier)
1822    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for frac_snow_nobio','','')
1823
1824    ALLOCATE (snowrho(kjpindex,nsnow),stat=ier)
1825    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snowrho','','')
1826
1827    ALLOCATE (snowheat(kjpindex,nsnow),stat=ier)
1828    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snowheat','','')
1829
1830    ALLOCATE (snowgrain(kjpindex,nsnow),stat=ier)
1831    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snowgrain','','')
1832
1833    ALLOCATE (snowtemp(kjpindex,nsnow),stat=ier)
1834    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snowtemp','','')
1835
1836    ALLOCATE (snowdz(kjpindex,nsnow),stat=ier)
1837    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for snowdz','','')
1838
1839    ALLOCATE (mc_layh(kjpindex, nslm),stat=ier)
1840    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for mc_layh','','')
1841
1842    ALLOCATE (mcl_layh(kjpindex, nslm),stat=ier)
1843    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for mcl_layh','','')
1844
1845    ALLOCATE (soilmoist(kjpindex, nslm),stat=ier)
1846    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for soilmoist','','')
1847
1848    ALLOCATE (mcl_layh_s(kjpindex, nslm, nstm),stat=ier)
1849    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for mcl_layh_s','','')
1850
1851    ALLOCATE (mc_layh_s(kjpindex, nslm, nstm),stat=ier)
1852    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for mc_layh_s','','')
1853
1854    ALLOCATE (soilmoist_s(kjpindex, nslm, nstm),stat=ier)
1855    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for soilmoist_s','','')
1856
1857    ALLOCATE(tdeep(kjpindex,ngrnd,nvm),stat=ier)
1858    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for tdeep','','')
1859
1860    ALLOCATE(hsdeep(kjpindex,ngrnd,nvm),stat=ier)
1861    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for hsdeep','','')
1862
1863    ALLOCATE(heat_Zimov(kjpindex,ngrnd,nvm),stat=ier)
1864    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for heat_Zimov','','')
1865
1866    ALLOCATE(sfluxCH4_deep(kjpindex),stat=ier)
1867    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for sfluxCH4_deep','','')
1868
1869    ALLOCATE(sfluxCO2_deep(kjpindex),stat=ier)
1870    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for sfluxCO2_deep','','')
1871
1872    ALLOCATE(som_total(kjpindex,ngrnd,nvm,nelements),stat=ier)
1873    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for som_total','','')
1874
1875    ALLOCATE (altmax(kjpindex,nvm),stat=ier)
1876    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for altmax','','')
1877    ALLOCATE(depth_organic_soil(kjpindex),stat=ier)
1878    IF (ier /= 0) CALL ipslerr_p(3,'sechiba_init','Pb in alloc for depth_organic_soil','','')
1879
1880    !! 1.6 Initialize indexing table for the vegetation fields.
1881    ! In SECHIBA we work on reduced grids but to store in the full 3D filed vegetation variable
1882    ! we need another index table : indexveg, indexsoil, indexnobio and indexgrnd
1883    DO ji = 1, kjpindex
1884       !
1885       DO jv = 1, nlai+1
1886          indexlai((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1887       ENDDO
1888       !
1889       DO jv = 1, nvm
1890          indexveg((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1891       ENDDO
1892       !     
1893       DO jv = 1, nstm
1894          indexsoil((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1895       ENDDO
1896       !     
1897       DO jv = 1, nnobio
1898          indexnobio((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1899       ENDDO
1900       !
1901       DO jv = 1, ngrnd
1902          indexgrnd((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1903       ENDDO
1904       !
1905       DO jv = 1, nsnow
1906          indexsnow((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij
1907       ENDDO
1908
1909       DO jv = 1, nslm
1910          indexnslm((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij
1911       ENDDO
1912
1913       DO jv = 1, nslm
1914          indexlayer((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1915       ENDDO
1916       !
1917       DO jv = 1, 2
1918          indexalb((jv-1)*kjpindex + ji) = INDEX(ji) + (jv-1)*kjpij + offset_omp - offset_mpi
1919       ENDDO
1920       !
1921    ENDDO
1922
1923!! 2. Read the default value that will be put into variable which are not in the restart file
1924    CALL ioget_expval(val_exp)
1925   
1926    IF (printlev_loc>=3) WRITE (numout,*) ' sechiba_init done '
1927
1928  END SUBROUTINE sechiba_init
1929 
1930
1931!! ==============================================================================================================================\n
1932!! SUBROUTINE   : sechiba_clear
1933!!
1934!>\BRIEF        Deallocate memory of sechiba's variables
1935!!
1936!! DESCRIPTION  : None
1937!!
1938!! RECENT CHANGE(S): None
1939!!
1940!! MAIN OUTPUT VARIABLE(S): None
1941!!
1942!! REFERENCE(S) : None
1943!!
1944!! FLOWCHART    : None
1945!! \n
1946!_ ================================================================================================================================
1947
1948  SUBROUTINE sechiba_clear()
1949
1950!! 1. Initialize first run
1951
1952    l_first_sechiba=.TRUE.
1953
1954!! 2. Deallocate dynamic variables of sechiba
1955
1956    IF ( ALLOCATED (indexveg)) DEALLOCATE (indexveg)
1957    IF ( ALLOCATED (indexlai)) DEALLOCATE (indexlai)
1958    IF ( ALLOCATED (indexsoil)) DEALLOCATE (indexsoil)
1959    IF ( ALLOCATED (indexnobio)) DEALLOCATE (indexnobio)
1960    IF ( ALLOCATED (indexsnow)) DEALLOCATE (indexsnow)
1961    IF ( ALLOCATED (indexgrnd)) DEALLOCATE (indexgrnd)
1962    IF ( ALLOCATED (indexlayer)) DEALLOCATE (indexlayer)
1963    IF ( ALLOCATED (indexnslm)) DEALLOCATE (indexnslm)
1964    IF ( ALLOCATED (indexalb)) DEALLOCATE (indexalb)
1965    IF ( ALLOCATED (flood_res)) DEALLOCATE (flood_res)
1966    IF ( ALLOCATED (flood_frac)) DEALLOCATE (flood_frac)
1967    IF ( ALLOCATED (snow)) DEALLOCATE (snow)
1968    IF ( ALLOCATED (snow_age)) DEALLOCATE (snow_age)
1969    IF ( ALLOCATED (drysoil_frac)) DEALLOCATE (drysoil_frac)
1970    IF ( ALLOCATED (evap_bare_lim)) DEALLOCATE (evap_bare_lim)
1971    IF ( ALLOCATED (evap_bare_lim_ns)) DEALLOCATE (evap_bare_lim_ns)
1972    IF ( ALLOCATED (evapot)) DEALLOCATE (evapot)
1973    IF ( ALLOCATED (evapot_corr)) DEALLOCATE (evapot_corr)
1974    IF ( ALLOCATED (humrel)) DEALLOCATE (humrel)
1975    IF ( ALLOCATED (vegstress)) DEALLOCATE (vegstress)
1976    IF ( ALLOCATED (soiltile)) DEALLOCATE (soiltile)
1977    IF ( ALLOCATED (fraclut)) DEALLOCATE (fraclut)
1978    IF ( ALLOCATED (nwdFraclut)) DEALLOCATE (nwdFraclut)
1979    IF ( ALLOCATED (njsc)) DEALLOCATE (njsc)
1980    IF ( ALLOCATED (reinf_slope)) DEALLOCATE (reinf_slope)
1981    IF ( ALLOCATED (vbeta1)) DEALLOCATE (vbeta1)
1982    IF ( ALLOCATED (vbeta4)) DEALLOCATE (vbeta4)
1983    IF ( ALLOCATED (vbeta5)) DEALLOCATE (vbeta5)
1984    IF ( ALLOCATED (soilcap)) DEALLOCATE (soilcap)
1985    IF ( ALLOCATED (soilflx)) DEALLOCATE (soilflx)
1986    IF ( ALLOCATED (temp_sol)) DEALLOCATE (temp_sol)
1987    IF ( ALLOCATED (qsurf)) DEALLOCATE (qsurf)
1988    IF ( ALLOCATED (qsintveg)) DEALLOCATE (qsintveg)
1989    IF ( ALLOCATED (vbeta2))  DEALLOCATE (vbeta2)
1990    IF ( ALLOCATED (vbeta3)) DEALLOCATE (vbeta3)
1991    IF ( ALLOCATED (vbeta3pot)) DEALLOCATE (vbeta3pot)
1992    IF ( ALLOCATED (gsmean)) DEALLOCATE (gsmean)
1993    IF ( ALLOCATED (cimean)) DEALLOCATE (cimean)
1994    IF ( ALLOCATED (gpp)) DEALLOCATE (gpp)
1995    IF ( ALLOCATED (temp_growth)) DEALLOCATE (temp_growth) 
1996    IF ( ALLOCATED (veget)) DEALLOCATE (veget)
1997    IF ( ALLOCATED (veget_max)) DEALLOCATE (veget_max)
1998    IF ( ALLOCATED (tot_bare_soil)) DEALLOCATE (tot_bare_soil)
1999    IF ( ALLOCATED (lai)) DEALLOCATE (lai)
2000    IF ( ALLOCATED (frac_age)) DEALLOCATE (frac_age)
2001    IF ( ALLOCATED (height)) DEALLOCATE (height)
2002    IF ( ALLOCATED (roughheight)) DEALLOCATE (roughheight)
2003    IF ( ALLOCATED (frac_nobio)) DEALLOCATE (frac_nobio)
2004    IF ( ALLOCATED (snow_nobio)) DEALLOCATE (snow_nobio)
2005    IF ( ALLOCATED (snow_nobio_age)) DEALLOCATE (snow_nobio_age)
2006    IF ( ALLOCATED (assim_param)) DEALLOCATE (assim_param)
2007    IF ( ALLOCATED (vevapflo)) DEALLOCATE (vevapflo)
2008    IF ( ALLOCATED (vevapsno)) DEALLOCATE (vevapsno)
2009    IF ( ALLOCATED (vevapnu)) DEALLOCATE (vevapnu)
2010    IF ( ALLOCATED (totfrac_nobio)) DEALLOCATE (totfrac_nobio)
2011    IF ( ALLOCATED (floodout)) DEALLOCATE (floodout)
2012    IF ( ALLOCATED (runoff)) DEALLOCATE (runoff)
2013    IF ( ALLOCATED (drainage)) DEALLOCATE (drainage)
2014    IF ( ALLOCATED (reinfiltration)) DEALLOCATE (reinfiltration)
2015    IF ( ALLOCATED (irrigation)) DEALLOCATE (irrigation)
2016    IF ( ALLOCATED (tot_melt)) DEALLOCATE (tot_melt)
2017    IF ( ALLOCATED (vbeta)) DEALLOCATE (vbeta)
2018    IF ( ALLOCATED (rau)) DEALLOCATE (rau)
2019    IF ( ALLOCATED (deadleaf_cover)) DEALLOCATE (deadleaf_cover)
2020    IF ( ALLOCATED (stempdiag)) DEALLOCATE (stempdiag)
2021    IF ( ALLOCATED (co2_flux)) DEALLOCATE (co2_flux)
2022    IF ( ALLOCATED (co2_to_bm)) DEALLOCATE (co2_to_bm)
2023    IF ( ALLOCATED (shumdiag)) DEALLOCATE (shumdiag)
2024    IF ( ALLOCATED (shumdiag_perma)) DEALLOCATE (shumdiag_perma)
2025    IF ( ALLOCATED (litterhumdiag)) DEALLOCATE (litterhumdiag)
2026    IF ( ALLOCATED (ptnlev1)) DEALLOCATE (ptnlev1)
2027    IF ( ALLOCATED (k_litt)) DEALLOCATE (k_litt)
2028    IF ( ALLOCATED (vevapwet)) DEALLOCATE (vevapwet)
2029    IF ( ALLOCATED (transpir)) DEALLOCATE (transpir)
2030    IF ( ALLOCATED (transpot)) DEALLOCATE (transpot)
2031    IF ( ALLOCATED (qsintmax)) DEALLOCATE (qsintmax)
2032    IF ( ALLOCATED (rveget)) DEALLOCATE (rveget)
2033    IF ( ALLOCATED (rstruct)) DEALLOCATE (rstruct)
2034    IF ( ALLOCATED (frac_snow_veg)) DEALLOCATE (frac_snow_veg)
2035    IF ( ALLOCATED (frac_snow_nobio)) DEALLOCATE (frac_snow_nobio)
2036    IF ( ALLOCATED (snowrho)) DEALLOCATE (snowrho)
2037    IF ( ALLOCATED (snowgrain)) DEALLOCATE (snowgrain)
2038    IF ( ALLOCATED (snowtemp)) DEALLOCATE (snowtemp)
2039    IF ( ALLOCATED (snowdz)) DEALLOCATE (snowdz)
2040    IF ( ALLOCATED (snowheat)) DEALLOCATE (snowheat)
2041    IF ( ALLOCATED (cgrnd_snow)) DEALLOCATE (cgrnd_snow)
2042    IF ( ALLOCATED (dgrnd_snow)) DEALLOCATE (dgrnd_snow)
2043    IF ( ALLOCATED (lambda_snow)) DEALLOCATE(lambda_snow) 
2044    IF ( ALLOCATED (gtemp)) DEALLOCATE (gtemp)
2045    IF ( ALLOCATED (pgflux)) DEALLOCATE (pgflux)
2046    IF ( ALLOCATED (mc_layh)) DEALLOCATE (mc_layh)
2047    IF ( ALLOCATED (mcl_layh)) DEALLOCATE (mcl_layh)
2048    IF ( ALLOCATED (soilmoist)) DEALLOCATE (soilmoist)
2049    IF ( ALLOCATED (mcl_layh_s)) DEALLOCATE (mc_layh_s)
2050    IF ( ALLOCATED (mc_layh_s)) DEALLOCATE (mc_layh_s)
2051    IF ( ALLOCATED (soilmoist_s)) DEALLOCATE (soilmoist_s)
2052    IF ( ALLOCATED (tdeep)) DEALLOCATE (tdeep)
2053    IF ( ALLOCATED (hsdeep)) DEALLOCATE (hsdeep)
2054    IF ( ALLOCATED (heat_Zimov)) DEALLOCATE (heat_Zimov)
2055    IF ( ALLOCATED (sfluxCH4_deep)) DEALLOCATE (sfluxCH4_deep)
2056    IF ( ALLOCATED (sfluxCO2_deep)) DEALLOCATE (sfluxCO2_deep)
2057    IF ( ALLOCATED (som_total)) DEALLOCATE (som_total)
2058    IF ( ALLOCATED (altmax)) DEALLOCATE (altmax)
2059
2060!! 3. Clear all allocated memory
2061
2062    CALL pft_parameters_clear
2063    CALL slowproc_clear 
2064    CALL diffuco_clear 
2065    CALL enerbil_clear 
2066    CALL hydrol_clear 
2067    CALL thermosoil_clear
2068    CALL condveg_clear 
2069    CALL routing_clear
2070
2071  END SUBROUTINE sechiba_clear
2072
2073
2074!! ==============================================================================================================================\n
2075!! SUBROUTINE   : sechiba_var_init
2076!!
2077!>\BRIEF        Calculate air density as a function of air temperature and
2078!! pressure for each terrestrial pixel.
2079!!
2080!! RECENT CHANGE(S): None
2081!!
2082!! MAIN OUTPUT VARIABLE(S): air density (::rau, kg m^{-3}).
2083!!
2084!! REFERENCE(S) : None
2085!!
2086!! FLOWCHART    : None
2087!! \n
2088!_ ================================================================================================================================
2089
2090  SUBROUTINE sechiba_var_init (kjpindex, rau, pb, temp_air) 
2091
2092!! 0.1 Input variables
2093
2094    INTEGER(i_std), INTENT (in)                    :: kjpindex        !! Domain size - terrestrial pixels only (unitless)
2095    REAL(r_std),DIMENSION (kjpindex), INTENT (in)  :: pb              !! Surface pressure (hPa)
2096    REAL(r_std),DIMENSION (kjpindex), INTENT (in)  :: temp_air        !! Air temperature (K)
2097   
2098!! 0.2 Output variables
2099
2100    REAL(r_std),DIMENSION (kjpindex), INTENT (out) :: rau             !! Air density @tex $(kg m^{-3})$ @endtex
2101
2102!! 0.3 Modified variables
2103
2104!! 0.4 Local variables
2105
2106    INTEGER(i_std)                                 :: ji              !! Indices (unitless)
2107!_ ================================================================================================================================
2108   
2109!! 1. Calculate intial air density (::rau)
2110   
2111    DO ji = 1,kjpindex
2112       rau(ji) = pa_par_hpa * pb(ji) / (cte_molr*temp_air(ji))
2113    END DO
2114
2115    IF (printlev_loc>=3) WRITE (numout,*) ' sechiba_var_init done '
2116
2117  END SUBROUTINE sechiba_var_init
2118
2119
2120!! ==============================================================================================================================\n
2121!! SUBROUTINE   : sechiba_end
2122!!
2123!>\BRIEF        Swap old for newly calculated soil temperature.
2124!!
2125!! RECENT CHANGE(S): None
2126!!
2127!! MAIN OUTPUT VARIABLE(S): soil temperature (::temp_sol; K)
2128!!
2129!! REFERENCE(S) : None
2130!!
2131!! FLOWCHART    : None
2132!! \n
2133!! ================================================================================================================================
2134
2135  SUBROUTINE sechiba_end (kjpindex, temp_sol_new, temp_sol)
2136                         
2137
2138!! 0.1 Input variables
2139
2140    INTEGER(i_std), INTENT (in)                       :: kjpindex           !! Domain size - terrestrial pixels only (unitless)
2141    REAL(r_std),DIMENSION (kjpindex), INTENT (in)     :: temp_sol_new       !! New soil temperature (K)
2142   
2143    !! 0.2 Output variables
2144    REAL(r_std),DIMENSION (kjpindex), INTENT (out)    :: temp_sol           !! Soil temperature (K)
2145
2146!_ ================================================================================================================================
2147   
2148!! 1. Swap temperature
2149
2150    temp_sol(:) = temp_sol_new(:)
2151   
2152    IF (printlev_loc>=3) WRITE (numout,*) ' sechiba_end done '
2153
2154  END SUBROUTINE sechiba_end
2155
2156!! ==============================================================================================================================\n
2157!! SUBROUTINE   : sechiba_interface_orchidee_inca
2158!!
2159!>\BRIEF        make the interface between surface and atmospheric chemistry
2160!!
2161!! DESCRIPTION  : This subroutine is called from INCA, the atmospheric chemistry model. It is used to transfer variables from ORCHIDEE to INCA.
2162!!
2163!! RECENT CHANGE(S): move from chemistry module to be more generic (feb - 2017)
2164!!
2165!! MAIN OUTPUT VARIABLE(S): emission COV to be transport by orchidee to inca in fields_out array
2166!!
2167!! REFERENCE(S) : None
2168!!
2169!! FLOWCHART    : None
2170!! \n
2171!! ================================================================================================================================
2172  SUBROUTINE sechiba_interface_orchidee_inca( &
2173       nvm_out, veget_max_out, veget_frac_out, lai_out, snow_out, &
2174       field_out_names, fields_out, field_in_names, fields_in)
2175
2176
2177    INTEGER, INTENT(out)                      :: nvm_out            !! Number of vegetation types
2178    REAL(r_std), DIMENSION (:,:), INTENT(out) :: veget_max_out      !! Max. fraction of vegetation type (LAI -> infty)
2179    REAL(r_std), DIMENSION (:,:), INTENT(out) :: veget_frac_out     !! Fraction of vegetation type (unitless, 0-1) 
2180    REAL(r_std), DIMENSION (:,:), INTENT(out) :: lai_out            !! Surface foliere
2181    REAL(r_std), DIMENSION (:)  , INTENT(out) :: snow_out           !! Snow mass [Kg/m^2]
2182
2183    !
2184    ! Optional arguments
2185    !
2186    ! Names and fields for emission variables : to be transport by Orchidee to Inca
2187    CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_out_names
2188    REAL(r_std),DIMENSION(:,:,:), OPTIONAL, INTENT(OUT) :: fields_out
2189    !
2190    ! Names and fields for deposit variables : to be transport from chemistry model by INCA to ORCHIDEE.
2191    CHARACTER(LEN=*),DIMENSION(:), OPTIONAL, INTENT(IN) :: field_in_names
2192    REAL(r_std),DIMENSION(:,:), OPTIONAL, INTENT(IN)    :: fields_in
2193
2194
2195    ! Variables always transmitted from sechiba to inca
2196    nvm_out = nvm 
2197    veget_max_out(:,:)  = veget_max(:,:) 
2198    veget_frac_out(:,:) = veget(:,:) 
2199    lai_out(:,:)  = lai(:,:) 
2200    snow_out(:)  = snow(:) 
2201
2202    ! Call chemistry_flux_interface if at least one of variables field_out_names or
2203    ! field_in_names is present in the argument list of sechiba_interface_orchidee_inca when called from inca.
2204    IF (PRESENT(field_out_names) .AND. .NOT. PRESENT(field_in_names)) THEN
2205       CALL chemistry_flux_interface(field_out_names=field_out_names, fields_out=fields_out)
2206    ELSE IF (.NOT. PRESENT(field_out_names) .AND. PRESENT(field_in_names)) THEN
2207       CALL chemistry_flux_interface(field_in_names=field_in_names, fields_in=fields_in)
2208    ELSE IF (PRESENT(field_out_names) .AND. PRESENT(field_in_names)) THEN
2209       CALL chemistry_flux_interface(field_out_names=field_out_names, fields_out=fields_out, &
2210            field_in_names=field_in_names, fields_in=fields_in)
2211    ENDIF
2212
2213  END SUBROUTINE sechiba_interface_orchidee_inca
2214
2215
2216END MODULE sechiba
2217
Note: See TracBrowser for help on using the repository browser.