Context Navigation

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

Last change on this file 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: