Context Navigation

source: tags/ORCHIDEE_1_9_6/ORCHIDEE/src_stomate/stomate.f90 @ 880

Last change on this file since 880 was 843, checked in by didier.solyga, 12 years ago
Improved parameters documentation and presentation. The presentation is the same for all paramters so a script can automatically build a orchidee.default file.
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 134.2 KB

Line
1	!< $HeadURL$
2	!< $Date$
3	!< $Author$
4	!< $Revision$
5	!IPSL (2006)
6	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
7	!
8	MODULE stomate
9	!---------------------------------------------------------------------
10	! Daily update of leaf area index etc.
11	!---------------------------------------------------------------------
12	USE netcdf
13	!-
14	USE ioipsl
15	USE defprec
16	USE grid
17	USE constantes
18	USE pft_parameters
19	USE stomate_io
20	USE stomate_data
21	USE stomate_season
22	USE stomate_lpj
23	USE stomate_assimtemp
24	USE stomate_litter
25	USE stomate_vmax
26	USE stomate_soilcarbon
27	USE stomate_resp
28	USE parallel
29	! USE Write_field_p
30	!-
31	IMPLICIT NONE
32	PRIVATE
33	PUBLIC stomate_main,stomate_clear,init_forcing,forcing_read
34	!
35	INTEGER,PARAMETER :: r_typ =nf90_real4
36	!-
37	! variables used inside stomate module : declaration and initialisation
38	!-
39	! total natural space (fraction of total space)
40	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: veget_cov_max
41	! new year total natural space (fraction of total space)
42	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: veget_cov_max_new
43	! carbon pool: active, slow, or passive, (gC/(m**2 of ground))
44	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: carbon
45	! density of individuals (1/(m**2 of ground))
46	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: ind
47	! daily moisture availability
48	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: humrel_daily
49	! daily litter humidity
50	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: litterhum_daily
51	! daily 2 meter temperatures (K)
52	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: t2m_daily
53	! daily minimum 2 meter temperatures (K)
54	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: t2m_min_daily
55	! daily surface temperatures (K)
56	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: tsurf_daily
57	! daily soil temperatures (K)
58	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: tsoil_daily
59	! daily soil humidity
60	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: soilhum_daily
61	! daily precipitations (mm/day) (for phenology)
62	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: precip_daily
63	! daily gross primary productivity (gC/(m**2 of ground)/day)
64	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: gpp_daily
65	! daily net primary productivity (gC/(m**2 of ground)/day)
66	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: npp_daily
67	! Turnover rates (gC/(m**2 of ground)/day)
68	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: turnover_daily
69	! Turnover rates (gC/(m*2 of ground)dtradia)
70	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: turnover_littercalc
71	! Probability of fire
72	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: fireindex
73	! Longer term total litter above the ground (gC/m**2 of ground)
74	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: firelitter
75	! "monthly" moisture availability
76	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: humrel_month
77	! "weekly" moisture availability
78	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: humrel_week
79	! "long term" 2 meter temperatures (K)
80	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: t2m_longterm
81	! "long term" 2 meter reference temperatures (K)
82	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: tlong_ref
83	! "monthly" 2 meter temperatures (K)
84	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: t2m_month
85	! "weekly" 2 meter temperatures (K)
86	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: t2m_week
87	! "monthly" soil temperatures (K)
88	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: tsoil_month
89	! "monthly" soil humidity
90	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: soilhum_month
91	! growing degree days, threshold -5 deg C (for phenology)
92	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: gdd_m5_dormance
93	! growing degree days, since midwinter (for phenology)
94	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: gdd_midwinter
95	! number of chilling days since leaves were lost (for phenology)
96	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: ncd_dormance
97	! number of growing days, threshold -5 deg C (for phenology)
98	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: ngd_minus5
99	! last year's maximum moisture availability
100	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: maxhumrel_lastyear
101	! this year's maximum moisture availability
102	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: maxhumrel_thisyear
103	! last year's minimum moisture availability
104	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: minhumrel_lastyear
105	! this year's minimum moisture availability
106	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: minhumrel_thisyear
107	! last year's maximum "weekly" GPP (gC/m**2 covered/day)
108	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: maxgppweek_lastyear
109	! this year's maximum "weekly" GPP (gC/m**2 covered/day)
110	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: maxgppweek_thisyear
111	! last year's annual GDD0
112	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: gdd0_lastyear
113	! this year's annual GDD0
114	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: gdd0_thisyear
115	! last year's annual precipitation (mm/year)
116	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: precip_lastyear
117	! this year's annual precipitation (mm/year)
118	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: precip_thisyear
119	! PFT exists (equivalent to veget > 0 for natural PFTs)
120	LOGICAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: PFTpresent
121	! "long term" net primary productivity (gC/(m**2 of ground)/year)
122	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: npp_longterm
123	! last year's maximum leaf mass, for each PFT (gC/(m**2 of ground))
124	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: lm_lastyearmax
125	! this year's maximum leaf mass, for each PFT (gC/(m**2 of ground))
126	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: lm_thisyearmax
127	! last year's maximum fpc for each natural PFT, on ground
128	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: maxfpc_lastyear
129	! this year's maximum fpc for each PFT, on ground (see stomate_season)
130	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: maxfpc_thisyear
131	! "long term" turnover rate (gC/(m**2 of ground)/year)
132	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: turnover_longterm
133	! "weekly" GPP (gC/day/(m**2 covered)
134	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: gpp_week
135	! biomass (gC/(m**2 of ground))
136	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: biomass
137	! is the plant senescent?
138	! (only for deciduous trees - carbohydrate reserve)
139	LOGICAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: senescence
140	! how many days ago was the beginning of the growing season
141	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: when_growthinit
142	! age (years). For trees, mean stand age.
143	! For grasses, ears since introduction of PFT
144	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: age
145	! Winter too cold? between 0 and 1
146	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: adapted
147	! Winter sufficiently cold? between 0 and 1
148	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: regenerate
149	! fraction of litter above the ground belonging to different PFTs
150	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: litterpart
151	! metabolic and structural litter, above and below ground
152	! (gC/(m**2 of ground))
153	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:,:):: litter
154	! dead leaves on ground, per PFT, metabolic and structural,
155	! in gC/(m**2 of ground)
156	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: dead_leaves
157	! black carbon on the ground (gC/(m**2 of total ground))
158	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: black_carbon
159	! ratio Lignine/Carbon in structural litter, above and below ground,
160	! (gC/(m**2 of ground))
161	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: lignin_struc
162	! carbon emitted into the atmosphere by fire (living and dead biomass)
163	! (in gC/m**2 of average ground/day)
164	!NV devient 2D
165	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: co2_fire
166	! co2 taken up (gC/(m**2 of total ground)/day)
167	!NV devient 2D
168	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: co2_to_bm_dgvm
169	! heterotrophic respiration (gC/day/m**2 of total ground)
170	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: resp_hetero_d
171	! maintenance respiration (gC/day/(m**2 of total ground))
172	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: resp_hetero_radia
173	! maintenance respiration (gC/day/(m**2 of total ground))
174	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: resp_maint_d
175	! growth respiration (gC/day/(m**2 of total ground))
176	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: resp_growth_d
177	! vegetation fractions (on ground) after last light competition
178	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: veget_lastlight
179	! is the PFT everywhere in the grid box or very localized (after its intoduction)
180	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: everywhere
181	! in order for this PFT to be introduced,
182	! does it have to be present in an adjacent grid box?
183	LOGICAL,ALLOCATABLE,SAVE,DIMENSION(:,:) :: need_adjacent
184	! leaf age (d)
185	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: leaf_age
186	! fraction of leaves in leaf age class
187	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: leaf_frac
188	! How much time ago was the PFT eliminated for the last time (y)
189	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: RIP_time
190	! duration of dormance (d)
191	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: time_lowgpp
192	! time elapsed since strongest moisture availability (d)
193	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: time_hum_min
194	! minimum moisture during dormance
195	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: hum_min_dormance
196	! time constant of probability of a leaf to be eaten by a herbivore (days)
197	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: herbivores
198	! npp total written for forcesoil...
199	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: npp_equil
200	! npp total ...
201	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: npp_tot
202	! moisture control of heterotrophic respiration
203	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: control_moist
204	! temperature control of heterotrophic respiration, above and below
205	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: control_temp
206	! quantity of carbon going into carbon pools from litter decomposition
207	! (gC/(m**2 of ground)/day)
208	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:,:):: soilcarbon_input
209	! quantity of carbon going into carbon pools from litter decomposition daily
210	! (gC/(m**2 of ground)/day)
211	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:):: soilcarbon_input_daily
212	! moisture control of heterotrophic respiration daily
213	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: control_moist_daily
214	! temperature control of heterotrophic respiration, above and below daily
215	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: control_temp_daily
216	! how many states were calculated for a given soil forcing time step
217	! turnover time of leaves
218	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: turnover_time
219	! daily total CO2 flux (gC/m**2/day)
220	!NV champs 2D
221	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: co2_flux_daily
222	! monthly total CO2 flux (gC/m**2)
223	!NV champs 2D
224	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: co2_flux_monthly
225	! conversion of biomass to litter (gC/(m**2 of ground))/day
226	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: bm_to_litter
227	! conversion of biomass to litter (gC/(m*2 of ground))dtradia
228	REAL(r_std), ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: bm_to_littercalc
229
230	INTEGER(i_std),ALLOCATABLE,SAVE,DIMENSION(:) :: nforce
231
232	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: harvest_above_monthly, cflux_prod_monthly
233
234	! "maximal" coverage fraction of a PFT (LAI -> infinity) on ground
235	REAL(r_std), ALLOCATABLE,SAVE,DIMENSION(:,:) :: fpc_max
236
237	! Date and EndOfYear, intialize and update in slowproc
238	! (Now managed in slowproc for land_use)
239	! time step of STOMATE in days
240	REAL(r_std),SAVE :: dt_days=zero ! Time step in days for stomate
241	! to check
242	REAL(r_std),SAVE :: day_counter=zero ! count each sechiba (dtradia) time step each day
243	! date (d)
244	INTEGER(i_std),SAVE :: date=0
245	! is it time for Stomate or update of LAI etc. ?
246	LOGICAL, SAVE :: do_slow=.FALSE.
247	! Do update of yearly variables ?
248	! This variable must be .TRUE. once a year
249	LOGICAL, SAVE :: EndOfYear=.FALSE.
250	! Do update of monthly variables ?
251	! This variable must be .TRUE. once a month
252	LOGICAL, SAVE :: EndOfMonth=.FALSE.
253	PUBLIC dt_days, day_counter, date, do_slow, EndOfYear
254
255
256	! forcing data in memory
257	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: clay_fm
258	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: humrel_daily_fm
259	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: litterhum_daily_fm
260	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: t2m_daily_fm
261	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: t2m_min_daily_fm
262	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: tsurf_daily_fm
263	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: tsoil_daily_fm
264	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: soilhum_daily_fm
265	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: precip_fm
266	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: gpp_daily_fm
267	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: veget_fm
268	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: veget_max_fm
269	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: lai_fm
270	PUBLIC clay_fm, humrel_daily_fm, litterhum_daily_fm, t2m_daily_fm, t2m_min_daily_fm, tsurf_daily_fm, tsoil_daily_fm, &
271	soilhum_daily_fm, precip_fm, gpp_daily_fm, veget_fm, veget_max_fm, lai_fm
272
273	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: clay_fm_g
274	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: humrel_daily_fm_g
275	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: litterhum_daily_fm_g
276	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: t2m_daily_fm_g
277	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: t2m_min_daily_fm_g
278	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: tsurf_daily_fm_g
279	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: tsoil_daily_fm_g
280	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: soilhum_daily_fm_g
281	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: precip_fm_g
282	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: gpp_daily_fm_g
283	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: veget_fm_g
284	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: veget_max_fm_g
285	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: lai_fm_g
286
287	INTEGER(i_std),ALLOCATABLE,SAVE,DIMENSION(:) :: isf
288	LOGICAL,ALLOCATABLE,SAVE,DIMENSION(:) :: nf_written
289	INTEGER(i_std),ALLOCATABLE,SAVE,DIMENSION(:) :: nf_cumul
290	PUBLIC isf, nf_written
291
292	! first call
293	LOGICAL,SAVE :: l_first_stomate = .TRUE.
294	! flag for cumul of forcing if teststomate
295	LOGICAL,SAVE :: cumul_forcing=.FALSE.
296	! flag for cumul of forcing if forcesoil
297	LOGICAL,SAVE :: cumul_Cforcing=.FALSE.
298	!
299	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: resp_maint_part_radia
300	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: resp_maint_part
301	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: resp_maint_radia
302
303	! Land cover change variables
304	! products remaining in the 10/100 year-turnover pool after the annual release for each compartment
305	! (10 or 100 + 1 : input from year of land cover change)
306	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: prod10
307	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: prod100
308	! annual release from the 10/100 year-turnover pool compartments
309	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: flux10
310	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: flux100
311	! release during first year following land cover change
312	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: convflux
313	! total annual release from the 10/100 year-turnover pool
314	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: cflux_prod10, cflux_prod100
315
316	! harvest above ground biomass for agriculture
317	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: harvest_above
318
319	! Carbon Mass total
320	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:) :: carb_mass_total
321
322	CONTAINS
323	!
324	!=
325	!
326	SUBROUTINE stomate_main &
327	& (kjit, kjpij, kjpindex, dtradia, dt_slow, &
328	& ldrestart_read, ldrestart_write, ldforcing_write, ldcarbon_write, &
329	& index, lalo, neighbours, resolution, contfrac, totfrac_nobio, clay, &
330	& t2m, t2m_min, temp_sol, stempdiag, &
331	& humrel, shumdiag, litterhumdiag, precip_rain, precip_snow, &
332	& gpp, deadleaf_cover, assim_param, &
333	& lai, height, veget, veget_max, qsintmax, &
334	& veget_max_new, totfrac_nobio_new, &
335	& hist_id, hist2_id, rest_id_stom, hist_id_stom, hist_id_stom_IPCC, &
336	& co2_flux, fco2_lu, resp_maint,resp_hetero,resp_growth)
337	!---------------------------------------------------------------------
338	!
339	IMPLICIT NONE
340	! 0 interface description
341	!
342	! 0.1 input
343	!
344	! 0.1.1 input scalar
345	!
346	! Time step number
347	INTEGER(i_std),INTENT(in) :: kjit
348	! Domain size
349	INTEGER(i_std),INTENT(in) :: kjpindex
350	! Total size of the un-compressed grid
351	INTEGER(i_std),INTENT(in) :: kjpij
352	! Time step of SECHIBA
353	REAL(r_std),INTENT(in) :: dtradia
354	! Time step of STOMATE
355	REAL(r_std),INTENT(in) :: dt_slow
356	! Logical for _restart_ file to read
357	LOGICAL,INTENT(in) :: ldrestart_read
358	! Logical for _restart_ file to write
359	LOGICAL,INTENT(in) :: ldrestart_write
360	! Logical for _forcing_ file to write
361	LOGICAL,INTENT(in) :: ldforcing_write
362	! Logical for _carbon_forcing_ file to write
363	LOGICAL,INTENT(in) :: ldcarbon_write
364	! SECHIBA's _history_ file identifier
365	INTEGER(i_std),INTENT(in) :: hist_id
366	! SECHIBA's _history_ file 2 identifier
367	INTEGER(i_std),INTENT(in) :: hist2_id
368	! STOMATE's _Restart_ file file identifier
369	INTEGER(i_std),INTENT(in) :: rest_id_stom
370	! STOMATE's _history_ file file identifier
371	INTEGER(i_std),INTENT(in) :: hist_id_stom
372	! STOMATE's IPCC _history_ file file identifier
373	INTEGER(i_std),INTENT(in) :: hist_id_stom_IPCC
374	!
375	! 0.1.2 input fields
376	!
377	! Indeces of the points on the map
378	INTEGER(i_std),DIMENSION(kjpindex),INTENT(in) :: index
379	! Geogr. coordinates (latitude,longitude) (degrees)
380	REAL(r_std),DIMENSION(kjpindex,2),INTENT(in) :: lalo
381	! neighoring grid points if land
382	INTEGER(i_std),DIMENSION(kjpindex,8),INTENT(in) :: neighbours
383	! size in x an y of the grid (m)
384	REAL(r_std),DIMENSION(kjpindex,2),INTENT(in) :: resolution
385	! fraction of continent in the grid
386	REAL(r_std),DIMENSION (kjpindex), INTENT (in) :: contfrac
387	! fraction of grid cell covered by lakes, land ice, cities, ...
388	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: totfrac_nobio
389	! clay fraction
390	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: clay
391	! Relative humidity ("moisture availability")
392	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: humrel
393	! 2 m air temperature (K)
394	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: t2m
395	! min. 2 m air temp. during forcing time step (K)
396	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: t2m_min
397	! Surface temperature
398	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: temp_sol
399	! Soil temperature
400	REAL(r_std),DIMENSION(kjpindex,nbdl),INTENT(in) :: stempdiag
401	! Relative soil moisture
402	REAL(r_std),DIMENSION(kjpindex,nbdl),INTENT(in) :: shumdiag
403	! Litter humidity
404	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: litterhumdiag
405	! Rain precipitation
406	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: precip_rain
407	! Snow precipitation
408	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: precip_snow
409	! GPP (gC/(m**2 of total ground)/time step)
410	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: gpp
411	! new "maximal" coverage fraction of a PFT (LAI -> infinity) on ground
412	! only if EndOfYear is activated
413	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: veget_max_new
414	! new fraction of grid cell covered by lakes, land ice, cities, ...
415	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: totfrac_nobio_new
416	!
417	! 0.2 output
418	!
419	! 0.2.1 output scalar
420	!
421	! 0.2.2 output fields
422	!
423	! CO2 flux in gC/m**2 of average ground/dt
424	!NV champs 2D
425	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out) :: co2_flux
426	REAL(r_std),DIMENSION(kjpindex),INTENT(out) :: fco2_lu
427	! autotrophic respiration in gC/m**2 of surface/dt
428	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out) :: resp_maint
429	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out) :: resp_growth
430	!NV champs2D
431	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out) :: resp_hetero
432	!
433	! 0.3 modified
434	!
435	! 0.3.1 modified scalar
436	! 0.3.2 modified fields
437	!
438	! Surface foliere
439	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(inout) :: lai
440	! Fraction of vegetation type from hydrological module. Takes into account ice etc.
441	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(inout) :: veget
442	! Maximum fraction of vegetation type from hydrological module. Takes into account ice etc.
443	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(inout) :: veget_max
444	! height (m)
445	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(inout) :: height
446	! min+max+opt temps & vmax for photosynthesis
447	REAL(r_std),DIMENSION(kjpindex,nvm,npco2),INTENT(inout) :: assim_param
448	! fraction of soil covered by dead leaves
449	REAL(r_std),DIMENSION(kjpindex),INTENT(inout) :: deadleaf_cover
450	! Maximum water on vegetation for interception
451	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(inout) :: qsintmax
452	!
453	! 0.4 local declaration
454	!
455	! 0.4.1 variables defined on nvm in SECHIBA, nvm in STOMATE/LPJ
456
457	! gross primary productivity (gC/(m**2 of total ground)/day)
458	REAL(r_std),DIMENSION(kjpindex,nvm) :: gpp_d
459	! fractional coverage: actually covered space, taking into account LAI (= grid scale fpc).
460	! Fraction of ground.
461	REAL(r_std),DIMENSION(kjpindex,nvm) :: veget_cov
462
463	! 0.4.2 other
464	!
465	! soil level used for LAI
466	INTEGER(i_std),SAVE :: lcanop
467	! counts time until next STOMATE time step
468	REAL(r_std) :: day_counter_read
469	! STOMATE time step read in restart file
470	REAL(r_std) :: dt_days_read
471	! STOMATE date
472	INTEGER(i_std) :: date_read
473	! Maximum STOMATE time step (days)
474	REAL(r_std),PARAMETER :: max_dt_days = 5.
475	! Writing frequency for history file (d)
476	REAL(r_std) :: hist_days
477	! precipitation (mm/day)
478	REAL(r_std),DIMENSION(kjpindex) :: precip
479	! Maximum rate of carboxylation
480	REAL(r_std),DIMENSION(kjpindex,nvm) :: vcmax
481	! Maximum rate of RUbp regeneration
482	REAL(r_std),DIMENSION(kjpindex,nvm) :: vjmax
483	! Min temperature for photosynthesis (deg C)
484	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_photo_min
485	! Opt temperature for photosynthesis (deg C)
486	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_photo_opt
487	! Max temperature for photosynthesis (deg C)
488	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_photo_max
489	!NV supprimé
490	! total autotrophic respiration (gC/day/(m**2 of total ground))
491	! REAL(r_std),DIMENSION(kjpindex) :: resp_auto_tot
492	!NV gpp_tot supprimé
493	! total photosynthesis (gC/day/(m**2 of total ground))
494	! REAL(r_std),DIMENSION(kjpindex) :: gpp_tot
495	! -- LOOP
496	REAL(r_std) :: net_co2_flux_monthly, net_co2_flux_monthly_sum
497	INTEGER :: ios
498	! for forcing file: "daily" gpp
499	REAL(r_std),DIMENSION(kjpindex,nvm) :: gpp_daily_x
500	! total "vegetation" cover
501	REAL(r_std),DIMENSION(kjpindex) :: cvegtot
502	INTEGER(i_std) :: ji, jv, i, j
503	REAL(r_std) :: trans_veg
504	REAL(r_std) :: tmp_day(1)
505	!
506
507	INTEGER(i_std) :: ier
508	! moisture control of heterotrophic respiration
509	REAL(r_std),DIMENSION(kjpindex,nlevs) :: control_moist_inst
510	! temperature control of heterotrophic respiration, above and below
511	REAL(r_std),DIMENSION(kjpindex,nlevs) :: control_temp_inst
512	! quantity of carbon going into carbon pools from litter decomposition
513	! (gC/(m**2 of ground)/day)
514	REAL(r_std),DIMENSION(kjpindex,ncarb,nvm) :: soilcarbon_input_inst
515	! time step of soil forcing file (in days)
516	REAL(r_std),SAVE :: dt_forcesoil
517	INTEGER(i_std),SAVE :: nparan ! Number of time steps per year for carbon spinup
518	INTEGER(i_std),SAVE :: nbyear
519	INTEGER(i_std),PARAMETER :: nparanmax=366 ! Number max of time steps per year for carbon spinup
520	REAL(r_std) :: sf_time
521	INTEGER(i_std),SAVE :: iatt
522	INTEGER(i_std),SAVE :: iatt_old=1
523	INTEGER(i_std) :: max_totsize, totsize_1step,totsize_tmp
524	REAL(r_std) :: xn
525	INTEGER(i_std),SAVE :: nsfm, nsft
526	INTEGER(i_std),SAVE :: iisf
527	!-
528	CHARACTER(LEN=100), SAVE :: forcing_name,Cforcing_name
529	INTEGER(i_std),SAVE :: Cforcing_id
530	INTEGER(i_std),PARAMETER :: ndm = 10
531	INTEGER(i_std) :: vid
532	INTEGER(i_std) :: nneigh,direct
533	INTEGER(i_std),DIMENSION(ndm) :: d_id
534
535
536	! root temperature (convolution of root and soil temperature profiles)
537	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_root
538	REAL(r_std),DIMENSION(kjpindex,nvm,nparts) :: coeff_maint
539	! temperature which is pertinent for maintenance respiration (K)
540	REAL(r_std),DIMENSION(kjpindex,nparts) :: t_maint_radia
541	! integration constant for root profile
542	REAL(r_std),DIMENSION(kjpindex) :: rpc
543	! long term annual mean temperature, C
544	REAL(r_std),DIMENSION(kjpindex) :: tl
545	! slope of maintenance respiration coefficient (1/K)
546	REAL(r_std),DIMENSION(kjpindex) :: slope
547	! root depth. This will, one day, be a prognostic variable.
548	! It will be calculated by
549	! STOMATE (save in restart file & give to hydrology module!),
550	! probably somewhere
551	! in the allocation routine. For the moment, it is prescribed.
552	REAL(r_std),DIMENSION(kjpindex,nvm) :: rprof
553	INTEGER(i_std) :: l,k
554	! litter heterotrophic respiration (gC/day/m**2 by PFT)
555	REAL(r_std),DIMENSION(kjpindex,nvm) :: resp_hetero_litter
556	! soil heterotrophic respiration (gC/day/m**2 by PFT)
557	REAL(r_std),DIMENSION(kjpindex,nvm) :: resp_hetero_soil
558	INTEGER(i_std) :: iyear
559
560	REAL(r_std),DIMENSION(nbp_glo) :: clay_g
561	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:,:) :: soilcarbon_input_g
562	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: control_moist_g
563	REAL(r_std),ALLOCATABLE,DIMENSION(:,:,:) :: control_temp_g
564	REAL(r_std),ALLOCATABLE,DIMENSION(:,:) :: npp_equil_g
565
566	REAL(r_std), DIMENSION(kjpindex) :: vartmp
567	REAL(r_std) :: net_cflux_prod_monthly_sum , net_cflux_prod_monthly_tot
568	REAL(r_std) :: net_harvest_above_monthly_sum, net_harvest_above_monthly_tot
569	REAL(r_std) :: net_biosp_prod_monthly_sum , net_biosp_prod_monthly_tot
570	!---------------------------------------------------------------------
571	! first of all: store time step in common value
572	itime = kjit
573
574	DO j=1,nvm
575	rprof(:,j) = 1./humcste(j)
576	ENDDO
577	!-
578	! 1 do initialisation
579	!-
580	resp_growth=zero
581	IF (l_first_stomate) THEN
582	IF (long_print) THEN
583	WRITE (numout,*) ' l_first_stomate : call stomate_init'
584	ENDIF
585	!
586	!
587	! 1.1 allocation, file definitions. Set flags.
588	!
589	CALL stomate_init (kjpij, kjpindex, index, ldforcing_write, lalo, &
590	rest_id_stom, hist_id_stom, hist_id_stom_IPCC)
591
592	!
593	! 1.2 read PFT data
594	!
595	CALL data
596	!
597	! 1.3 Initial conditions
598	!
599	! 1.3.1 read STOMATE's start file
600	!
601	co2_flux(:,:) = zero
602	fco2_lu(:) = zero
603	!
604	CALL readstart &
605	& (kjpindex, index, lalo, resolution, &
606	& day_counter_read, dt_days_read, date_read, &
607	& ind, adapted, regenerate, &
608	& humrel_daily, litterhum_daily, &
609	& t2m_daily, t2m_min_daily, tsurf_daily, tsoil_daily, &
610	& soilhum_daily, precip_daily, &
611	& gpp_daily, npp_daily, turnover_daily, &
612	& humrel_month, humrel_week, &
613	& t2m_longterm, tlong_ref, t2m_month, t2m_week, &
614	& tsoil_month, soilhum_month, fireindex, firelitter, &
615	& maxhumrel_lastyear, maxhumrel_thisyear, &
616	& minhumrel_lastyear, minhumrel_thisyear, &
617	& maxgppweek_lastyear, maxgppweek_thisyear, &
618	& gdd0_lastyear, gdd0_thisyear, &
619	& precip_lastyear, precip_thisyear, &
620	& gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
621	& PFTpresent, npp_longterm, lm_lastyearmax, lm_thisyearmax, &
622	& maxfpc_lastyear, maxfpc_thisyear, &
623	& turnover_longterm, gpp_week, biomass, resp_maint_part, &
624	& leaf_age, leaf_frac, &
625	& senescence, when_growthinit, age, &
626	& resp_hetero_d, resp_maint_d, resp_growth_d, co2_fire, co2_to_bm_dgvm, &
627	& veget_lastlight, everywhere, need_adjacent, RIP_time, time_lowgpp, &
628	& time_hum_min, hum_min_dormance, &
629	& litterpart, litter, dead_leaves, &
630	& carbon, black_carbon, lignin_struc,turnover_time,&
631	& prod10,prod100,flux10, flux100, &
632	& convflux, cflux_prod10, cflux_prod100, bm_to_litter, carb_mass_total)
633
634	! 1.4 read the boundary conditions
635	!
636	CALL readbc (kjpindex, lalo, resolution, tlong_ref)
637	!
638	! 1.5 check time step
639	!
640	! 1.5.1 allow STOMATE's time step to change
641	! although this is dangerous
642	!
643	IF (dt_days /= dt_days_read) THEN
644	WRITE(numout,*) 'slow_processes: STOMATE time step changes:', &
645	dt_days_read,' -> ',dt_days
646	ENDIF
647
648	! 1.5.2 time step has to be a multiple of a full day
649
650	IF ( ( dt_days-REAL(NINT(dt_days),r_std) ) > min_stomate ) THEN
651	WRITE(numout,*) 'slow_processes: STOMATE time step is wrong:', &
652	& dt_days,' days.'
653	STOP
654	ENDIF
655
656	! 1.5.3 upper limit to STOMATE's time step
657
658	IF ( dt_days > max_dt_days ) THEN
659	WRITE(numout,*) 'slow_processes: STOMATE time step too large:', &
660	& dt_days,' days.'
661	STOP
662	ENDIF
663
664	! 1.5.4 STOMATE time step must not be less than the forcing time step
665
666	IF ( dtradia > dt_days*one_day ) THEN
667	WRITE(numout,*) &
668	& 'slow_processes: STOMATE time step smaller than forcing time step.'
669	STOP
670	ENDIF
671
672	! 1.5.5 For teststomate : test day_counter
673	IF ( abs(day_counter - day_counter_read) > min_stomate ) THEN
674	WRITE(numout,*) 'slow_processes: STOMATE day counter changes:', &
675	day_counter_read,' -> ',day_counter
676	ENDIF
677
678	! 1.5.6 date check
679	IF (date /= date_read) THEN
680	WRITE(numout,*) 'slow_processes: STOMATE date changes:', &
681	date_read,' -> ',date
682	ENDIF
683
684	! 1.5.6 some more messages
685
686	WRITE(numout,*) 'slow_processes, STOMATE time step (d): ', dt_days
687
688	!
689	! 1.6 write forcing file for stomate?
690	!
691	IF (ldforcing_write) THEN
692
693	!Config Key = STOMATE_FORCING_NAME
694	!Config Desc = Name of STOMATE's forcing file
695	!Config If = OK_STOMATE
696	!Config Def = NONE
697	!Config Help = Name that will be given
698	!Config to STOMATE's offline forcing file
699	!Config Units = [FILE]
700	!-
701	forcing_name = stomate_forcing_name ! compatibilité avec driver Nicolas
702	CALL getin_p('STOMATE_FORCING_NAME',forcing_name)
703
704	IF ( TRIM(forcing_name) /= 'NONE' ) THEN
705	IF (is_root_prc) CALL SYSTEM ('rm -f '//TRIM(forcing_name))
706	WRITE(numout,*) 'writing a forcing file for STOMATE.'
707
708	!Config Key = STOMATE_FORCING_MEMSIZE
709	!Config Desc = Size of STOMATE forcing data in memory
710	!Config If = OK_STOMATE
711	!Config Def = 50
712	!Config Help = This variable determines how many
713	!Config forcing states will be kept in memory.
714	!Config Must be a compromise between memory
715	!Config use and frequeny of disk access.
716	!Config Units = [MegaBytes]
717
718	max_totsize = 50
719	CALL getin_p('STOMATE_FORCING_MEMSIZE', max_totsize)
720	max_totsize = max_totsize*1000000
721
722	totsize_1step = &
723	& SIZE(clay)*KIND(clay) &
724	& +SIZE(humrel_daily)*KIND(humrel_daily) &
725	& +SIZE(litterhum_daily)*KIND(litterhum_daily) &
726	& +SIZE(t2m_daily)*KIND(t2m_daily) &
727	& +SIZE(t2m_min_daily)*KIND(t2m_min_daily) &
728	& +SIZE(tsurf_daily)*KIND(tsurf_daily) &
729	& +SIZE(tsoil_daily)*KIND(tsoil_daily) &
730	& +SIZE(soilhum_daily)*KIND(soilhum_daily) &
731	& +SIZE(precip_daily)*KIND(precip_daily) &
732	& +SIZE(gpp_daily_x)*KIND(gpp_daily_x) &
733	& +SIZE(veget)*KIND(veget) &
734	& +SIZE(veget_max)*KIND(veget_max) &
735	& +SIZE(lai)*KIND(lai)
736
737	CALL reduce_sum(totsize_1step,totsize_tmp)
738	CALL bcast(totsize_tmp)
739	totsize_1step=totsize_tmp
740	! total number of forcing steps
741	nsft = INT(one_year/(dt_slow/one_day))
742
743	! number of forcing steps in memory
744	nsfm = MIN(nsft, &
745	& MAX(1,NINT( REAL(max_totsize,r_std) &
746	& /REAL(totsize_1step,r_std))))
747
748	CALL init_forcing (kjpindex,nsfm,nsft)
749
750	isf(:) = (/ (i,i=1,nsfm) /)
751	nf_written(:) = .FALSE.
752	nf_cumul(:) = 0
753
754	iisf = 0
755
756	!-
757	IF (is_root_prc) THEN
758	ier = NF90_CREATE (TRIM(forcing_name),NF90_SHARE,forcing_id)
759	ier = NF90_PUT_ATT (forcing_id,NF90_GLOBAL,'dtradia',dtradia)
760	ier = NF90_PUT_ATT (forcing_id,NF90_GLOBAL,'dt_slow',dt_slow)
761	ier = NF90_PUT_ATT (forcing_id,NF90_GLOBAL, &
762	& 'nsft',REAL(nsft,r_std))
763	ier = NF90_PUT_ATT (forcing_id,NF90_GLOBAL, &
764	& 'kjpij',REAL(iim_g*jjm_g,r_std))
765	ier = NF90_PUT_ATT (forcing_id,NF90_GLOBAL, &
766	& 'kjpindex',REAL(nbp_glo,r_std))
767	!-
768	ier = NF90_DEF_DIM (forcing_id,'points',nbp_glo,d_id(1))
769	ier = NF90_DEF_DIM (forcing_id,'layers',nbdl,d_id(2))
770	ier = NF90_DEF_DIM (forcing_id,'pft',nvm,d_id(3))
771	direct=2
772	ier = NF90_DEF_DIM (forcing_id,'direction',direct,d_id(4))
773	nneigh=8
774	ier = NF90_DEF_DIM (forcing_id,'nneigh',nneigh,d_id(5))
775	ier = NF90_DEF_DIM (forcing_id,'time',nsft,d_id(6))
776	ier = NF90_DEF_DIM (forcing_id,'nbparts',nparts,d_id(7))
777	!-
778	ier = NF90_DEF_VAR (forcing_id,'points', r_typ,d_id(1),vid)
779	ier = NF90_DEF_VAR (forcing_id,'layers', r_typ,d_id(2),vid)
780	ier = NF90_DEF_VAR (forcing_id,'pft', r_typ,d_id(3),vid)
781	ier = NF90_DEF_VAR (forcing_id,'direction', r_typ,d_id(4),vid)
782	ier = NF90_DEF_VAR (forcing_id,'nneigh', r_typ,d_id(5),vid)
783	ier = NF90_DEF_VAR (forcing_id,'time', r_typ,d_id(6),vid)
784	ier = NF90_DEF_VAR (forcing_id,'nbparts', r_typ,d_id(7),vid)
785	ier = NF90_DEF_VAR (forcing_id,'index', r_typ,d_id(1),vid)
786	ier = NF90_DEF_VAR (forcing_id,'contfrac', r_typ,d_id(1),vid)
787	ier = NF90_DEF_VAR (forcing_id,'lalo', &
788	& r_typ,(/ d_id(1),d_id(4) /),vid)
789	ier = NF90_DEF_VAR (forcing_id,'neighbours', &
790	& r_typ,(/ d_id(1),d_id(5) /),vid)
791	ier = NF90_DEF_VAR (forcing_id,'resolution', &
792	& r_typ,(/ d_id(1),d_id(4) /),vid)
793	ier = NF90_DEF_VAR (forcing_id,'clay', &
794	& r_typ,(/ d_id(1),d_id(6) /),vid)
795	ier = NF90_DEF_VAR (forcing_id,'humrel', &
796	& r_typ,(/ d_id(1),d_id(3),d_id(6) /),vid)
797	ier = NF90_DEF_VAR (forcing_id,'litterhum', &
798	& r_typ,(/ d_id(1),d_id(6) /),vid)
799	ier = NF90_DEF_VAR (forcing_id,'t2m', &
800	& r_typ,(/ d_id(1),d_id(6) /),vid)
801	ier = NF90_DEF_VAR (forcing_id,'t2m_min', &
802	& r_typ,(/ d_id(1),d_id(6) /),vid)
803	ier = NF90_DEF_VAR (forcing_id,'tsurf', &
804	& r_typ,(/ d_id(1),d_id(6) /),vid)
805	ier = NF90_DEF_VAR (forcing_id,'tsoil', &
806	& r_typ,(/ d_id(1),d_id(2),d_id(6) /),vid)
807	ier = NF90_DEF_VAR (forcing_id,'soilhum', &
808	& r_typ,(/ d_id(1),d_id(2),d_id(6) /),vid)
809	ier = NF90_DEF_VAR (forcing_id,'precip', &
810	& r_typ,(/ d_id(1),d_id(6) /),vid)
811	ier = NF90_DEF_VAR (forcing_id,'gpp', &
812	& r_typ,(/ d_id(1),d_id(3),d_id(6) /),vid)
813	ier = NF90_DEF_VAR (forcing_id,'veget', &
814	& r_typ,(/ d_id(1),d_id(3),d_id(6) /),vid)
815	ier = NF90_DEF_VAR (forcing_id,'veget_max', &
816	& r_typ,(/ d_id(1),d_id(3),d_id(6) /),vid)
817	ier = NF90_DEF_VAR (forcing_id,'lai', &
818	& r_typ,(/ d_id(1),d_id(3),d_id(6) /),vid)
819	ier = NF90_ENDDEF (forcing_id)
820	!-
821	ier = NF90_INQ_VARID (forcing_id,'points',vid)
822	ier = NF90_PUT_VAR (forcing_id,vid, &
823	& (/(REAL(i,r_std),i=1,nbp_glo) /))
824	ier = NF90_INQ_VARID (forcing_id,'layers',vid)
825	ier = NF90_PUT_VAR (forcing_id,vid,(/(REAL(i,r_std),i=1,nbdl)/))
826	ier = NF90_INQ_VARID (forcing_id,'pft',vid)
827	ier = NF90_PUT_VAR (forcing_id,vid,(/(REAL(i,r_std),i=1,nvm)/))
828	ier = NF90_INQ_VARID (forcing_id,'direction',vid)
829	ier = NF90_PUT_VAR (forcing_id,vid,(/(REAL(i,r_std),i=1,2)/))
830	ier = NF90_INQ_VARID (forcing_id,'nneigh',vid)
831	ier = NF90_PUT_VAR (forcing_id,vid,(/(REAL(i,r_std),i=1,8)/))
832	ier = NF90_INQ_VARID (forcing_id,'time',vid)
833	ier = NF90_PUT_VAR (forcing_id,vid,(/(REAL(i,r_std),i=1,nsft)/))
834	ier = NF90_INQ_VARID (forcing_id,'nbparts',vid)
835	ier = NF90_PUT_VAR (forcing_id,vid,(/(REAL(i,r_std),i=1,nparts)/))
836	ier = NF90_INQ_VARID (forcing_id,'index',vid)
837	ier = NF90_PUT_VAR (forcing_id,vid,REAL(index_g,r_std))
838	ier = NF90_INQ_VARID (forcing_id,'contfrac',vid)
839	ier = NF90_PUT_VAR (forcing_id,vid,REAL(contfrac_g,r_std))
840	ier = NF90_INQ_VARID (forcing_id,'lalo',vid)
841	ier = NF90_PUT_VAR (forcing_id,vid,lalo_g)
842	!ym attention a neighbours, a modifier plus tard
843	ier = NF90_INQ_VARID (forcing_id,'neighbours',vid)
844	ier = NF90_PUT_VAR (forcing_id,vid,REAL(neighbours_g,r_std))
845	ier = NF90_INQ_VARID (forcing_id,'resolution',vid)
846	ier = NF90_PUT_VAR (forcing_id,vid,resolution_g)
847	ENDIF
848	ENDIF
849	ENDIF
850	!
851	! 1.7 write forcing file for the soil?
852	!
853	IF (ldcarbon_write) THEN
854	!
855	!Config Key = STOMATE_CFORCING_NAME
856	!Config Desc = Name of STOMATE's carbon forcing file
857	!Config If = OK_STOMATE
858	!Config Def = NONE
859	!Config Help = Name that will be given to STOMATE's carbon
860	!Config offline forcing file
861	!Config Units = [FILE]
862	!-
863	Cforcing_name = stomate_Cforcing_name ! compatibilité avec driver Nicolas
864	CALL getin_p('STOMATE_CFORCING_NAME',Cforcing_name)
865
866	IF ( TRIM(Cforcing_name) /= 'NONE' ) THEN
867	IF (is_root_prc) CALL SYSTEM ('rm -f '//TRIM(Cforcing_name))
868	!
869	! time step of forcesoil
870	!
871	!Config Key = FORCESOIL_STEP_PER_YEAR
872	!Config Desc = Number of time steps per year for carbon spinup.
873	!Config If = OK_STOMATE
874	!Config Def = 365
875	!Config Help = Number of time steps per year for carbon spinup.
876	!Config Units = [days, months, year]
877	nparan = 365
878	CALL getin_p('FORCESOIL_STEP_PER_YEAR', nparan)
879
880	IF ( nparan < 1 ) nparan = 1
881
882	!Config Key = FORCESOIL_NB_YEAR
883	!Config Desc = Number of years saved for carbon spinup.
884	!Config If = OK_STOMATE
885	!Config Def = 1
886	!Config Help = Number of years saved for carbon spinup. If internal parameter cumul_Cforcing is TRUE in stomate.f90
887	!Config Then this parameter is forced to one.
888	!Config Units = [years]
889	nbyear=1
890	CALL getin_p('FORCESOIL_NB_YEAR', nbyear)
891	IF ( cumul_Cforcing ) THEN
892	CALL ipslerr (1,'stomate', &
893	& 'Internal parameter cumul_Cforcing is TRUE in stomate.f90', &
894	& 'Parameter FORCESOIL_NB_YEAR is forced to 1.', &
895	& 'Then variable nbyear is forced to 1.')
896	nbyear=1
897	ENDIF
898
899	dt_forcesoil = zero
900	nparan = nparan+1
901	DO WHILE (dt_forcesoil < dt_slow/one_day)
902	nparan = nparan-1
903	IF (nparan < 1) THEN
904	STOP 'Problem 1 with number of soil forcing time steps.'
905	ENDIF
906	dt_forcesoil = one_year/REAL(nparan,r_std)
907	ENDDO
908
909	IF ( nparan > nparanmax ) THEN
910	STOP 'Problem 2 with number of soil forcing time steps.'
911	ENDIF
912
913	WRITE(numout,*) 'time step of soil forcing (d): ',dt_forcesoil
914
915	ALLOCATE( nforce(nparan*nbyear))
916	nforce(:) = 0
917
918	ALLOCATE(control_moist(kjpindex,nlevs,nparan*nbyear))
919	ALLOCATE(npp_equil(kjpindex,nparan*nbyear))
920	ALLOCATE(npp_tot(kjpindex))
921	ALLOCATE(control_temp(kjpindex,nlevs,nparan*nbyear))
922	ALLOCATE(soilcarbon_input(kjpindex,ncarb,nvm,nparan*nbyear))
923
924	npp_equil(:,:) = zero
925	npp_tot(:) = zero
926	control_moist(:,:,:) = zero
927	control_temp(:,:,:) = zero
928	soilcarbon_input(:,:,:,:) = zero
929	ENDIF
930	ENDIF
931	!
932	! 1.8 calculate STOMATE's vegetation fractions
933	! from veget, veget_max
934	!
935	DO j=1,nvm
936	WHERE ((1.-totfrac_nobio(:)) > min_sechiba)
937	veget_cov(:,j)=veget(:,j)/( 1.-totfrac_nobio(:) )
938	veget_cov_max(:,j)=veget_max(:,j)/( 1.-totfrac_nobio(:) )
939	ELSEWHERE
940	veget_cov(:,j) = zero
941	veget_cov_max(:,j) = zero
942	ENDWHERE
943	ENDDO
944	IF ( EndOfYear ) THEN
945	DO j=1,nvm
946	WHERE ((1.-totfrac_nobio_new(:)) > min_sechiba)
947	veget_cov_max_new(:,j)=veget_max_new(:,j)/( 1.-totfrac_nobio_new(:) )
948	ELSEWHERE
949	veget_cov_max_new(:,j) = zero
950	ENDWHERE
951	ENDDO
952	ENDIF
953	!
954	! 1.9 initialize some variables
955	! STOMATE diagnoses some variables for SECHIBA :
956	! assim_param, deadleaf_cover, etc.
957	! These variables can be recalculated easily
958	! from STOMATE's prognostic variables.
959	! height is saved in Sechiba.
960	!
961	IF (control%ok_stomate) THEN
962	CALL stomate_var_init &
963	& (kjpindex, veget_cov, veget_cov_max, leaf_age, leaf_frac, &
964	& tlong_ref, t2m_month, dead_leaves, &
965	& veget, lai, qsintmax, deadleaf_cover, assim_param)
966	harvest_above = zero
967	ENDIF
968	!
969	! 1.10 reset flag
970	!
971	l_first_stomate = .FALSE.
972	!
973	! 1.11 return
974	!
975	RETURN
976	ENDIF ! first call
977	IF (bavard >= 4) THEN
978	WRITE(numout,*) 'DATE ',date,' ymds', year, month, day, sec, '-- stp --', itime, do_slow
979	ENDIF
980	!-
981	! 2 prepares restart file for the next simulation
982	!-
983	IF (ldrestart_write) THEN
984	IF (long_print) THEN
985	WRITE (numout,*) &
986	& ' we have to complete restart file with STOMATE variables'
987	ENDIF
988	CALL writerestart &
989	& (kjpindex, index, &
990	& day_counter, dt_days, date, &
991	& ind, adapted, regenerate, &
992	& humrel_daily, litterhum_daily, &
993	& t2m_daily, t2m_min_daily, tsurf_daily, tsoil_daily, &
994	& soilhum_daily, precip_daily, &
995	& gpp_daily, npp_daily, turnover_daily, &
996	& humrel_month, humrel_week, &
997	& t2m_longterm, tlong_ref, t2m_month, t2m_week, &
998	& tsoil_month, soilhum_month, fireindex, firelitter, &
999	& maxhumrel_lastyear, maxhumrel_thisyear, &
1000	& minhumrel_lastyear, minhumrel_thisyear, &
1001	& maxgppweek_lastyear, maxgppweek_thisyear, &
1002	& gdd0_lastyear, gdd0_thisyear, &
1003	& precip_lastyear, precip_thisyear, &
1004	& gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
1005	& PFTpresent, npp_longterm, lm_lastyearmax, lm_thisyearmax, &
1006	& maxfpc_lastyear, maxfpc_thisyear, &
1007	& turnover_longterm, gpp_week, biomass, resp_maint_part, &
1008	& leaf_age, leaf_frac, &
1009	& senescence, when_growthinit, age, &
1010	& resp_hetero_d, resp_maint_d, resp_growth_d, co2_fire, co2_to_bm_dgvm, &
1011	& veget_lastlight, everywhere, need_adjacent, &
1012	& RIP_time, time_lowgpp, &
1013	& time_hum_min, hum_min_dormance, &
1014	& litterpart, litter, dead_leaves, &
1015	& carbon, black_carbon, lignin_struc,turnover_time,&
1016	& prod10,prod100,flux10, flux100, &
1017	& convflux, cflux_prod10, cflux_prod100, bm_to_litter, carb_mass_total)
1018
1019	IF (ldforcing_write .AND. TRIM(forcing_name) /= 'NONE' ) THEN
1020	CALL forcing_write(forcing_id,1,iisf)
1021	!
1022	IF (is_root_prc) ier = NF90_CLOSE (forcing_id)
1023	forcing_id=-1
1024	ENDIF
1025
1026	IF (ldcarbon_write .AND. TRIM(Cforcing_name) /= 'NONE' ) THEN
1027	WRITE(numout,*) &
1028	& 'stomate: writing the forcing file for carbon spinup'
1029	!
1030	DO iatt = 1, nparan*nbyear
1031	IF ( nforce(iatt) > 0 ) THEN
1032	soilcarbon_input(:,:,:,iatt) = &
1033	& soilcarbon_input(:,:,:,iatt)/REAL(nforce(iatt),r_std)
1034	control_moist(:,:,iatt) = &
1035	& control_moist(:,:,iatt)/REAL(nforce(iatt),r_std)
1036	control_temp(:,:,iatt) = &
1037	& control_temp(:,:,iatt)/REAL(nforce(iatt),r_std)
1038	npp_equil(:,iatt) = &
1039	& npp_equil(:,iatt)/REAL(nforce(iatt),r_std)
1040	ELSE
1041	WRITE(numout,*) &
1042	& 'We have no soil carbon forcing data for this time step:', &
1043	& iatt
1044	WRITE(numout,*) ' -> we set them to zero'
1045	soilcarbon_input(:,:,:,iatt) = zero
1046	control_moist(:,:,iatt) = zero
1047	control_temp(:,:,iatt) = zero
1048	npp_equil(:,iatt) = zero
1049	ENDIF
1050	ENDDO
1051
1052	IF (is_root_prc) THEN
1053	ALLOCATE(soilcarbon_input_g(nbp_glo,ncarb,nvm,nparan*nbyear))
1054	ALLOCATE(control_moist_g(nbp_glo,nlevs,nparan*nbyear))
1055	ALLOCATE(control_temp_g(nbp_glo,nlevs,nparan*nbyear))
1056	ALLOCATE(npp_equil_g(nbp_glo,nparan*nbyear))
1057	ENDIF
1058
1059	CALL gather(clay,clay_g)
1060	CALL gather(soilcarbon_input,soilcarbon_input_g)
1061	CALL gather(control_moist,control_moist_g)
1062	CALL gather(control_temp,control_temp_g)
1063	CALL gather(npp_equil,npp_equil_g)
1064
1065	!-
1066	IF (is_root_prc) THEN
1067	ier = NF90_CREATE (TRIM(Cforcing_name),NF90_WRITE,Cforcing_id)
1068	ier = NF90_PUT_ATT (Cforcing_id,NF90_GLOBAL, &
1069	& 'kjpindex',REAL(nbp_glo,r_std))
1070	ier = NF90_PUT_ATT (Cforcing_id,NF90_GLOBAL, &
1071	& 'nparan',REAL(nparan,r_std))
1072	ier = NF90_PUT_ATT (Cforcing_id,NF90_GLOBAL, &
1073	& 'nbyear',REAL(nbyear,r_std))
1074	ier = NF90_DEF_DIM (Cforcing_id,'points',nbp_glo,d_id(1))
1075	ier = NF90_DEF_DIM (Cforcing_id,'carbtype',ncarb,d_id(2))
1076	ier = NF90_DEF_DIM (Cforcing_id,'vegtype',nvm,d_id(3))
1077	ier = NF90_DEF_DIM (Cforcing_id,'level',nlevs,d_id(4))
1078	ier = NF90_DEF_DIM (Cforcing_id,'time_step',nparan*nbyear,d_id(5))
1079	!-
1080	ier = NF90_DEF_VAR (Cforcing_id,'points', r_typ,d_id(1),vid)
1081	ier = NF90_DEF_VAR (Cforcing_id,'carbtype', r_typ,d_id(2),vid)
1082	ier = NF90_DEF_VAR (Cforcing_id,'vegtype', r_typ,d_id(3),vid)
1083	ier = NF90_DEF_VAR (Cforcing_id,'level', r_typ,d_id(4),vid)
1084	ier = NF90_DEF_VAR (Cforcing_id,'time_step', r_typ,d_id(5),vid)
1085	ier = NF90_DEF_VAR (Cforcing_id,'index', r_typ,d_id(1),vid)
1086	ier = NF90_DEF_VAR (Cforcing_id,'clay', r_typ,d_id(1),vid)
1087	ier = NF90_DEF_VAR (Cforcing_id,'soilcarbon_input',r_typ, &
1088	& (/ d_id(1),d_id(2),d_id(3),d_id(5) /),vid)
1089	ier = NF90_DEF_VAR (Cforcing_id,'control_moist',r_typ, &
1090	& (/ d_id(1),d_id(4),d_id(5) /),vid)
1091	ier = NF90_DEF_VAR (Cforcing_id,'control_temp',r_typ, &
1092	& (/ d_id(1),d_id(4),d_id(5) /),vid)
1093	ier = NF90_DEF_VAR (Cforcing_id,'npp_equil',r_typ, &
1094	& (/ d_id(1),d_id(5) /),vid)
1095	ier = NF90_ENDDEF (Cforcing_id)
1096	!-
1097	ier = NF90_INQ_VARID (Cforcing_id,'points',vid)
1098	ier = NF90_PUT_VAR (Cforcing_id,vid, &
1099	& (/(REAL(i,r_std),i=1,nbp_glo)/))
1100	ier = NF90_INQ_VARID (Cforcing_id,'carbtype',vid)
1101	ier = NF90_PUT_VAR (Cforcing_id,vid, &
1102	& (/(REAL(i,r_std),i=1,ncarb)/))
1103	ier = NF90_INQ_VARID (Cforcing_id,'vegtype',vid)
1104	ier = NF90_PUT_VAR (Cforcing_id,vid, &
1105	& (/(REAL(i,r_std),i=1,nvm)/))
1106	ier = NF90_INQ_VARID (Cforcing_id,'level',vid)
1107	ier = NF90_PUT_VAR (Cforcing_id,vid, &
1108	& (/(REAL(i,r_std),i=1,nlevs)/))
1109	ier = NF90_INQ_VARID (Cforcing_id,'time_step',vid)
1110	ier = NF90_PUT_VAR (Cforcing_id,vid, &
1111	& (/(REAL(i,r_std),i=1,nparan*nbyear)/))
1112	ier = NF90_INQ_VARID (Cforcing_id,'index',vid)
1113	ier = NF90_PUT_VAR (Cforcing_id,vid, REAL(index_g,r_std) )
1114	ier = NF90_INQ_VARID (Cforcing_id,'clay',vid)
1115	ier = NF90_PUT_VAR (Cforcing_id,vid, clay_g )
1116	ier = NF90_INQ_VARID (Cforcing_id,'soilcarbon_input',vid)
1117	ier = NF90_PUT_VAR (Cforcing_id,vid, soilcarbon_input_g )
1118	ier = NF90_INQ_VARID (Cforcing_id,'control_moist',vid)
1119	ier = NF90_PUT_VAR (Cforcing_id,vid, control_moist_g )
1120	ier = NF90_INQ_VARID (Cforcing_id,'control_temp',vid)
1121	ier = NF90_PUT_VAR (Cforcing_id,vid, control_temp_g )
1122	ier = NF90_INQ_VARID (Cforcing_id,'npp_equil',vid)
1123	ier = NF90_PUT_VAR (Cforcing_id,vid, npp_equil_g )
1124	!-
1125	ier = NF90_CLOSE (Cforcing_id)
1126	Cforcing_id = -1
1127	ENDIF
1128
1129	IF (is_root_prc) THEN
1130	DEALLOCATE(soilcarbon_input_g)
1131	DEALLOCATE(control_moist_g)
1132	DEALLOCATE(control_temp_g)
1133	DEALLOCATE(npp_equil_g)
1134	ENDIF
1135	ENDIF
1136	RETURN
1137	ENDIF ! write restart-files
1138	!
1139	! 4 Special treatment for some input arrays.
1140	!
1141	!
1142	! 4.1 Sum of liquid and solid precipitation
1143	!
1144	precip = ( precip_rain+precip_snow )*one_day/dtradia
1145	!
1146	! 4.2 Copy.
1147	!
1148	! 4.2.1 calculate STOMATE's vegetation fractions
1149	! from veget and veget_max
1150	!
1151	DO j=1,nvm
1152	WHERE ((1.-totfrac_nobio(:)) > min_sechiba)
1153	veget_cov(:,j)=veget(:,j)/( 1.-totfrac_nobio(:) )
1154	veget_cov_max(:,j)=veget_max(:,j)/( 1.-totfrac_nobio(:) )
1155	ELSEWHERE
1156	veget_cov(:,j) = zero
1157	veget_cov_max(:,j) = zero
1158	ENDWHERE
1159	ENDDO
1160	IF ( EndOfYear ) THEN
1161	DO j=1,nvm
1162	WHERE ((1.-totfrac_nobio(:)) > min_sechiba)
1163	veget_cov_max_new(:,j)=veget_max_new(:,j)/( 1.-totfrac_nobio(:) )
1164	ELSEWHERE
1165	veget_cov_max_new(:,j) = zero
1166	ENDWHERE
1167	ENDDO
1168	ENDIF
1169	gpp_d(:,1) = zero
1170	DO j = 2,nvm
1171	! 4.2.2 GPP
1172	! gpp in gC/m**2 of total ground/day
1173	WHERE (veget_cov_max(:,j) > min_stomate)
1174	gpp_d(:,j) = gpp(:,j)/ veget_cov_max(:,j)* one_day/dtradia
1175	ELSEWHERE
1176	gpp_d(:,j) = zero
1177	ENDWHERE
1178	ENDDO
1179
1180	IF ( day == 1 .AND. sec .LT. dtradia ) THEN
1181	EndOfMonth=.TRUE.
1182	ELSE
1183	EndOfMonth=.FALSE.
1184	ENDIF
1185	!
1186	! 5 "daily" variables
1187	! Note: If dt_days /= 1, then xx_daily are not daily variables,
1188	! but that is not really a problem.
1189	!
1190	!
1191	! 5.1 accumulate instantaneous variables
1192	! and eventually calculate daily mean value
1193	!
1194
1195	CALL stomate_accu (kjpindex, nvm, dt_slow, dtradia, &
1196	& do_slow, humrel, humrel_daily)
1197	CALL stomate_accu (kjpindex, 1, dt_slow, dtradia, &
1198	& do_slow, litterhumdiag, litterhum_daily)
1199	CALL stomate_accu (kjpindex, 1, dt_slow, dtradia, &
1200	& do_slow, t2m, t2m_daily)
1201	CALL stomate_accu (kjpindex, 1, dt_slow, dtradia, &
1202	& do_slow, temp_sol, tsurf_daily)
1203	CALL stomate_accu (kjpindex, nbdl, dt_slow, dtradia, &
1204	& do_slow, stempdiag, tsoil_daily)
1205	CALL stomate_accu (kjpindex, nbdl, dt_slow, dtradia, &
1206	& do_slow, shumdiag, soilhum_daily)
1207	CALL stomate_accu (kjpindex, 1, dt_slow, dtradia, &
1208	& do_slow, precip, precip_daily)
1209	CALL stomate_accu (kjpindex, nvm, dt_slow, dtradia, &
1210	& do_slow, gpp_d, gpp_daily)
1211
1212	!
1213	! 5.2 daily minimum temperature
1214	!
1215	t2m_min_daily(:) = MIN( t2m_min(:), t2m_min_daily(:) )
1216	!
1217	! 5.3 calculate respiration (NV 14/5/2002)
1218	!
1219	! 5.3.1 calculate maintenance respiration
1220	!
1221	!NV maintenant lai est passé comme argument car on n'a plus les problemes de nat/agri!
1222	!donc plus à etre calculé ici....
1223	CALL maint_respiration &
1224	& (kjpindex,dtradia,lai,t2m,tlong_ref,stempdiag,height,veget_cov_max, &
1225	& rprof,biomass,resp_maint_part_radia)
1226
1227	resp_maint_radia(:,:) = zero
1228	DO j=2,nvm
1229	DO k= 1, nparts
1230	resp_maint_radia(:,j) = resp_maint_radia(:,j) &
1231	& +resp_maint_part_radia(:,j,k)
1232	ENDDO
1233	ENDDO
1234	resp_maint_part(:,:,:)= resp_maint_part(:,:,:) &
1235	& +resp_maint_part_radia(:,:,:)
1236	!
1237	! 5.3.2 the whole litter stuff:
1238	! litter update, lignin content, PFT parts, litter decay,
1239	! litter heterotrophic respiration, dead leaf soil cover.
1240	! No vertical discretisation in the soil for litter decay.
1241	!
1242	turnover_littercalc = turnover_daily * dtradia/one_day
1243	bm_to_littercalc = bm_to_litter*dtradia/one_day
1244
1245	CALL littercalc (kjpindex, dtradia/one_day, &
1246	turnover_littercalc, bm_to_littercalc, &
1247	veget_cov_max, temp_sol, stempdiag, shumdiag, litterhumdiag, &
1248	litterpart, litter, dead_leaves, lignin_struc, &
1249	deadleaf_cover, resp_hetero_litter, &
1250	soilcarbon_input_inst, control_temp_inst, control_moist_inst)
1251	resp_hetero_litter=resp_hetero_litter*dtradia/one_day
1252	!
1253	! 5.3.3 soil carbon dynamics: heterotrophic respiration from the soil.
1254	! For the moment, no vertical discretisation.
1255	! We might later introduce a vertical discretisation.
1256	! However, in that case, we would have to treat the vertical
1257	! exchanges of carbon between the different levels.
1258	!
1259	CALL soilcarbon (kjpindex, dtradia/one_day, clay, &
1260	soilcarbon_input_inst, control_temp_inst, control_moist_inst, &
1261	carbon, resp_hetero_soil)
1262	resp_hetero_soil=resp_hetero_soil*dtradia/one_day
1263	resp_hetero_radia = resp_hetero_litter+resp_hetero_soil
1264	resp_hetero_d = resp_hetero_d + resp_hetero_radia
1265	CALL stomate_accu (kjpindex, nlevs, dt_slow, dtradia, &
1266	& do_slow, control_moist_inst, control_moist_daily)
1267	CALL stomate_accu (kjpindex, nlevs, dt_slow, dtradia, &
1268	& do_slow, control_temp_inst, control_temp_daily)
1269	DO i=1,ncarb
1270	CALL stomate_accu (kjpindex, nvm, dt_slow, dtradia, &
1271	& do_slow, soilcarbon_input_inst(:,i,:), soilcarbon_input_daily(:,i,:))
1272	ENDDO
1273	!
1274	! 6 Daily processes
1275	!
1276	IF (do_slow) THEN
1277
1278	! 6.0 update lai
1279	IF (control%ok_pheno) THEN ! lai from stomate
1280	lai(:,ibare_sechiba) = zero
1281	DO j = 2, nvm
1282	lai(:,j) = biomass(:,j,ileaf)*sla(j)
1283	ENDDO
1284	ELSE
1285	CALL setlai(kjpindex,lai) ! lai prescribed
1286	ENDIF
1287
1288	!
1289	! 6.1 total natural space
1290	!
1291	!
1292	! 6.2 Calculate longer-term "meteorological" and biological parameters
1293	!
1294	CALL season &
1295	& (kjpindex, dt_days, EndOfYear, &
1296	& veget_cov, veget_cov_max, &
1297	& humrel_daily, t2m_daily, tsoil_daily, soilhum_daily, &
1298	& precip_daily, npp_daily, biomass, &
1299	& turnover_daily, gpp_daily, when_growthinit, &
1300	& maxhumrel_lastyear, maxhumrel_thisyear, &
1301	& minhumrel_lastyear, minhumrel_thisyear, &
1302	& maxgppweek_lastyear, maxgppweek_thisyear, &
1303	& gdd0_lastyear, gdd0_thisyear, &
1304	& precip_lastyear, precip_thisyear, &
1305	& lm_lastyearmax, lm_thisyearmax, &
1306	& maxfpc_lastyear, maxfpc_thisyear, &
1307	& humrel_month, humrel_week, t2m_longterm, &
1308	& tlong_ref, t2m_month, t2m_week, tsoil_month, soilhum_month, &
1309	& npp_longterm, turnover_longterm, gpp_week, &
1310	& gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
1311	& time_lowgpp, time_hum_min, hum_min_dormance, herbivores)
1312	!
1313	! 6.3 STOMATE: allocation, phenology, etc.
1314	!
1315	IF (control%ok_stomate) THEN
1316
1317	! 6.3.1 call stomate
1318
1319	CALL StomateLpj &
1320	& (kjpindex, dt_days, EndOfYear, EndOfMonth, &
1321	& neighbours, resolution, &
1322	& clay, herbivores, &
1323	& tsurf_daily, tsoil_daily, t2m_daily, t2m_min_daily, &
1324	& litterhum_daily, soilhum_daily, &
1325	& maxhumrel_lastyear, minhumrel_lastyear, &
1326	& gdd0_lastyear, precip_lastyear, &
1327	& humrel_month, humrel_week, tlong_ref, t2m_month, t2m_week, &
1328	& tsoil_month, soilhum_month, &
1329	& gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
1330	& turnover_longterm, gpp_daily, time_lowgpp, &
1331	& time_hum_min, maxfpc_lastyear, resp_maint_part,&
1332	& PFTpresent, age, fireindex, firelitter, &
1333	& leaf_age, leaf_frac, biomass, ind, adapted, regenerate, &
1334	& senescence, when_growthinit, litterpart, litter, &
1335	& dead_leaves, carbon, black_carbon, lignin_struc, &
1336	& veget_cov_max, veget_cov, npp_longterm, lm_lastyearmax, &
1337	& veget_lastlight, everywhere, need_adjacent, RIP_time, &
1338	& lai, rprof,npp_daily, turnover_daily, turnover_time,&
1339	& control_moist_inst, control_temp_inst, soilcarbon_input_inst, &
1340	& co2_to_bm_dgvm, co2_fire, &
1341	& resp_hetero_d, resp_maint_d, resp_growth_d, &
1342	& height, deadleaf_cover, vcmax, vjmax, &
1343	& t_photo_min, t_photo_opt, t_photo_max,bm_to_litter,&
1344	& prod10, prod100, flux10, flux100, veget_cov_max_new,&
1345	& convflux, cflux_prod10, cflux_prod100, harvest_above, carb_mass_total, lcchange,&
1346	& fpc_max)
1347
1348	!
1349	! fco2_lu --> luccarb
1350	fco2_lu(:) = convflux(:) &
1351	& + cflux_prod10(:) &
1352	& + cflux_prod100(:) &
1353	& + harvest_above(:)
1354	!
1355	! 6.4 output: transform from dimension nvm to dimension nvm
1356	! Several Stomate-PFTs may correspond
1357	! to a single Sechiba-PFT (see 4.2).
1358	! We sum up the vegetation cover over these Stomate-PFTs.
1359	! Mean LAI, height, and Vmax is calculated
1360	! by ponderating with (maximum) vegetation cover.
1361	!
1362
1363	! 6.4.1 calculate veget, veget_max,
1364	! from veget_cov and veget_cov_max
1365	veget(:,:) = zero
1366	veget_max(:,:) = zero
1367	!
1368	DO j = 1, nvm
1369	veget(:,j) = veget(:,j) + &
1370	& veget_cov(:,j) * ( 1.-totfrac_nobio(:) )
1371	veget_max(:,j) = veget_max(:,j) + &
1372	& veget_cov_max(:,j) * ( 1.-totfrac_nobio(:) )
1373	ENDDO
1374	!
1375	! 6.4.2 photosynthesis parameters
1376
1377	assim_param(:,:,ivcmax) = zero
1378	assim_param(:,:,ivjmax) = zero
1379	assim_param(:,:,itmin) = zero
1380	assim_param(:,:,itopt) = zero
1381	assim_param(:,:,itmax) = zero
1382
1383	DO j = 2,nvm
1384	assim_param(:,j,ivcmax)=vcmax(:,j)
1385	ENDDO
1386
1387	DO j = 2, nvm
1388	assim_param(:,j,ivjmax)=vjmax(:,j)
1389	ENDDO
1390
1391	DO j = 2, nvm
1392	assim_param(:,j,itmin)=t_photo_min(:,j)
1393	ENDDO
1394
1395	DO j = 2, nvm
1396	assim_param(:,j,itopt)=t_photo_opt(:,j)
1397	ENDDO
1398
1399	DO j = 2, nvm
1400	assim_param(:,j,itmax)=t_photo_max(:,j)
1401	ENDDO
1402
1403	!
1404	! 6.5 update forcing variables for soil carbon
1405	!
1406	IF (ldcarbon_write .AND. TRIM(Cforcing_name) /= 'NONE') THEN
1407	npp_tot(:)=0
1408
1409	DO j=2,nvm
1410	npp_tot(:)=npp_tot(:) + npp_daily(:,j)
1411	ENDDO
1412	sf_time = MODULO(REAL(date,r_std)-1,one_year*REAL(nbyear,r_std))
1413	iatt=FLOOR(sf_time/dt_forcesoil) + 1
1414	IF (iatt == 0) iatt = iatt_old + 1
1415
1416	IF ((iatt<iatt_old) .and. (.not. cumul_Cforcing)) THEN
1417	nforce(:)=0
1418	soilcarbon_input(:,:,:,:) = zero
1419	control_moist(:,:,:) = zero
1420	control_temp(:,:,:) = zero
1421	npp_equil(:,:) = zero
1422	ENDIF
1423	iatt_old=iatt
1424
1425
1426	nforce(iatt) = nforce(iatt)+1
1427	soilcarbon_input(:,:,:,iatt) = soilcarbon_input(:,:,:,iatt)+ soilcarbon_input_daily(:,:,:)
1428	control_moist(:,:,iatt) = control_moist(:,:,iatt)+ control_moist_daily(:,:)
1429	control_temp(:,:,iatt) = control_temp(:,:,iatt)+ control_temp_daily(:,:)
1430	npp_equil(:,iatt) = npp_equil(:,iatt)+npp_tot(:)
1431	ENDIF
1432	!
1433	! 6.6 updates qsintmax
1434	!
1435	DO j=2,nvm
1436	qsintmax(:,j) = qsintcst * veget(:,j) * lai(:,j)
1437	ENDDO
1438	ENDIF
1439	!
1440	! 6.7 write forcing file?
1441	! ldforcing_write should only be .TRUE.
1442	! if STOMATE is run in coupled mode.
1443	! In stand-alone mode, the forcing file is read!
1444	!
1445	IF ( ldforcing_write .AND. TRIM(forcing_name) /= 'NONE' ) THEN
1446
1447	gpp_daily_x(:,:) = zero
1448	!gpp needs to be multiplied by coverage for forcing (see above)
1449	DO j = 2, nvm
1450	gpp_daily_x(:,j) = gpp_daily_x(:,j) + &
1451	& gpp_daily(:,j) * dt_slow / one_day * veget_cov_max(:,j)
1452	ENDDO
1453	!
1454	! bare soil moisture availability has not been treated
1455	! in STOMATE (doesn't matter)
1456	!
1457	humrel_daily(:,ibare_sechiba) = humrel(:,ibare_sechiba)
1458
1459	! next forcing step in memory
1460	iisf = iisf+1
1461
1462	! how many times have we treated this forcing state
1463	xn = REAL(nf_cumul(isf(iisf)),r_std)
1464
1465	! cumulate. be careful :
1466	! precipitation is multiplied by dt_slow/one_day
1467	IF (cumul_forcing) THEN
1468	clay_fm(:,iisf) = (xn*clay_fm(:,iisf)+clay(:))/(xn+1.)
1469	humrel_daily_fm(:,:,iisf) = &
1470	& (xn*humrel_daily_fm(:,:,iisf) + humrel_daily(:,:))/(xn+1.)
1471	litterhum_daily_fm(:,iisf) = &
1472	& (xn*litterhum_daily_fm(:,iisf)+litterhum_daily(:))/(xn+1.)
1473	t2m_daily_fm(:,iisf) = &
1474	& (xn*t2m_daily_fm(:,iisf)+t2m_daily(:))/(xn+1.)
1475	t2m_min_daily_fm(:,iisf) = &
1476	& (xn*t2m_min_daily_fm(:,iisf)+t2m_min_daily(:))/(xn+1.)
1477	tsurf_daily_fm(:,iisf) = &
1478	& (xn*tsurf_daily_fm(:,iisf)+tsurf_daily(:))/(xn+1.)
1479	tsoil_daily_fm(:,:,iisf) = &
1480	& (xn*tsoil_daily_fm(:,:,iisf)+tsoil_daily(:,:))/(xn+1.)
1481	soilhum_daily_fm(:,:,iisf) = &
1482	& (xn*soilhum_daily_fm(:,:,iisf)+soilhum_daily(:,:))/(xn+1.)
1483	precip_fm(:,iisf) = &
1484	& (xnprecip_fm(:,iisf)+precip_daily(:)dt_slow/one_day)/(xn+1.)
1485	gpp_daily_fm(:,:,iisf) = &
1486	& (xn*gpp_daily_fm(:,:,iisf) + gpp_daily_x(:,:))/(xn+1.)
1487	veget_fm(:,:,iisf) = &
1488	& (xn*veget_fm(:,:,iisf) + veget(:,:) )/(xn+1.)
1489	veget_max_fm(:,:,iisf) = &
1490	& (xn*veget_max_fm(:,:,iisf) + veget_max(:,:) )/(xn+1.)
1491	lai_fm(:,:,iisf) = &
1492	& (xn*lai_fm(:,:,iisf) + lai(:,:) )/(xn+1.)
1493	ELSE
1494	clay_fm(:,iisf) = clay(:)
1495	humrel_daily_fm(:,:,iisf) = humrel_daily(:,:)
1496	litterhum_daily_fm(:,iisf) = litterhum_daily(:)
1497	t2m_daily_fm(:,iisf) = t2m_daily(:)
1498	t2m_min_daily_fm(:,iisf) =t2m_min_daily(:)
1499	tsurf_daily_fm(:,iisf) = tsurf_daily(:)
1500	tsoil_daily_fm(:,:,iisf) =tsoil_daily(:,:)
1501	soilhum_daily_fm(:,:,iisf) =soilhum_daily(:,:)
1502	precip_fm(:,iisf) = precip_daily(:)
1503	gpp_daily_fm(:,:,iisf) =gpp_daily_x(:,:)
1504	veget_fm(:,:,iisf) = veget(:,:)
1505	veget_max_fm(:,:,iisf) =veget_max(:,:)
1506	lai_fm(:,:,iisf) =lai(:,:)
1507	ENDIF
1508	nf_cumul(isf(iisf)) = nf_cumul(isf(iisf))+1
1509
1510	! do we have to write the forcing states?
1511	IF (iisf == nsfm) THEN
1512
1513	! write these forcing states
1514	CALL forcing_write(forcing_id,1,nsfm)
1515	! determine which forcing states must be read
1516	isf(1) = isf(nsfm)+1
1517	IF ( isf(1) > nsft ) isf(1) = 1
1518	DO iisf = 2, nsfm
1519	isf(iisf) = isf(iisf-1)+1
1520	IF (isf(iisf) > nsft) isf(iisf) = 1
1521	ENDDO
1522
1523	! read them
1524	!for debug only! CALL forcing_read(forcing_id,nsfm)
1525
1526	iisf = 0
1527
1528	ENDIF
1529
1530	ENDIF
1531
1532
1533	! 6.8 compute daily co2_flux
1534
1535	! CO2 flux in gC/m**2/sec
1536	! (positive towards the atmosphere) is sum of:
1537	! 1/ heterotrophic respiration from ground
1538	! 2/ maintenance respiration from the plants
1539	! 3/ growth respiration from the plants
1540	! 4/ co2 created by fire
1541	! 5/ - co2 taken up in the DGVM to establish saplings.
1542	! 6/ - co2 taken up by photosyntyhesis
1543
1544	!NV co2_flux devient un champs dépendant du PFT
1545	co2_flux_daily(:,:)= &
1546	& resp_maint_d(:,:) + resp_growth_d(:,:) + resp_hetero_d(:,:) + &
1547	& co2_fire(:,:) - co2_to_bm_dgvm(:,:) - gpp_daily(:,:)
1548
1549	IF ( hist_id_stom_IPCC > 0 ) THEN
1550	vartmp(:)=SUM(co2_flux_dailyveget_cov_max,dim=2)/1e3/one_daycontfrac
1551	CALL histwrite (hist_id_stom_IPCC, "nep", itime, &
1552	vartmp, kjpindex, hori_index)
1553	ENDIF
1554	!
1555	co2_flux_monthly(:,:) = co2_flux_monthly(:,:) + co2_flux_daily(:,:)
1556	! Monthly Cumulative fluxes of fluc and harvest
1557	harvest_above_monthly(:) = harvest_above_monthly(:) + harvest_above(:)
1558	cflux_prod_monthly(:) = cflux_prod_monthly(:) + convflux(:) + cflux_prod10(:) + cflux_prod100(:)
1559	IF ( EndOfMonth ) THEN
1560	IF ( control%ok_stomate ) THEN
1561	CALL histwrite (hist_id_stomate, 'CO2FLUX', itime, &
1562	co2_flux_monthly, kjpindex*nvm, horipft_index)
1563	ENDIF
1564	!MM
1565	! Si on supprimer le cumul par mois,
1566	! il ne faut pas oublier cette modif resolution(:,1)resolution(:,2)contfrac(:)
1567	DO j=2, nvm
1568	co2_flux_monthly(:,j) = co2_flux_monthly(:,j)* &
1569	resolution(:,1)resolution(:,2)contfrac(:)
1570	ENDDO
1571	!NV modif calcul du net_co2_flux
1572	net_co2_flux_monthly = zero
1573	DO ji=1,kjpindex
1574	DO j=2,nvm
1575	net_co2_flux_monthly = net_co2_flux_monthly + &
1576	& co2_flux_monthly(ji,j)*veget_cov_max(ji,j)
1577	ENDDO
1578	ENDDO
1579	! Total ( land) Cumulative fluxes of fluc and harvest
1580	net_cflux_prod_monthly_sum=&
1581	& SUM(cflux_prod_monthly(:)resolution(:,1)resolution(:,2)contfrac(:))1e-15
1582	CALL reduce_sum(net_cflux_prod_monthly_sum,net_cflux_prod_monthly_tot)
1583	CALL bcast(net_cflux_prod_monthly_tot)
1584
1585	net_harvest_above_monthly_sum=&
1586	& SUM(harvest_above_monthly(:)resolution(:,1)resolution(:,2)contfrac(:))1e-15
1587	CALL reduce_sum(net_harvest_above_monthly_sum,net_harvest_above_monthly_tot)
1588	CALL bcast(net_harvest_above_monthly_tot)
1589
1590	net_co2_flux_monthly = net_co2_flux_monthly*1e-15
1591	CALL reduce_sum(net_co2_flux_monthly,net_co2_flux_monthly_sum)
1592	CALL bcast(net_co2_flux_monthly_sum)
1593
1594	WRITE(numout,9010) 'GLOBAL net_cflux_prod_monthly (Peta gC/month) = ',net_cflux_prod_monthly_tot
1595	WRITE(numout,9010) 'GLOBAL net_harvest_above_monthly (Peta gC/month) = ',net_harvest_above_monthly_tot
1596	WRITE(numout,9010) 'GLOBAL net_co2_flux_monthly (Peta gC/month) = ',net_co2_flux_monthly_sum
1597
1598	! Calculation of net biospheric production
1599	net_biosp_prod_monthly_tot = &
1600	& ( net_co2_flux_monthly_sum + net_cflux_prod_monthly_tot + net_harvest_above_monthly_tot )
1601	WRITE(numout,9010) 'GLOBAL net_biosp_prod_monthly (Peta gC/month) = ',net_biosp_prod_monthly_tot
1602
1603	9010 FORMAT(A52,F17.14)
1604	!!$ IF ( control%ok_stomate ) THEN
1605	!!$ vartmp(:)=net_co2_flux_monthly_sum
1606	!!$ CALL histwrite (hist_id_stomate, 'CO2FLUX_MONTHLY_SUM', itime, &
1607	!!$ vartmp, kjpindex, hori_index )
1608	!!$ ENDIF
1609	!!$ IF (is_root_prc) THEN
1610	!!$ OPEN( unit=39, &
1611	!!$ file="stomate_co2flux.data", &
1612	!!$ action="write", &
1613	!!$ position="append", &
1614	!!$ iostat=ios )
1615	!!$ IF ( ios /= 0 ) THEN
1616	!!$ STOP "Erreur lors de la lecture/ecriture du fichier stomate_co2flux.data"
1617	!!$ ELSE
1618	!!$ WRITE(numout,*)
1619	!!$ WRITE(numout,*) "Ecriture du fichier stomate_co2flux.data"
1620	!!$ WRITE(numout,*)
1621	!!$ END IF
1622	!!$ WRITE(39,*) net_co2_flux_monthly_sum
1623	!!$ CLOSE( unit=39 )
1624	!!$ ENDIF
1625	co2_flux_monthly(:,:) = zero
1626	harvest_above_monthly(:) = zero
1627	cflux_prod_monthly(:) = zero
1628	ENDIF
1629	!
1630	! 6.9 reset daily variables
1631	!
1632	humrel_daily(:,:) = zero
1633	litterhum_daily(:) = zero
1634	t2m_daily(:) = zero
1635	t2m_min_daily(:) = large_value
1636	tsurf_daily(:) = zero
1637	tsoil_daily(:,:) = zero
1638	soilhum_daily(:,:) = zero
1639	precip_daily(:) = zero
1640	gpp_daily(:,:) = zero
1641	resp_maint_part(:,:,:)=zero
1642	resp_hetero_d=zero
1643	IF (bavard >= 3) THEN
1644	WRITE(numout,*) 'stomate_main: daily processes done'
1645	ENDIF
1646
1647	ENDIF ! daily processes?
1648	!
1649	! 7 Outputs from Stomate
1650	! co2_flux is assigned a value only if STOMATE is activated.
1651	! Otherwise, the calling hydrological module must do this itself.
1652	!
1653	IF ( control%ok_stomate ) THEN
1654
1655	resp_maint(:,:) = resp_maint_radia(:,:)*veget_cov_max(:,:)
1656	resp_maint(:,ibare_sechiba) = zero
1657	resp_growth(:,:)= resp_growth_d(:,:)veget_cov_max(:,:)dtradia/one_day
1658	!
1659	resp_hetero(:,:) = resp_hetero_radia(:,:)*veget_cov_max(:,:)
1660
1661	! CO2 flux in gC/m**2/sec (positive towards the atmosphere) is sum of:
1662	! 1/ heterotrophic respiration from ground
1663	! 2/ maintenance respiration from the plants
1664	! 3/ growth respiration from the plants
1665	! 4/ co2 created by fire
1666	! 5/ - co2 taken up in the DGVM to establish saplings.
1667	! 6/ - co2 taken up by photosyntyhesis
1668
1669	co2_flux(:,:) = resp_hetero(:,:) + resp_maint(:,:) + resp_growth(:,:) &
1670	& + (co2_fire(:,:)-co2_to_bm_dgvm(:,:))*veget_cov_max(:,:)/one_day &
1671	& - gpp(:,:)
1672
1673	ENDIF
1674	!
1675	! 8 message
1676	!
1677	IF ( (bavard >= 2).AND.EndOfYear.AND.do_slow) THEN
1678	WRITE(numout,*) 'stomate: EndOfYear'
1679	ENDIF
1680	IF (bavard >= 4) WRITE(numout,*) 'Leaving stomate_main'
1681	IF (long_print) WRITE (numout,*) ' stomate_main done '
1682	!--------------------------
1683	END SUBROUTINE stomate_main
1684	!
1685	!=
1686	!
1687	SUBROUTINE stomate_init &
1688	& (kjpij, kjpindex, index, ldforcing_write, lalo, &
1689	& rest_id_stom, hist_id_stom, hist_id_stom_IPCC)
1690	!---------------------------------------------------------------------
1691	! interface description
1692	! input scalar
1693	! Total size of the un-compressed grid
1694	INTEGER(i_std),INTENT(in) :: kjpij
1695	! Domain size
1696	INTEGER(i_std),INTENT(in) :: kjpindex
1697	! Logical for _forcing_ file to write
1698	LOGICAL,INTENT(in) :: ldforcing_write
1699	! Geogr. coordinates (latitude,longitude) (degrees)
1700	REAL(r_std),DIMENSION(kjpindex,2),INTENT(in) :: lalo
1701	! STOMATE's _Restart_ file file identifier
1702	INTEGER(i_std),INTENT(in) :: rest_id_stom
1703	! STOMATE's _history_ file file identifier
1704	INTEGER(i_std),INTENT(in) :: hist_id_stom
1705	! STOMATE's IPCC _history_ file file identifier
1706	INTEGER(i_std),INTENT(in) :: hist_id_stom_IPCC
1707	! input fields
1708	! Indeces of the points on the map
1709	INTEGER(i_std),DIMENSION(kjpindex),INTENT(in) :: index
1710	! local declaration
1711	REAL(r_std) :: tmp_day(1)
1712	! soil depth taken for canopy
1713	REAL(r_std) :: zcanop
1714	! soil depths at diagnostic levels
1715	REAL(r_std),DIMENSION(nbdl) :: zsoil
1716	! Index
1717	INTEGER(i_std) :: l
1718	! allocation error
1719	LOGICAL :: l_error
1720	INTEGER(i_std) :: ier
1721	! indices
1722	INTEGER(i_std) :: ji,j,ipd
1723	!---------------------------------------------------------------------
1724	!
1725	! 1 online diagnostics
1726	! (by default, "bavard" is set to 1 in stomate_constants)
1727	!
1728	!Config Key = BAVARD
1729	!Config Desc = level of online diagnostics in STOMATE (0-4)
1730	!Config If = OK_STOMATE
1731	!Config Def = 1
1732	!Config Help = With this variable, you can determine
1733	!Config how much online information STOMATE
1734	!Config gives during the run. 0 means
1735	!Config virtually no info.
1736	!Config Units = [-]
1737	!
1738	bavard = 1
1739	CALL getin_p('BAVARD', bavard)
1740
1741	IF ( kjpindex > 0 ) THEN
1742	!
1743	!Config Key = STOMATE_DIAGPT
1744	!Config Desc = Index of grid point for online diagnostics
1745	!Config If = OK_STOMATE
1746	!Config Def = 1
1747	!Config Help = This is the index of the grid point which
1748	! will be used for online diagnostics.
1749	!Config Units = [-]
1750	ipd = 1
1751	CALL getin_p('STOMATE_DIAGPT',ipd)
1752	ipd = MIN( ipd, kjpindex )
1753	WRITE(numout,*) 'Stomate: '
1754	WRITE(numout,*) ' Index of grid point for online diagnostics: ',ipd
1755	WRITE(numout,*) ' Lon, lat:',lalo(ipd,2),lalo(ipd,1)
1756	WRITE(numout,*) ' Index of this point on GCM grid: ',index(ipd)
1757	!
1758	ENDIF
1759	!
1760	! 2 check consistency of flags
1761	!
1762	IF ( ( .NOT. control%ok_stomate ) .AND. control%ok_dgvm ) THEN
1763	WRITE(numout,*) 'Cannot do dynamical vegetation without STOMATE.'
1764	WRITE(numout,*) 'We stop.'
1765	STOP
1766	ENDIF
1767
1768	IF ((.NOT.control%ok_co2).AND.control%ok_stomate) THEN
1769	WRITE(numout,*) 'Cannot call STOMATE without GPP.'
1770	WRITE(numout,*) 'We stop.'
1771	STOP
1772	ENDIF
1773
1774	IF ( ( .NOT. control%ok_co2 ) .AND. ldforcing_write ) THEN
1775	WRITE(numout,*) &
1776	& 'Cannot write forcing file if photosynthesis is not activated'
1777	WRITE(numout,*) 'We stop.'
1778	STOP
1779	ENDIF
1780	!
1781	! 3 messages
1782	!
1783	WRITE(numout,*) 'stomate: first call'
1784	WRITE(numout,*) ' Photosynthesis: ', control%ok_co2
1785	WRITE(numout,*) ' STOMATE: ', control%ok_stomate
1786	WRITE(numout,*) ' LPJ: ', control%ok_dgvm
1787	!
1788	! 4 allocation of STOMATE's variables
1789	!
1790	l_error = .FALSE.
1791	ALLOCATE(veget_cov_max(kjpindex,nvm),stat=ier)
1792	l_error = l_error .OR. (ier /= 0)
1793	ALLOCATE(veget_cov_max_new(kjpindex,nvm),stat=ier)
1794	l_error = l_error .OR. (ier /= 0)
1795	ALLOCATE(ind(kjpindex,nvm),stat=ier)
1796	l_error = l_error .OR. (ier /= 0)
1797	ALLOCATE(adapted(kjpindex,nvm),stat=ier)
1798	l_error = l_error .OR. (ier /= 0)
1799	ALLOCATE(regenerate(kjpindex,nvm),stat=ier)
1800	l_error = l_error .OR. (ier /= 0)
1801	ALLOCATE(humrel_daily(kjpindex,nvm),stat=ier)
1802	l_error = l_error .OR. (ier /= 0)
1803	ALLOCATE(litterhum_daily(kjpindex),stat=ier)
1804	l_error = l_error .OR. (ier /= 0)
1805	ALLOCATE(t2m_daily(kjpindex),stat=ier)
1806	l_error = l_error .OR. (ier /= 0)
1807	ALLOCATE(t2m_min_daily(kjpindex),stat=ier)
1808	l_error = l_error .OR. (ier /= 0)
1809	ALLOCATE(tsurf_daily(kjpindex),stat=ier)
1810	l_error = l_error .OR. (ier /= 0)
1811	ALLOCATE(tsoil_daily(kjpindex,nbdl),stat=ier)
1812	l_error = l_error .OR. (ier /= 0)
1813	ALLOCATE(soilhum_daily(kjpindex,nbdl),stat=ier)
1814	l_error = l_error .OR. (ier /= 0)
1815	ALLOCATE(precip_daily(kjpindex),stat=ier)
1816	l_error = l_error .OR. (ier /= 0)
1817	ALLOCATE(gpp_daily(kjpindex,nvm),stat=ier)
1818	l_error = l_error .OR. (ier /= 0)
1819	ALLOCATE(npp_daily(kjpindex,nvm),stat=ier)
1820	l_error = l_error .OR. (ier /= 0)
1821	ALLOCATE(turnover_daily(kjpindex,nvm,nparts),stat=ier)
1822	l_error = l_error .OR. (ier /= 0)
1823	ALLOCATE(turnover_littercalc(kjpindex,nvm,nparts),stat=ier)
1824	l_error = l_error .OR. (ier /= 0)
1825	ALLOCATE(humrel_month(kjpindex,nvm),stat=ier)
1826	l_error = l_error .OR. (ier /= 0)
1827	ALLOCATE(humrel_week(kjpindex,nvm),stat=ier)
1828	l_error = l_error .OR. (ier /= 0)
1829	ALLOCATE(t2m_longterm(kjpindex),stat=ier)
1830	l_error = l_error .OR. (ier /= 0)
1831	ALLOCATE(tlong_ref(kjpindex),stat=ier)
1832	l_error = l_error .OR. (ier /= 0)
1833	ALLOCATE(t2m_month(kjpindex),stat=ier)
1834	l_error = l_error .OR. (ier /= 0)
1835	ALLOCATE(t2m_week(kjpindex),stat=ier)
1836	l_error = l_error .OR. (ier /= 0)
1837	ALLOCATE(tsoil_month(kjpindex,nbdl),stat=ier)
1838	l_error = l_error .OR. (ier /= 0)
1839	ALLOCATE(soilhum_month(kjpindex,nbdl),stat=ier)
1840	l_error = l_error .OR. (ier /= 0)
1841	ALLOCATE(fireindex(kjpindex,nvm),stat=ier)
1842	l_error = l_error .OR. (ier /= 0)
1843	ALLOCATE(firelitter(kjpindex,nvm),stat=ier)
1844	l_error = l_error .OR. (ier /= 0)
1845	ALLOCATE(maxhumrel_lastyear(kjpindex,nvm),stat=ier)
1846	l_error = l_error .OR. (ier /= 0)
1847	ALLOCATE(maxhumrel_thisyear(kjpindex,nvm),stat=ier)
1848	l_error = l_error .OR. (ier /= 0)
1849	ALLOCATE(minhumrel_lastyear(kjpindex,nvm),stat=ier)
1850	l_error = l_error .OR. (ier /= 0)
1851	ALLOCATE(minhumrel_thisyear(kjpindex,nvm),stat=ier)
1852	l_error = l_error .OR. (ier /= 0)
1853	ALLOCATE(maxgppweek_lastyear(kjpindex,nvm),stat=ier)
1854	l_error = l_error .OR. (ier /= 0)
1855	ALLOCATE(maxgppweek_thisyear(kjpindex,nvm),stat=ier)
1856	l_error = l_error .OR. (ier /= 0)
1857	ALLOCATE(gdd0_lastyear(kjpindex),stat=ier)
1858	l_error = l_error .OR. (ier /= 0)
1859	ALLOCATE(gdd0_thisyear(kjpindex),stat=ier)
1860	l_error = l_error .OR. (ier /= 0)
1861	ALLOCATE(precip_lastyear(kjpindex),stat=ier)
1862	l_error = l_error .OR. (ier /= 0)
1863	ALLOCATE(precip_thisyear(kjpindex),stat=ier)
1864	l_error = l_error .OR. (ier /= 0)
1865	ALLOCATE(gdd_m5_dormance(kjpindex,nvm),stat=ier)
1866	l_error = l_error .OR. (ier /= 0)
1867	ALLOCATE(gdd_midwinter(kjpindex,nvm),stat=ier)
1868	l_error = l_error .OR. (ier /= 0)
1869	ALLOCATE(ncd_dormance(kjpindex,nvm),stat=ier)
1870	l_error = l_error .OR. (ier /= 0)
1871	ALLOCATE(ngd_minus5(kjpindex,nvm),stat=ier)
1872	l_error = l_error .OR. (ier /= 0)
1873	ALLOCATE(PFTpresent(kjpindex,nvm),stat=ier)
1874	l_error = l_error .OR. (ier /= 0)
1875	ALLOCATE(npp_longterm(kjpindex,nvm),stat=ier)
1876	l_error = l_error .OR. (ier /= 0)
1877	ALLOCATE(lm_lastyearmax(kjpindex,nvm),stat=ier)
1878	l_error = l_error .OR. (ier /= 0)
1879	ALLOCATE(lm_thisyearmax(kjpindex,nvm),stat=ier)
1880	l_error = l_error .OR. (ier /= 0)
1881	ALLOCATE(maxfpc_lastyear(kjpindex,nvm),stat=ier)
1882	l_error = l_error .OR. (ier /= 0)
1883	ALLOCATE(maxfpc_thisyear(kjpindex,nvm),stat=ier)
1884	l_error = l_error .OR. (ier /= 0)
1885	ALLOCATE(turnover_longterm(kjpindex,nvm,nparts),stat=ier)
1886	l_error = l_error .OR. (ier /= 0)
1887	ALLOCATE(gpp_week(kjpindex,nvm),stat=ier)
1888	l_error = l_error .OR. (ier /= 0)
1889	ALLOCATE(biomass(kjpindex,nvm,nparts),stat=ier)
1890	l_error = l_error .OR. (ier /= 0)
1891	ALLOCATE(senescence(kjpindex,nvm),stat=ier)
1892	l_error = l_error .OR. (ier /= 0)
1893	ALLOCATE(when_growthinit(kjpindex,nvm),stat=ier)
1894	l_error = l_error .OR. (ier /= 0)
1895	ALLOCATE(age(kjpindex,nvm),stat=ier)
1896	l_error = l_error .OR. (ier /= 0)
1897	ALLOCATE(resp_hetero_d(kjpindex,nvm),stat=ier)
1898	l_error = l_error .OR. (ier /= 0)
1899	ALLOCATE(resp_hetero_radia(kjpindex,nvm),stat=ier)
1900	l_error = l_error .OR. (ier /= 0)
1901	ALLOCATE(resp_maint_d(kjpindex,nvm),stat=ier)
1902	l_error = l_error .OR. (ier /= 0)
1903	ALLOCATE(resp_growth_d(kjpindex,nvm),stat=ier)
1904	l_error = l_error .OR. (ier /= 0)
1905	!NV devient 2D
1906	ALLOCATE(co2_fire(kjpindex,nvm),stat=ier)
1907	l_error = l_error .OR. (ier /= 0)
1908	!NV devient 2D
1909	ALLOCATE(co2_to_bm_dgvm(kjpindex,nvm),stat=ier)
1910	l_error = l_error .OR. (ier /= 0)
1911	ALLOCATE(veget_lastlight(kjpindex,nvm),stat=ier)
1912	l_error = l_error .OR. (ier /= 0)
1913	ALLOCATE(everywhere(kjpindex,nvm),stat=ier)
1914	l_error = l_error .OR. (ier /= 0)
1915	ALLOCATE(need_adjacent(kjpindex,nvm),stat=ier)
1916	l_error = l_error .OR. (ier /= 0)
1917	ALLOCATE(leaf_age(kjpindex,nvm,nleafages),stat=ier)
1918	l_error = l_error .OR. (ier /= 0)
1919	ALLOCATE(leaf_frac(kjpindex,nvm,nleafages),stat=ier)
1920	l_error = l_error .OR. (ier /= 0)
1921	ALLOCATE(RIP_time(kjpindex,nvm),stat=ier)
1922	l_error = l_error .OR. (ier /= 0)
1923	ALLOCATE(time_lowgpp(kjpindex,nvm),stat=ier)
1924	l_error = l_error .OR. (ier /= 0)
1925	ALLOCATE(time_hum_min(kjpindex,nvm),stat=ier)
1926	l_error = l_error .OR. (ier /= 0)
1927	ALLOCATE(hum_min_dormance(kjpindex,nvm),stat=ier)
1928	l_error = l_error .OR. (ier /= 0)
1929	ALLOCATE(litterpart(kjpindex,nvm,nlitt),stat=ier)
1930	l_error = l_error .OR. (ier /= 0)
1931	ALLOCATE(litter(kjpindex,nlitt,nvm,nlevs),stat=ier)
1932	l_error = l_error .OR. (ier /= 0)
1933	ALLOCATE(dead_leaves(kjpindex,nvm,nlitt),stat=ier)
1934	l_error = l_error .OR. (ier /= 0)
1935	ALLOCATE(carbon(kjpindex,ncarb,nvm),stat=ier)
1936	l_error = l_error .OR. (ier /= 0)
1937	ALLOCATE(black_carbon(kjpindex),stat=ier)
1938	l_error = l_error .OR. (ier /= 0)
1939	ALLOCATE(lignin_struc(kjpindex,nvm,nlevs),stat=ier)
1940	l_error = l_error .OR. (ier /= 0)
1941	ALLOCATE(turnover_time(kjpindex,nvm),stat=ier)
1942	l_error = l_error .OR. (ier /= 0)
1943	ALLOCATE(co2_flux_daily(kjpindex,nvm),stat=ier)
1944	l_error = l_error .OR. (ier /= 0)
1945	ALLOCATE(co2_flux_monthly(kjpindex,nvm),stat=ier)
1946	l_error = l_error .OR. (ier /= 0)
1947	ALLOCATE (cflux_prod_monthly(kjpindex), stat=ier)
1948	l_error = l_error .OR. (ier.NE.0)
1949	ALLOCATE (harvest_above_monthly(kjpindex), stat=ier)
1950	l_error = l_error .OR. (ier.NE.0)
1951	ALLOCATE(bm_to_litter(kjpindex,nvm,nparts),stat=ier)
1952	l_error = l_error .OR. (ier /= 0)
1953	ALLOCATE(bm_to_littercalc(kjpindex,nvm,nparts),stat=ier)
1954	l_error = l_error .OR. (ier /= 0)
1955	ALLOCATE(herbivores(kjpindex,nvm),stat=ier)
1956	l_error = l_error .OR. (ier /= 0)
1957	ALLOCATE(hori_index(kjpindex),stat=ier)
1958	l_error = l_error .OR. (ier /= 0)
1959	ALLOCATE(horipft_index(kjpindex*nvm),stat=ier)
1960	l_error = l_error .OR. (ier /= 0)
1961	ALLOCATE(resp_maint_part_radia(kjpindex,nvm,nparts),stat=ier)
1962	l_error = l_error .OR. (ier /= 0)
1963	ALLOCATE(resp_maint_radia(kjpindex,nvm),stat=ier)
1964	l_error = l_error .OR. (ier /= 0)
1965	ALLOCATE(resp_maint_part(kjpindex,nvm,nparts),stat=ier)
1966	l_error = l_error .OR. (ier /= 0)
1967	resp_maint_part(:,:,:)=zero
1968	ALLOCATE (horip10_index(kjpindex*10), stat=ier)
1969	l_error = l_error .OR. (ier.NE.0)
1970	ALLOCATE (horip100_index(kjpindex*100), stat=ier)
1971	l_error = l_error .OR. (ier.NE.0)
1972	ALLOCATE (horip11_index(kjpindex*11), stat=ier)
1973	l_error = l_error .OR. (ier.NE.0)
1974	ALLOCATE (horip101_index(kjpindex*101), stat=ier)
1975	l_error = l_error .OR. (ier.NE.0)
1976	ALLOCATE (prod10(kjpindex,0:10), stat=ier)
1977	l_error = l_error .OR. (ier.NE.0)
1978	ALLOCATE (prod100(kjpindex,0:100), stat=ier)
1979	l_error = l_error .OR. (ier.NE.0)
1980	ALLOCATE (flux10(kjpindex,10), stat=ier)
1981	l_error = l_error .OR. (ier.NE.0)
1982	ALLOCATE (flux100(kjpindex,100), stat=ier)
1983	l_error = l_error .OR. (ier.NE.0)
1984	ALLOCATE (convflux(kjpindex), stat=ier)
1985	l_error = l_error .OR. (ier.NE.0)
1986	ALLOCATE (cflux_prod10(kjpindex), stat=ier)
1987	l_error = l_error .OR. (ier.NE.0)
1988	ALLOCATE (cflux_prod100(kjpindex), stat=ier)
1989	l_error = l_error .OR. (ier.NE.0)
1990	ALLOCATE (harvest_above(kjpindex), stat=ier)
1991	l_error = l_error .OR. (ier.NE.0)
1992	ALLOCATE (carb_mass_total(kjpindex), stat=ier)
1993	l_error = l_error .OR. (ier.NE.0)
1994	ALLOCATE (soilcarbon_input_daily(kjpindex,ncarb,nvm), stat=ier)
1995	l_error = l_error .OR. (ier.NE.0)
1996	ALLOCATE (control_temp_daily(kjpindex,nlevs), stat=ier)
1997	l_error = l_error .OR. (ier.NE.0)
1998	ALLOCATE (control_moist_daily(kjpindex,nlevs), stat=ier)
1999	l_error = l_error .OR. (ier.NE.0)
2000	!
2001	ALLOCATE (fpc_max(kjpindex,nvm), stat=ier)
2002	l_error = l_error .OR. (ier.NE.0)
2003	!
2004	IF (l_error) THEN
2005	STOP 'stomate_init: error in memory allocation'
2006	ENDIF
2007	!
2008	! 5 file definitions: stored in common variables
2009	!
2010	hist_id_stomate = hist_id_stom
2011	hist_id_stomate_IPCC = hist_id_stom_IPCC
2012	rest_id_stomate = rest_id_stom
2013	hori_index(:) = index(:)
2014
2015	! Get the indexing table for the vegetation fields.
2016	! In STOMATE we work on
2017	! reduced grids but to store in the full 3D filed vegetation variable
2018	! we need another index table : indexpft
2019
2020	DO j = 1, nvm
2021	DO ji = 1, kjpindex
2022	horipft_index((j-1)kjpindex+ji) = index(ji)+(j-1)kjpij
2023	ENDDO
2024	ENDDO
2025
2026	! indexing tables added for land cover change fields
2027	DO j = 1, 10
2028	DO ji = 1, kjpindex
2029	horip10_index((j-1)kjpindex+ji) = index(ji)+(j-1)kjpij
2030	ENDDO
2031	ENDDO
2032
2033	DO j = 1, 100
2034	DO ji = 1, kjpindex
2035	horip100_index((j-1)kjpindex+ji) = index(ji)+(j-1)kjpij
2036	ENDDO
2037	ENDDO
2038
2039	DO j = 1, 11
2040	DO ji = 1, kjpindex
2041	horip11_index((j-1)kjpindex+ji) = index(ji)+(j-1)kjpij
2042	ENDDO
2043	ENDDO
2044
2045	DO j = 1, 101
2046	DO ji = 1, kjpindex
2047	horip101_index((j-1)kjpindex+ji) = index(ji)+(j-1)kjpij
2048	ENDDO
2049	ENDDO
2050
2051	!
2052	! 6 some global initializations
2053	!
2054	! edit shilong
2055	! bm_to_litter(:,:,:) = zero
2056	turnover_daily(:,:,:) = zero
2057	resp_hetero_d(:,:) = zero
2058	co2_flux_daily(:,:) = zero
2059	co2_flux_monthly(:,:) = zero
2060	cflux_prod_monthly(:) = zero
2061	harvest_above_monthly(:) = zero
2062	control_moist_daily(:,:) = zero
2063	control_temp_daily(:,:) = zero
2064	soilcarbon_input_daily(:,:,:) = zero
2065
2066	! initialisation of land cover change variables
2067	prod10(:,:) = zero
2068	prod100(:,:) = zero
2069	flux10(:,:) = zero
2070	flux100(:,:) = zero
2071	convflux(:) = zero
2072	cflux_prod10(:) = zero
2073	cflux_prod100(:)= zero
2074
2075	fpc_max(:,:)=zero
2076	!--------------------------
2077	END SUBROUTINE stomate_init
2078	!
2079	!=
2080	!
2081	SUBROUTINE stomate_clear
2082	!---------------------------------------------------------------------
2083	! 1. Deallocate all dynamics variables
2084	IF (ALLOCATED(veget_cov_max)) DEALLOCATE(veget_cov_max)
2085	IF (ALLOCATED(veget_cov_max_new)) DEALLOCATE(veget_cov_max_new)
2086	IF (ALLOCATED(ind)) DEALLOCATE(ind)
2087	IF (ALLOCATED(adapted)) DEALLOCATE(adapted)
2088	IF (ALLOCATED(regenerate)) DEALLOCATE(regenerate)
2089	IF (ALLOCATED(humrel_daily)) DEALLOCATE(humrel_daily)
2090	IF (ALLOCATED(litterhum_daily)) DEALLOCATE(litterhum_daily)
2091	IF (ALLOCATED(t2m_daily)) DEALLOCATE(t2m_daily)
2092	IF (ALLOCATED(t2m_min_daily)) DEALLOCATE(t2m_min_daily)
2093	IF (ALLOCATED(tsurf_daily)) DEALLOCATE(tsurf_daily)
2094	IF (ALLOCATED(tsoil_daily)) DEALLOCATE(tsoil_daily)
2095	IF (ALLOCATED(soilhum_daily)) DEALLOCATE(soilhum_daily)
2096	IF (ALLOCATED(precip_daily)) DEALLOCATE(precip_daily)
2097	IF (ALLOCATED(gpp_daily)) DEALLOCATE(gpp_daily)
2098	IF (ALLOCATED(npp_daily)) DEALLOCATE(npp_daily)
2099	IF (ALLOCATED(turnover_daily)) DEALLOCATE(turnover_daily)
2100	IF (ALLOCATED(turnover_littercalc)) DEALLOCATE(turnover_littercalc)
2101	IF (ALLOCATED(humrel_month)) DEALLOCATE(humrel_month)
2102	IF (ALLOCATED(humrel_week)) DEALLOCATE(humrel_week)
2103	IF (ALLOCATED(t2m_longterm)) DEALLOCATE(t2m_longterm)
2104	IF (ALLOCATED(tlong_ref)) DEALLOCATE(tlong_ref)
2105	IF (ALLOCATED(t2m_month)) DEALLOCATE(t2m_month)
2106	IF (ALLOCATED(t2m_week)) DEALLOCATE(t2m_week)
2107	IF (ALLOCATED(tsoil_month)) DEALLOCATE(tsoil_month)
2108	IF (ALLOCATED(soilhum_month)) DEALLOCATE(soilhum_month)
2109	IF (ALLOCATED(fireindex)) DEALLOCATE(fireindex)
2110	IF (ALLOCATED(firelitter)) DEALLOCATE(firelitter)
2111	IF (ALLOCATED(maxhumrel_lastyear)) DEALLOCATE(maxhumrel_lastyear)
2112	IF (ALLOCATED(maxhumrel_thisyear)) DEALLOCATE(maxhumrel_thisyear)
2113	IF (ALLOCATED(minhumrel_lastyear)) DEALLOCATE(minhumrel_lastyear)
2114	IF (ALLOCATED(minhumrel_thisyear)) DEALLOCATE(minhumrel_thisyear)
2115	IF (ALLOCATED(maxgppweek_lastyear)) DEALLOCATE(maxgppweek_lastyear)
2116	IF (ALLOCATED(maxgppweek_thisyear)) DEALLOCATE(maxgppweek_thisyear)
2117	IF (ALLOCATED(gdd0_lastyear)) DEALLOCATE(gdd0_lastyear)
2118	IF (ALLOCATED(gdd0_thisyear)) DEALLOCATE(gdd0_thisyear)
2119	IF (ALLOCATED(precip_lastyear)) DEALLOCATE(precip_lastyear)
2120	IF (ALLOCATED(precip_thisyear)) DEALLOCATE(precip_thisyear)
2121	IF (ALLOCATED(gdd_m5_dormance)) DEALLOCATE(gdd_m5_dormance)
2122	IF (ALLOCATED(gdd_midwinter)) DEALLOCATE(gdd_midwinter)
2123	IF (ALLOCATED(ncd_dormance)) DEALLOCATE(ncd_dormance)
2124	IF (ALLOCATED(ngd_minus5)) DEALLOCATE(ngd_minus5)
2125	IF (ALLOCATED(PFTpresent)) DEALLOCATE(PFTpresent)
2126	IF (ALLOCATED(npp_longterm)) DEALLOCATE(npp_longterm)
2127	IF (ALLOCATED(lm_lastyearmax)) DEALLOCATE(lm_lastyearmax)
2128	IF (ALLOCATED(lm_thisyearmax)) DEALLOCATE(lm_thisyearmax)
2129	IF (ALLOCATED(maxfpc_lastyear)) DEALLOCATE(maxfpc_lastyear)
2130	IF (ALLOCATED(maxfpc_thisyear)) DEALLOCATE(maxfpc_thisyear)
2131	IF (ALLOCATED(turnover_longterm)) DEALLOCATE(turnover_longterm)
2132	IF (ALLOCATED(gpp_week)) DEALLOCATE(gpp_week)
2133	IF (ALLOCATED(biomass)) DEALLOCATE(biomass)
2134	IF (ALLOCATED(senescence)) DEALLOCATE(senescence)
2135	IF (ALLOCATED(when_growthinit)) DEALLOCATE(when_growthinit)
2136	IF (ALLOCATED(age)) DEALLOCATE(age)
2137	IF (ALLOCATED(resp_hetero_d)) DEALLOCATE(resp_hetero_d)
2138	IF (ALLOCATED(resp_hetero_radia)) DEALLOCATE(resp_hetero_radia)
2139	IF (ALLOCATED(resp_maint_d)) DEALLOCATE(resp_maint_d)
2140	IF (ALLOCATED(resp_growth_d)) DEALLOCATE(resp_growth_d)
2141	IF (ALLOCATED(co2_fire)) DEALLOCATE(co2_fire)
2142	IF (ALLOCATED(co2_to_bm_dgvm)) DEALLOCATE(co2_to_bm_dgvm)
2143	IF (ALLOCATED(veget_lastlight)) DEALLOCATE(veget_lastlight)
2144	IF (ALLOCATED(everywhere)) DEALLOCATE(everywhere)
2145	IF (ALLOCATED(need_adjacent)) DEALLOCATE(need_adjacent)
2146	IF (ALLOCATED(leaf_age)) DEALLOCATE(leaf_age)
2147	IF (ALLOCATED(leaf_frac)) DEALLOCATE(leaf_frac)
2148	IF (ALLOCATED(RIP_time)) DEALLOCATE(RIP_time)
2149	IF (ALLOCATED(time_lowgpp)) DEALLOCATE(time_lowgpp)
2150	IF (ALLOCATED(time_hum_min)) DEALLOCATE(time_hum_min)
2151	IF (ALLOCATED(hum_min_dormance)) DEALLOCATE(hum_min_dormance)
2152	IF (ALLOCATED(litterpart)) DEALLOCATE(litterpart)
2153	IF (ALLOCATED(litter)) DEALLOCATE(litter)
2154	IF (ALLOCATED(dead_leaves)) DEALLOCATE(dead_leaves)
2155	IF (ALLOCATED(carbon)) DEALLOCATE(carbon)
2156	IF (ALLOCATED(black_carbon)) DEALLOCATE(black_carbon)
2157	IF (ALLOCATED(lignin_struc)) DEALLOCATE(lignin_struc)
2158	IF (ALLOCATED(turnover_time)) DEALLOCATE(turnover_time)
2159	IF (ALLOCATED(co2_flux_daily)) DEALLOCATE(co2_flux_daily)
2160	IF (ALLOCATED(co2_flux_monthly)) DEALLOCATE(co2_flux_monthly)
2161	IF (ALLOCATED(harvest_above_monthly)) DEALLOCATE (harvest_above_monthly)
2162	IF (ALLOCATED(cflux_prod_monthly)) DEALLOCATE (cflux_prod_monthly)
2163	IF (ALLOCATED(bm_to_litter)) DEALLOCATE(bm_to_litter)
2164	IF (ALLOCATED(bm_to_littercalc)) DEALLOCATE(bm_to_littercalc)
2165	IF (ALLOCATED(herbivores)) DEALLOCATE(herbivores)
2166	IF (ALLOCATED(resp_maint_part_radia)) DEALLOCATE(resp_maint_part_radia)
2167	IF (ALLOCATED(resp_maint_radia)) DEALLOCATE(resp_maint_radia)
2168	IF (ALLOCATED(resp_maint_part)) DEALLOCATE(resp_maint_part)
2169	IF (ALLOCATED(hori_index)) DEALLOCATE(hori_index)
2170	IF (ALLOCATED(horipft_index)) DEALLOCATE(horipft_index)
2171	IF (ALLOCATED(clay_fm)) DEALLOCATE(clay_fm)
2172	IF (ALLOCATED(humrel_daily_fm)) DEALLOCATE(humrel_daily_fm)
2173	IF (ALLOCATED(litterhum_daily_fm)) DEALLOCATE(litterhum_daily_fm)
2174	IF (ALLOCATED(t2m_daily_fm)) DEALLOCATE(t2m_daily_fm)
2175	IF (ALLOCATED(t2m_min_daily_fm)) DEALLOCATE(t2m_min_daily_fm)
2176	IF (ALLOCATED(tsurf_daily_fm)) DEALLOCATE(tsurf_daily_fm)
2177	IF (ALLOCATED(tsoil_daily_fm)) DEALLOCATE(tsoil_daily_fm)
2178	IF (ALLOCATED(soilhum_daily_fm)) DEALLOCATE(soilhum_daily_fm)
2179	IF (ALLOCATED(precip_fm)) DEALLOCATE(precip_fm)
2180	IF (ALLOCATED(gpp_daily_fm)) DEALLOCATE(gpp_daily_fm)
2181	IF (ALLOCATED(veget_fm)) DEALLOCATE(veget_fm)
2182	IF (ALLOCATED(veget_max_fm)) DEALLOCATE(veget_max_fm)
2183	IF (ALLOCATED(lai_fm)) DEALLOCATE(lai_fm)
2184
2185	IF (is_root_prc) THEN
2186	IF (ALLOCATED(clay_fm_g)) DEALLOCATE(clay_fm_g)
2187	IF (ALLOCATED(humrel_daily_fm_g)) DEALLOCATE(humrel_daily_fm_g)
2188	IF (ALLOCATED(litterhum_daily_fm_g)) DEALLOCATE(litterhum_daily_fm_g)
2189	IF (ALLOCATED(t2m_daily_fm_g)) DEALLOCATE(t2m_daily_fm_g)
2190	IF (ALLOCATED(t2m_min_daily_fm_g)) DEALLOCATE(t2m_min_daily_fm_g)
2191	IF (ALLOCATED(tsurf_daily_fm_g)) DEALLOCATE(tsurf_daily_fm_g)
2192	IF (ALLOCATED(tsoil_daily_fm_g)) DEALLOCATE(tsoil_daily_fm_g)
2193	IF (ALLOCATED(soilhum_daily_fm_g)) DEALLOCATE(soilhum_daily_fm_g)
2194	IF (ALLOCATED(precip_fm_g)) DEALLOCATE(precip_fm_g)
2195	IF (ALLOCATED(gpp_daily_fm_g)) DEALLOCATE(gpp_daily_fm_g)
2196	IF (ALLOCATED(veget_fm_g)) DEALLOCATE(veget_fm_g)
2197	IF (ALLOCATED(veget_max_fm_g)) DEALLOCATE(veget_max_fm_g)
2198	IF (ALLOCATED(lai_fm_g)) DEALLOCATE(lai_fm_g)
2199	ENDIF
2200
2201	IF (ALLOCATED(isf)) DEALLOCATE(isf)
2202	IF (ALLOCATED(nf_written)) DEALLOCATE(nf_written)
2203	IF (ALLOCATED(nf_cumul)) DEALLOCATE(nf_cumul)
2204	IF (ALLOCATED(nforce)) DEALLOCATE(nforce)
2205	IF (ALLOCATED(control_moist)) DEALLOCATE(control_moist)
2206	IF (ALLOCATED(control_temp)) DEALLOCATE(control_temp)
2207	IF (ALLOCATED(soilcarbon_input)) DEALLOCATE(soilcarbon_input)
2208	IF ( ALLOCATED (horip10_index)) DEALLOCATE (horip10_index)
2209	IF ( ALLOCATED (horip100_index)) DEALLOCATE (horip100_index)
2210	IF ( ALLOCATED (horip11_index)) DEALLOCATE (horip11_index)
2211	IF ( ALLOCATED (horip101_index)) DEALLOCATE (horip101_index)
2212	IF ( ALLOCATED (prod10)) DEALLOCATE (prod10)
2213	IF ( ALLOCATED (prod100)) DEALLOCATE (prod100)
2214	IF ( ALLOCATED (flux10)) DEALLOCATE (flux10)
2215	IF ( ALLOCATED (flux100)) DEALLOCATE (flux100)
2216	IF ( ALLOCATED (convflux)) DEALLOCATE (convflux)
2217	IF ( ALLOCATED (cflux_prod10)) DEALLOCATE (cflux_prod10)
2218	IF ( ALLOCATED (cflux_prod100)) DEALLOCATE (cflux_prod100)
2219	IF ( ALLOCATED (harvest_above)) DEALLOCATE (harvest_above)
2220	IF ( ALLOCATED (soilcarbon_input_daily)) DEALLOCATE (soilcarbon_input_daily)
2221	IF ( ALLOCATED (control_temp_daily)) DEALLOCATE (control_temp_daily)
2222	IF ( ALLOCATED (control_moist_daily)) DEALLOCATE (control_moist_daily)
2223
2224	IF ( ALLOCATED (fpc_max)) DEALLOCATE (fpc_max)
2225
2226	! 2. reset l_first
2227	l_first_stomate=.TRUE.
2228	! 3. call to clear functions
2229	CALL get_reftemp_clear
2230	CALL season_clear
2231	CALL stomatelpj_clear
2232	CALL littercalc_clear
2233	CALL vmax_clear
2234	!---------------------------
2235	END SUBROUTINE stomate_clear
2236	!
2237	!=
2238	!
2239	SUBROUTINE stomate_var_init &
2240	& (kjpindex, veget_cov, veget_cov_max, leaf_age, leaf_frac, &
2241	& tlong_ref, t2m_month, dead_leaves, &
2242	& veget, lai, qsintmax, deadleaf_cover, assim_param)
2243
2244	!---------------------------------------------------------------------
2245	! this subroutine outputs values of assim_param etc.
2246	! only if ok_stomate = .TRUE.
2247	! otherwise,the calling procedure must do it itself.
2248	!
2249	! interface description
2250	! input scalar
2251	! Domain size
2252	INTEGER(i_std),INTENT(in) :: kjpindex
2253	! input fields
2254	! fractional coverage: actually covered space
2255	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: veget_cov
2256	! fractional coverage: maximum covered space
2257	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: veget_cov_max
2258	! "long term" 2 meter reference temperatures (K)
2259	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: tlong_ref
2260	! "monthly" 2 meter temperatures (K)
2261	REAL(r_std),DIMENSION(kjpindex),INTENT(in) :: t2m_month
2262	! dead leaves on ground, per PFT, metabolic and structural,
2263	! in gC/(m**2 of ground)
2264	REAL(r_std),DIMENSION(kjpindex,nvm,nlitt),INTENT(in) :: dead_leaves
2265	! Fraction of vegetation type
2266	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: veget
2267	! Surface foliere
2268	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(in) :: lai
2269	! modified fields (actually NOT modified)
2270	! leaf age (d)
2271	REAL(r_std),DIMENSION(kjpindex,nvm,nleafages),INTENT(inout) :: leaf_age
2272	! fraction of leaves in leaf age class
2273	REAL(r_std),DIMENSION(kjpindex,nvm,nleafages),INTENT(inout) :: leaf_frac
2274	! output scalar
2275	! output fields
2276	! Maximum water on vegetation for interception
2277	REAL(r_std),DIMENSION(kjpindex,nvm),INTENT(out) :: qsintmax
2278	! fraction of soil covered by dead leaves
2279	REAL(r_std),DIMENSION(kjpindex), INTENT (out) :: deadleaf_cover
2280	! min+max+opt temps & vmax for photosynthesis
2281	REAL(r_std),DIMENSION(kjpindex,nvm,npco2),INTENT(out) :: assim_param
2282
2283	!-
2284	! local declaration
2285	!-
2286	! dummy time step, must be zero
2287	REAL(r_std),PARAMETER :: dt_0 = zero
2288	REAL(r_std),DIMENSION(kjpindex,nvm) :: vcmax
2289	REAL(r_std),DIMENSION(kjpindex,nvm) :: vjmax
2290	! Min temperature for photosynthesis (deg C)
2291	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_photo_min
2292	! Opt temperature for photosynthesis (deg C)
2293	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_photo_opt
2294	! Max temperature for photosynthesis (deg C)
2295	REAL(r_std),DIMENSION(kjpindex,nvm) :: t_photo_max
2296	! Index
2297	INTEGER(i_std) :: j
2298	!---------------------------------------------------------------------
2299
2300	IF (control%ok_stomate) THEN
2301	!
2302	! 1 photosynthesis parameters
2303	!
2304	!
2305	! 1.1 vcmax
2306	! only if STOMATE is activated
2307	!
2308	CALL vmax (kjpindex, dt_0, leaf_age, leaf_frac, vcmax, vjmax)
2309	!
2310	! 1.2 assimilation temperatures
2311	!
2312	CALL assim_temp(kjpindex, tlong_ref, t2m_month, &
2313	& t_photo_min, t_photo_opt, t_photo_max)
2314	!
2315	! 1.3 transform into nvm vegetation types
2316	!
2317	assim_param(:,:,ivcmax) = zero
2318
2319	DO j = 2, nvm
2320	assim_param(:,j,ivcmax)=vcmax(:,j)
2321	ENDDO
2322
2323	assim_param(:,:,ivjmax) = zero
2324
2325	DO j = 2, nvm
2326	assim_param(:,j,ivjmax)=vjmax(:,j)
2327	ENDDO
2328
2329	assim_param(:,:,itmin) = zero
2330
2331	DO j = 2, nvm
2332	assim_param(:,j,itmin)=t_photo_min(:,j)
2333	ENDDO
2334
2335	assim_param(:,:,itopt) = zero
2336
2337	DO j = 2, nvm
2338	assim_param(:,j,itopt)=t_photo_opt(:,j)
2339	ENDDO
2340
2341	assim_param(:,:,itmax) = zero
2342
2343	DO j = 2, nvm
2344	assim_param(:,j,itmax)=t_photo_max(:,j)
2345	ENDDO
2346
2347	!
2348	! 2 dead leaf cover
2349	!
2350	! edit shilong
2351	! CALL deadleaf (kjpindex, dead_leaves, deadleaf_cover)
2352	CALL deadleaf (kjpindex, veget_cov_max, dead_leaves, deadleaf_cover)
2353	!
2354	! 3 qsintmax
2355	!
2356	qsintmax(:,ibare_sechiba) = zero
2357	DO j=2,nvm
2358	qsintmax(:,j) = qsintcstveget(:,j)lai(:,j)
2359	ENDDO
2360
2361	ENDIF ! ok_stomate = .TRUE.
2362	!--------------------------------
2363	END SUBROUTINE stomate_var_init
2364	!
2365	!=
2366	!
2367	SUBROUTINE stomate_accu &
2368	& (npts, n_dim2, dt_tot, dt, ldmean, field_in, field_out)
2369	!---------------------------------------------------------------------
2370	!
2371	! 0 declarations
2372	!
2373	! 0.1 input
2374	!
2375	! Domain size
2376	INTEGER(i_std),INTENT(in) :: npts
2377	! 2nd dimension (1 or nvm)
2378	INTEGER(i_std),INTENT(in) :: n_dim2
2379	! Time step of STOMATE (days)
2380	REAL(r_std),INTENT(in) :: dt_tot
2381	! Time step in days
2382	REAL(r_std),INTENT(in) :: dt
2383	! Calculate mean ?
2384	LOGICAL,INTENT(in) :: ldmean
2385	! Daily field
2386	REAL(r_std),DIMENSION(npts,n_dim2),INTENT(in) :: field_in
2387	!
2388	! 0.2 modified field
2389	!
2390	! Annual field
2391	REAL(r_std),DIMENSION(npts,n_dim2),INTENT(inout) :: field_out
2392	!---------------------------------------------------------------------
2393	!
2394	! 1 accumulation
2395	!
2396	field_out(:,:) = field_out(:,:)+field_in(:,:)*dt
2397	!
2398	! 2 mean fields
2399	!
2400	IF (ldmean) THEN
2401	field_out(:,:) = field_out(:,:)/dt_tot
2402	ENDIF
2403	!---------------------------------------------------------------------
2404	END SUBROUTINE stomate_accu
2405
2406	SUBROUTINE init_forcing (kjpindex,nsfm,nsft)
2407	!---------------------------------------------------------------------
2408	INTEGER(i_std),INTENT(in) :: kjpindex
2409	INTEGER(i_std),INTENT(in) :: nsfm
2410	INTEGER(i_std),INTENT(in) :: nsft
2411	!
2412	LOGICAL :: l_error
2413	INTEGER(i_std) :: ier
2414	!---------------------------------------------------------------------
2415	l_error = .FALSE.
2416	!
2417	ALLOCATE(clay_fm(kjpindex,nsfm),stat=ier)
2418	l_error = l_error .OR. (ier /= 0)
2419	IF (l_error) THEN
2420	WRITE(numout,*) 'Problem with memory allocation: forcing variables clay_fm ',kjpindex,nsfm
2421	STOP 'init_forcing'
2422	ENDIF
2423	ALLOCATE(humrel_daily_fm(kjpindex,nvm,nsfm),stat=ier)
2424	l_error = l_error .OR. (ier /= 0)
2425	IF (l_error) THEN
2426	WRITE(numout,*) 'Problem with memory allocation: forcing variables humrel_daily_fm ',kjpindex,nvm,nsfm
2427	STOP 'init_forcing'
2428	ENDIF
2429	ALLOCATE(litterhum_daily_fm(kjpindex,nsfm),stat=ier)
2430	l_error = l_error .OR. (ier /= 0)
2431	IF (l_error) THEN
2432	WRITE(numout,*) 'Problem with memory allocation: forcing variables litterhum_daily_fm ',kjpindex,nsfm
2433	STOP 'init_forcing'
2434	ENDIF
2435	ALLOCATE(t2m_daily_fm(kjpindex,nsfm),stat=ier)
2436	l_error = l_error .OR. (ier /= 0)
2437	IF (l_error) THEN
2438	WRITE(numout,*) 'Problem with memory allocation: forcing variables t2m_daily_fm ',kjpindex,nsfm
2439	STOP 'init_forcing'
2440	ENDIF
2441	ALLOCATE(t2m_min_daily_fm(kjpindex,nsfm),stat=ier)
2442	l_error = l_error .OR. (ier /= 0)
2443	IF (l_error) THEN
2444	WRITE(numout,*) 'Problem with memory allocation: forcing variables t2m_min_daily_fm ',kjpindex,nsfm
2445	STOP 'init_forcing'
2446	ENDIF
2447	ALLOCATE(tsurf_daily_fm(kjpindex,nsfm),stat=ier)
2448	l_error = l_error .OR. (ier /= 0)
2449	IF (l_error) THEN
2450	WRITE(numout,*) 'Problem with memory allocation: forcing variables tsurf_daily_fm ',kjpindex,nsfm
2451	STOP 'init_forcing'
2452	ENDIF
2453	ALLOCATE(tsoil_daily_fm(kjpindex,nbdl,nsfm),stat=ier)
2454	l_error = l_error .OR. (ier /= 0)
2455	IF (l_error) THEN
2456	WRITE(numout,*) 'Problem with memory allocation: forcing variables tsoil_daily_fm ',kjpindex,nbdl,nsfm
2457	STOP 'init_forcing'
2458	ENDIF
2459	ALLOCATE(soilhum_daily_fm(kjpindex,nbdl,nsfm),stat=ier)
2460	l_error = l_error .OR. (ier /= 0)
2461	IF (l_error) THEN
2462	WRITE(numout,*) 'Problem with memory allocation: forcing variables soilhum_daily_fm ',kjpindex,nbdl,nsfm
2463	STOP 'init_forcing'
2464	ENDIF
2465	ALLOCATE(precip_fm(kjpindex,nsfm),stat=ier)
2466	l_error = l_error .OR. (ier /= 0)
2467	IF (l_error) THEN
2468	WRITE(numout,*) 'Problem with memory allocation: forcing variables precip_fm ',kjpindex,nsfm
2469	STOP 'init_forcing'
2470	ENDIF
2471	ALLOCATE(gpp_daily_fm(kjpindex,nvm,nsfm),stat=ier)
2472	l_error = l_error .OR. (ier /= 0)
2473	IF (l_error) THEN
2474	WRITE(numout,*) 'Problem with memory allocation: forcing variables gpp_daily_fm ',kjpindex,nvm,nsfm
2475	STOP 'init_forcing'
2476	ENDIF
2477	ALLOCATE(veget_fm(kjpindex,nvm,nsfm),stat=ier)
2478	l_error = l_error .OR. (ier /= 0)
2479	IF (l_error) THEN
2480	WRITE(numout,*) 'Problem with memory allocation: forcing variables veget_fm ',kjpindex,nvm,nsfm
2481	STOP 'init_forcing'
2482	ENDIF
2483	ALLOCATE(veget_max_fm(kjpindex,nvm,nsfm),stat=ier)
2484	l_error = l_error .OR. (ier /= 0)
2485	IF (l_error) THEN
2486	WRITE(numout,*) 'Problem with memory allocation: forcing variables veget_max_fm ',kjpindex,nvm,nsfm
2487	STOP 'init_forcing'
2488	ENDIF
2489	ALLOCATE(lai_fm(kjpindex,nvm,nsfm),stat=ier)
2490	l_error = l_error .OR. (ier /= 0)
2491	IF (l_error) THEN
2492	WRITE(numout,*) 'Problem with memory allocation: forcing variables lai_fm ',kjpindex,nvm,nsfm
2493	STOP 'init_forcing'
2494	ENDIF
2495	ALLOCATE(isf(nsfm),stat=ier)
2496	l_error = l_error .OR. (ier /= 0)
2497	IF (l_error) THEN
2498	WRITE(numout,*) 'Problem with memory allocation: forcing variables isf ',nsfm
2499	STOP 'init_forcing'
2500	ENDIF
2501	ALLOCATE(nf_written(nsft),stat=ier)
2502	l_error = l_error .OR. (ier /= 0)
2503	IF (l_error) THEN
2504	WRITE(numout,*) 'Problem with memory allocation: forcing variables nf_written ',nsft
2505	STOP 'init_forcing'
2506	ENDIF
2507	ALLOCATE(nf_cumul(nsft),stat=ier)
2508	l_error = l_error .OR. (ier /= 0)
2509	IF (l_error) THEN
2510	WRITE(numout,*) 'Problem with memory allocation: forcing variables nf_cumul ',nsft
2511	STOP 'init_forcing'
2512	ENDIF
2513
2514	IF (is_root_prc) THEN
2515	ALLOCATE(clay_fm_g(nbp_glo,nsfm),stat=ier)
2516	l_error = l_error .OR. (ier /= 0)
2517	IF (l_error) THEN
2518	WRITE(numout,*) 'Problem with memory allocation: forcing variables clay_fm_g ',nbp_glo,nsfm
2519	STOP 'init_forcing'
2520	ENDIF
2521	ALLOCATE(humrel_daily_fm_g(nbp_glo,nvm,nsfm),stat=ier)
2522	l_error = l_error .OR. (ier /= 0)
2523	IF (l_error) THEN
2524	WRITE(numout,*) 'Problem with memory allocation: forcing variables humrel_daily_fm_g ',nbp_glo,nvm,nsfm
2525	STOP 'init_forcing'
2526	ENDIF
2527	ALLOCATE(litterhum_daily_fm_g(nbp_glo,nsfm),stat=ier)
2528	l_error = l_error .OR. (ier /= 0)
2529	IF (l_error) THEN
2530	WRITE(numout,*) 'Problem with memory allocation: forcing variables litterhum_daily_fm_g ',nbp_glo,nsfm
2531	STOP 'init_forcing'
2532	ENDIF
2533	ALLOCATE(t2m_daily_fm_g(nbp_glo,nsfm),stat=ier)
2534	l_error = l_error .OR. (ier /= 0)
2535	IF (l_error) THEN
2536	WRITE(numout,*) 'Problem with memory allocation: forcing variables t2m_daily_fm_g ',nbp_glo,nsfm
2537	STOP 'init_forcing'
2538	ENDIF
2539	ALLOCATE(t2m_min_daily_fm_g(nbp_glo,nsfm),stat=ier)
2540	l_error = l_error .OR. (ier /= 0)
2541	IF (l_error) THEN
2542	WRITE(numout,*) 'Problem with memory allocation: forcing variables t2m_min_daily_fm_g ',nbp_glo,nsfm
2543	STOP 'init_forcing'
2544	ENDIF
2545	ALLOCATE(tsurf_daily_fm_g(nbp_glo,nsfm),stat=ier)
2546	l_error = l_error .OR. (ier /= 0)
2547	IF (l_error) THEN
2548	WRITE(numout,*) 'Problem with memory allocation: forcing variables tsurf_daily_fm_g ',nbp_glo,nsfm
2549	STOP 'init_forcing'
2550	ENDIF
2551	ALLOCATE(tsoil_daily_fm_g(nbp_glo,nbdl,nsfm),stat=ier)
2552	l_error = l_error .OR. (ier /= 0)
2553	IF (l_error) THEN
2554	WRITE(numout,*) 'Problem with memory allocation: forcing variables tsoil_daily_fm_g ',nbp_glo,nbdl,nsfm
2555	STOP 'init_forcing'
2556	ENDIF
2557	ALLOCATE(soilhum_daily_fm_g(nbp_glo,nbdl,nsfm),stat=ier)
2558	l_error = l_error .OR. (ier /= 0)
2559	IF (l_error) THEN
2560	WRITE(numout,*) 'Problem with memory allocation: forcing variables soilhum_daily_fm_g ',nbp_glo,nbdl,nsfm
2561	STOP 'init_forcing'
2562	ENDIF
2563	ALLOCATE(precip_fm_g(nbp_glo,nsfm),stat=ier)
2564	l_error = l_error .OR. (ier /= 0)
2565	IF (l_error) THEN
2566	WRITE(numout,*) 'Problem with memory allocation: forcing variables precip_fm_g ',nbp_glo,nsfm
2567	STOP 'init_forcing'
2568	ENDIF
2569	ALLOCATE(gpp_daily_fm_g(nbp_glo,nvm,nsfm),stat=ier)
2570	l_error = l_error .OR. (ier /= 0)
2571	IF (l_error) THEN
2572	WRITE(numout,*) 'Problem with memory allocation: forcing variables gpp_daily_fm_g ',nbp_glo,nvm,nsfm
2573	STOP 'init_forcing'
2574	ENDIF
2575	ALLOCATE(veget_fm_g(nbp_glo,nvm,nsfm),stat=ier)
2576	l_error = l_error .OR. (ier /= 0)
2577	IF (l_error) THEN
2578	WRITE(numout,*) 'Problem with memory allocation: forcing variables veget_fm_g ',nbp_glo,nvm,nsfm
2579	STOP 'init_forcing'
2580	ENDIF
2581	ALLOCATE(veget_max_fm_g(nbp_glo,nvm,nsfm),stat=ier)
2582	l_error = l_error .OR. (ier /= 0)
2583	IF (l_error) THEN
2584	WRITE(numout,*) 'Problem with memory allocation: forcing variables veget_max_fm_g ',nbp_glo,nvm,nsfm
2585	STOP 'init_forcing'
2586	ENDIF
2587	ALLOCATE(lai_fm_g(nbp_glo,nvm,nsfm),stat=ier)
2588	l_error = l_error .OR. (ier /= 0)
2589	IF (l_error) THEN
2590	WRITE(numout,*) 'Problem with memory allocation: forcing variables lai_fm_g ',nbp_glo,nvm,nsfm
2591	STOP 'init_forcing'
2592	ENDIF
2593	ELSE
2594	ALLOCATE(clay_fm_g(0,nsfm),stat=ier)
2595	ALLOCATE(humrel_daily_fm_g(0,nvm,nsfm),stat=ier)
2596	ALLOCATE(litterhum_daily_fm_g(0,nsfm),stat=ier)
2597	ALLOCATE(t2m_daily_fm_g(0,nsfm),stat=ier)
2598	ALLOCATE(t2m_min_daily_fm_g(0,nsfm),stat=ier)
2599	ALLOCATE(tsurf_daily_fm_g(0,nsfm),stat=ier)
2600	ALLOCATE(tsoil_daily_fm_g(0,nbdl,nsfm),stat=ier)
2601	ALLOCATE(soilhum_daily_fm_g(0,nbdl,nsfm),stat=ier)
2602	ALLOCATE(precip_fm_g(0,nsfm),stat=ier)
2603	ALLOCATE(gpp_daily_fm_g(0,nvm,nsfm),stat=ier)
2604	ALLOCATE(veget_fm_g(0,nvm,nsfm),stat=ier)
2605	ALLOCATE(veget_max_fm_g(0,nvm,nsfm),stat=ier)
2606	ALLOCATE(lai_fm_g(0,nvm,nsfm),stat=ier)
2607	ENDIF
2608	!
2609	IF (l_error) THEN
2610	WRITE(numout,*) 'Problem with memory allocation: forcing variables'
2611	STOP 'init_forcing'
2612	ENDIF
2613	!
2614	CALL forcing_zero
2615	!--------------------------
2616	END SUBROUTINE init_forcing
2617	!
2618	!=
2619	!
2620	SUBROUTINE forcing_zero
2621	!---------------------------------------------------------------------
2622	clay_fm(:,:) = zero
2623	humrel_daily_fm(:,:,:) = zero
2624	litterhum_daily_fm(:,:) = zero
2625	t2m_daily_fm(:,:) = zero
2626	t2m_min_daily_fm(:,:) = zero
2627	tsurf_daily_fm(:,:) = zero
2628	tsoil_daily_fm(:,:,:) = zero
2629	soilhum_daily_fm(:,:,:) = zero
2630	precip_fm(:,:) = zero
2631	gpp_daily_fm(:,:,:) = zero
2632	veget_fm(:,:,:) = zero
2633	veget_max_fm(:,:,:) = zero
2634	lai_fm(:,:,:) = zero
2635	!--------------------------------
2636	END SUBROUTINE forcing_zero
2637	!
2638	!=
2639	!
2640	SUBROUTINE forcing_write(forcing_id,ibeg,iend)
2641	!---------------------------------------------------------------------
2642	INTEGER(i_std),INTENT(in) :: forcing_id
2643	INTEGER(i_std),INTENT(in) :: ibeg, iend
2644	!
2645	INTEGER(i_std) :: iisf, iblocks, nblocks
2646	INTEGER(i_std) :: ier
2647	INTEGER(i_std),DIMENSION(0:2) :: ifirst, ilast
2648	INTEGER(i_std),PARAMETER :: ndm = 10
2649	INTEGER(i_std),DIMENSION(ndm) :: start, count_force
2650	INTEGER(i_std) :: ndim, vid
2651	!---------------------------------------------------------------------
2652	!
2653	! determine blocks of forcing states that are contiguous in memory
2654	!
2655	nblocks = 0
2656	ifirst(:) = 1
2657	ilast(:) = 1
2658	!
2659	DO iisf = ibeg, iend
2660	IF ( (nblocks /= 0) &
2661	& .AND.(isf(iisf) == isf(ilast(nblocks))+1)) THEN
2662	! element is contiguous with last element found
2663	ilast(nblocks) = iisf
2664	ELSE
2665	! found first element of new block
2666	nblocks = nblocks+1
2667	IF (nblocks > 2) STOP 'Problem in forcing_write'
2668	ifirst(nblocks) = iisf
2669	ilast(nblocks) = iisf
2670	ENDIF
2671	ENDDO
2672	!
2673	CALL gather(clay_fm,clay_fm_g)
2674	CALL gather(humrel_daily_fm,humrel_daily_fm_g)
2675	CALL gather(litterhum_daily_fm,litterhum_daily_fm_g)
2676	CALL gather(t2m_daily_fm,t2m_daily_fm_g)
2677	CALL gather(t2m_min_daily_fm,t2m_min_daily_fm_g)
2678	CALL gather(tsurf_daily_fm,tsurf_daily_fm_g)
2679	CALL gather(tsoil_daily_fm,tsoil_daily_fm_g)
2680	CALL gather(soilhum_daily_fm,soilhum_daily_fm_g)
2681	CALL gather(precip_fm,precip_fm_g)
2682	CALL gather(gpp_daily_fm,gpp_daily_fm_g)
2683	CALL gather(veget_fm,veget_fm_g)
2684	CALL gather(veget_max_fm,veget_max_fm_g)
2685	CALL gather(lai_fm,lai_fm_g)
2686	IF (is_root_prc) THEN
2687	DO iblocks = 1, nblocks
2688	IF (ifirst(iblocks) /= ilast(iblocks)) THEN
2689	ndim = 2
2690	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2691	count_force(1:ndim) = SHAPE(clay_fm_g)
2692	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2693	ier = NF90_INQ_VARID (forcing_id,'clay',vid)
2694	ier = NF90_PUT_VAR (forcing_id,vid, &
2695	& clay_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2696	& start=start(1:ndim), count=count_force(1:ndim))
2697	ndim = 3;
2698	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2699	count_force(1:ndim) = SHAPE(humrel_daily_fm_g)
2700	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2701	ier = NF90_INQ_VARID (forcing_id,'humrel',vid)
2702	ier = NF90_PUT_VAR (forcing_id, vid, &
2703	& humrel_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2704	& start=start(1:ndim), count=count_force(1:ndim))
2705	ndim = 2;
2706	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2707	count_force(1:ndim) = SHAPE(litterhum_daily_fm_g)
2708	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2709	ier = NF90_INQ_VARID (forcing_id,'litterhum',vid)
2710	ier = NF90_PUT_VAR (forcing_id, vid, &
2711	& litterhum_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2712	& start=start(1:ndim), count=count_force(1:ndim))
2713	ndim = 2;
2714	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2715	count_force(1:ndim) = SHAPE(t2m_daily_fm_g)
2716	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2717	ier = NF90_INQ_VARID (forcing_id,'t2m',vid)
2718	ier = NF90_PUT_VAR (forcing_id, vid, &
2719	& t2m_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2720	& start=start(1:ndim), count=count_force(1:ndim))
2721	ndim = 2;
2722	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2723	count_force(1:ndim) = SHAPE(t2m_min_daily_fm_g)
2724	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2725	ier = NF90_INQ_VARID (forcing_id,'t2m_min',vid)
2726	ier = NF90_PUT_VAR (forcing_id, vid, &
2727	& t2m_min_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2728	& start=start(1:ndim), count=count_force(1:ndim))
2729	ndim = 2;
2730	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2731	count_force(1:ndim) = SHAPE(tsurf_daily_fm_g)
2732	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2733	ier = NF90_INQ_VARID (forcing_id,'tsurf',vid)
2734	ier = NF90_PUT_VAR (forcing_id, vid, &
2735	& tsurf_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2736	& start=start(1:ndim), count=count_force(1:ndim))
2737	ndim = 3;
2738	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2739	count_force(1:ndim) = SHAPE(tsoil_daily_fm_g)
2740	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2741	ier = NF90_INQ_VARID (forcing_id,'tsoil',vid)
2742	ier = NF90_PUT_VAR (forcing_id, vid, &
2743	& tsoil_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2744	& start=start(1:ndim), count=count_force(1:ndim))
2745	ndim = 3;
2746	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2747	count_force(1:ndim) = SHAPE(soilhum_daily_fm_g)
2748	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2749	ier = NF90_INQ_VARID (forcing_id,'soilhum',vid)
2750	ier = NF90_PUT_VAR (forcing_id, vid, &
2751	& soilhum_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2752	& start=start(1:ndim), count=count_force(1:ndim))
2753	ndim = 2;
2754	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2755	count_force(1:ndim) = SHAPE(precip_fm_g)
2756	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2757	ier = NF90_INQ_VARID (forcing_id,'precip',vid)
2758	ier = NF90_PUT_VAR (forcing_id, vid, &
2759	& precip_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2760	& start=start(1:ndim), count=count_force(1:ndim))
2761	ndim = 3;
2762	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2763	count_force(1:ndim) = SHAPE(gpp_daily_fm_g)
2764	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2765	ier = NF90_INQ_VARID (forcing_id,'gpp',vid)
2766	ier = NF90_PUT_VAR (forcing_id, vid, &
2767	& gpp_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2768	& start=start(1:ndim), count=count_force(1:ndim))
2769	ndim = 3;
2770	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2771	count_force(1:ndim) = SHAPE(veget_fm_g)
2772	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2773	ier = NF90_INQ_VARID (forcing_id,'veget',vid)
2774	ier = NF90_PUT_VAR (forcing_id, vid, &
2775	& veget_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2776	& start=start(1:ndim), count=count_force(1:ndim))
2777	ndim = 3;
2778	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2779	count_force(1:ndim) = SHAPE(veget_max_fm_g)
2780	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2781	ier = NF90_INQ_VARID (forcing_id,'veget_max',vid)
2782	ier = NF90_PUT_VAR (forcing_id, vid, &
2783	& veget_max_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2784	& start=start(1:ndim), count=count_force(1:ndim))
2785	ndim = 3;
2786	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2787	count_force(1:ndim) = SHAPE(lai_fm_g)
2788	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2789	ier = NF90_INQ_VARID (forcing_id,'lai',vid)
2790	ier = NF90_PUT_VAR (forcing_id, vid, &
2791	& lai_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2792	& start=start(1:ndim), count=count_force(1:ndim))
2793	ENDIF
2794	ENDDO
2795	ENDIF
2796	nf_written(isf(:)) = .TRUE.
2797	!---------------------------
2798	END SUBROUTINE forcing_write
2799	!
2800	!=
2801	!
2802	SUBROUTINE forcing_read(forcing_id,nsfm)
2803	!---------------------------------------------------------------------
2804	INTEGER(i_std),INTENT(in) :: forcing_id
2805	INTEGER(i_std),INTENT(in) :: nsfm
2806	!
2807	INTEGER(i_std) :: iisf, iblocks, nblocks
2808	INTEGER(i_std) :: ier
2809	LOGICAL :: a_er
2810	INTEGER(i_std),DIMENSION(0:2) :: ifirst, ilast
2811	INTEGER(i_std),PARAMETER :: ndm = 10
2812	INTEGER(i_std),DIMENSION(ndm) :: start, count_force
2813	INTEGER(i_std) :: ndim, vid
2814
2815	LOGICAL, PARAMETER :: check=.FALSE.
2816
2817	IF (check) WRITE(numout,*) "forcing_read "
2818	!---------------------------------------------------------------------
2819	!
2820	! set to zero if the corresponding forcing state
2821	! has not yet been written into the file
2822	!
2823	DO iisf = 1, nsfm
2824	IF (.NOT.nf_written(isf(iisf))) THEN
2825	clay_fm(:,iisf) = zero
2826	humrel_daily_fm(:,:,iisf) = zero
2827	litterhum_daily_fm(:,iisf) = zero
2828	t2m_daily_fm(:,iisf) = zero
2829	t2m_min_daily_fm(:,iisf) = zero
2830	tsurf_daily_fm(:,iisf) = zero
2831	tsoil_daily_fm(:,:,iisf) = zero
2832	soilhum_daily_fm(:,:,iisf) = zero
2833	precip_fm(:,iisf) = zero
2834	gpp_daily_fm(:,:,iisf) = zero
2835	veget_fm(:,:,iisf) = zero
2836	veget_max_fm(:,:,iisf) = zero
2837	lai_fm(:,:,iisf) = zero
2838	ENDIF
2839	ENDDO
2840	!
2841	! determine blocks of forcing states that are contiguous in memory
2842	!
2843	nblocks = 0
2844	ifirst(:) = 1
2845	ilast(:) = 1
2846	!
2847	DO iisf = 1, nsfm
2848	IF (nf_written(isf(iisf))) THEN
2849	IF ( (nblocks /= 0) &
2850	& .AND.(isf(iisf) == isf(ilast(nblocks))+1)) THEN
2851	! element is contiguous with last element found
2852	ilast(nblocks) = iisf
2853	ELSE
2854	! found first element of new block
2855	nblocks = nblocks+1
2856	IF (nblocks > 2) STOP 'Problem in forcing_read'
2857	!
2858	ifirst(nblocks) = iisf
2859	ilast(nblocks) = iisf
2860	ENDIF
2861	ENDIF
2862	ENDDO
2863	IF (check) WRITE(numout,*) "forcing_read nblocks, ifirst, ilast",nblocks, ifirst, ilast
2864	!
2865	IF (is_root_prc) THEN
2866	DO iblocks = 1, nblocks
2867	IF (check) WRITE(numout,*) "forcing_read iblocks, ifirst(iblocks), ilast(iblocks)",iblocks, &
2868	ifirst(iblocks), ilast(iblocks)
2869	IF (ifirst(iblocks) /= ilast(iblocks)) THEN
2870	a_er=.FALSE.
2871	ndim = 2;
2872	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2873	count_force(1:ndim) = SHAPE(clay_fm_g)
2874	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2875	ier = NF90_INQ_VARID (forcing_id,'clay',vid)
2876	a_er = a_er.OR.(ier.NE.0)
2877	ier = NF90_GET_VAR (forcing_id, vid, &
2878	& clay_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2879	& start=start(1:ndim), count=count_force(1:ndim))
2880	a_er = a_er.OR.(ier.NE.0)
2881	!---------
2882	ndim = 3;
2883	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2884	count_force(1:ndim) = SHAPE(humrel_daily_fm_g)
2885	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2886	ier = NF90_INQ_VARID (forcing_id,'humrel',vid)
2887	a_er = a_er.OR.(ier.NE.0)
2888	ier = NF90_GET_VAR (forcing_id, vid, &
2889	& humrel_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2890	& start=start(1:ndim), count=count_force(1:ndim))
2891	a_er = a_er.OR.(ier.NE.0)
2892	!---------
2893	ndim = 2;
2894	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2895	count_force(1:ndim) = SHAPE(litterhum_daily_fm_g)
2896	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2897	ier = NF90_INQ_VARID (forcing_id,'litterhum',vid)
2898	a_er = a_er.OR.(ier.NE.0)
2899	ier = NF90_GET_VAR (forcing_id, vid, &
2900	& litterhum_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2901	& start=start(1:ndim), count=count_force(1:ndim))
2902	a_er = a_er.OR.(ier.NE.0)
2903	!---------
2904	ndim = 2;
2905	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2906	count_force(1:ndim) = SHAPE(t2m_daily_fm_g)
2907	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2908	ier = NF90_INQ_VARID (forcing_id,'t2m',vid)
2909	a_er = a_er.OR.(ier.NE.0)
2910	ier = NF90_GET_VAR (forcing_id, vid, &
2911	& t2m_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2912	& start=start(1:ndim), count=count_force(1:ndim))
2913	a_er = a_er.OR.(ier.NE.0)
2914	!---------
2915	ndim = 2;
2916	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2917	count_force(1:ndim) = SHAPE(t2m_min_daily_fm_g)
2918	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2919	ier = NF90_INQ_VARID (forcing_id,'t2m_min',vid)
2920	a_er = a_er.OR.(ier.NE.0)
2921	ier = NF90_GET_VAR (forcing_id, vid, &
2922	& t2m_min_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2923	& start=start(1:ndim), count=count_force(1:ndim))
2924	a_er = a_er.OR.(ier.NE.0)
2925	!---------
2926	ndim = 2;
2927	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2928	count_force(1:ndim) = SHAPE(tsurf_daily_fm_g)
2929	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2930	ier = NF90_INQ_VARID (forcing_id,'tsurf',vid)
2931	a_er = a_er.OR.(ier.NE.0)
2932	ier = NF90_GET_VAR (forcing_id, vid, &
2933	& tsurf_daily_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2934	& start=start(1:ndim), count=count_force(1:ndim))
2935	a_er = a_er.OR.(ier.NE.0)
2936	!---------
2937	ndim = 3;
2938	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2939	count_force(1:ndim) = SHAPE(tsoil_daily_fm_g)
2940	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2941	ier = NF90_INQ_VARID (forcing_id,'tsoil',vid)
2942	a_er = a_er.OR.(ier.NE.0)
2943	ier = NF90_GET_VAR (forcing_id, vid, &
2944	& tsoil_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2945	& start=start(1:ndim), count=count_force(1:ndim))
2946	a_er = a_er.OR.(ier.NE.0)
2947	!---------
2948	ndim = 3;
2949	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2950	count_force(1:ndim) = SHAPE(soilhum_daily_fm_g)
2951	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2952	ier = NF90_INQ_VARID (forcing_id,'soilhum',vid)
2953	a_er = a_er.OR.(ier.NE.0)
2954	ier = NF90_GET_VAR (forcing_id, vid, &
2955	& soilhum_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2956	& start=start(1:ndim), count=count_force(1:ndim))
2957	a_er = a_er.OR.(ier.NE.0)
2958	!---------
2959	ndim = 2;
2960	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2961	count_force(1:ndim) = SHAPE(precip_fm_g)
2962	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2963	ier = NF90_INQ_VARID (forcing_id,'precip',vid)
2964	a_er = a_er.OR.(ier.NE.0)
2965	ier = NF90_GET_VAR (forcing_id, vid, &
2966	& precip_fm_g(:,ifirst(iblocks):ilast(iblocks)), &
2967	& start=start(1:ndim), count=count_force(1:ndim))
2968	a_er = a_er.OR.(ier.NE.0)
2969	!---------
2970	ndim = 3;
2971	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2972	count_force(1:ndim) = SHAPE(gpp_daily_fm_g)
2973	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2974	ier = NF90_INQ_VARID (forcing_id,'gpp',vid)
2975	a_er = a_er.OR.(ier.NE.0)
2976	ier = NF90_GET_VAR (forcing_id, vid, &
2977	& gpp_daily_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2978	& start=start(1:ndim), count=count_force(1:ndim))
2979	a_er = a_er.OR.(ier.NE.0)
2980	!---------
2981	ndim = 3;
2982	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2983	count_force(1:ndim) = SHAPE(veget_fm_g)
2984	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2985	ier = NF90_INQ_VARID (forcing_id,'veget',vid)
2986	a_er = a_er.OR.(ier.NE.0)
2987	ier = NF90_GET_VAR (forcing_id, vid, &
2988	& veget_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
2989	& start=start(1:ndim), count=count_force(1:ndim))
2990	a_er = a_er.OR.(ier.NE.0)
2991	!---------
2992	ndim = 3;
2993	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
2994	count_force(1:ndim) = SHAPE(veget_max_fm_g)
2995	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
2996	ier = NF90_INQ_VARID (forcing_id,'veget_max',vid)
2997	a_er = a_er.OR.(ier.NE.0)
2998	ier = NF90_GET_VAR (forcing_id, vid, &
2999	& veget_max_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
3000	& start=start(1:ndim), count=count_force(1:ndim))
3001	a_er = a_er.OR.(ier.NE.0)
3002	!---------
3003	ndim = 3;
3004	start(1:ndim) = 1; start(ndim) = isf(ifirst(iblocks));
3005	count_force(1:ndim) = SHAPE(lai_fm_g)
3006	count_force(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
3007	ier = NF90_INQ_VARID (forcing_id,'lai',vid)
3008	a_er = a_er.OR.(ier.NE.0)
3009	ier = NF90_GET_VAR (forcing_id, vid, &
3010	& lai_fm_g(:,:,ifirst(iblocks):ilast(iblocks)), &
3011	& start=start(1:ndim), count=count_force(1:ndim))
3012	a_er = a_er.OR.(ier.NE.0)
3013	IF (a_er) THEN
3014	CALL ipslerr (3,'forcing_read', &
3015	& 'PROBLEM when read forcing file', &
3016	& '','')
3017	ENDIF
3018	ENDIF
3019	ENDDO
3020	ENDIF
3021	CALL scatter(clay_fm_g,clay_fm)
3022	CALL scatter(humrel_daily_fm_g,humrel_daily_fm)
3023	CALL scatter(litterhum_daily_fm_g,litterhum_daily_fm)
3024	CALL scatter(t2m_daily_fm_g,t2m_daily_fm)
3025	CALL scatter(t2m_min_daily_fm_g,t2m_min_daily_fm)
3026	CALL scatter(tsurf_daily_fm_g,tsurf_daily_fm)
3027	CALL scatter(tsoil_daily_fm_g,tsoil_daily_fm)
3028	CALL scatter(soilhum_daily_fm_g,soilhum_daily_fm)
3029	CALL scatter(precip_fm_g,precip_fm)
3030	CALL scatter(gpp_daily_fm_g,gpp_daily_fm)
3031	CALL scatter(veget_fm_g,veget_fm)
3032	CALL scatter(veget_max_fm_g,veget_max_fm)
3033	CALL scatter(lai_fm_g,lai_fm)
3034	!--------------------------
3035	END SUBROUTINE forcing_read
3036	!
3037	!=
3038	!
3039	SUBROUTINE setlai(npts,lai)
3040	!---------------------------------------------------------------------
3041	! routine to force the lai in STOMATE (for assimilation procedures)
3042	! for the moment setlai only gives the lai from stomate,
3043	! this routine must be written in the future
3044	!
3045	! 0 declarations
3046	!
3047	! 0.1 input
3048	!
3049	! Domain size
3050	INTEGER(i_std),INTENT(in) :: npts
3051	! 0.3 output
3052	REAL(r_std),DIMENSION(npts,nvm),INTENT(out) :: lai
3053	! 0.4 local definitions
3054	INTEGER(i_std) :: j
3055	!---------------------------------------------------------------------
3056	lai(:,ibare_sechiba) = zero
3057	DO j=2,nvm
3058	lai(:,j) = biomass(:,j,ileaf)*sla(j)
3059	ENDDO
3060	!--------------------
3061	END SUBROUTINE setlai
3062	!
3063	!
3064	!-----------------
3065	END MODULE stomate

Note: See TracBrowser for help on using the repository browser.

Download in other formats: