Context Navigation

source: branches/publications/ORCHIDEE_CAN_r2290/src_stomate/stomate_phenology.f90 @ 7474

Last change on this file since 7474 was 2275, checked in by sebastiaan.luyssaert, 10 years ago
DEV: works on a single pixel. Experimenting with adaptation to waterstress. Committed to migrate the code to Curie. Do not use for regional simulations. Updates will follow.
Property svn:keywords set to `HeadURL Date Author Revision`
File size: 132.4 KB

Line
1	! =================================================================================================================================
2	! MODULE : stomate_phenology
3	!
4	! CONTACT : orchidee-help _at_ ipsl.jussieu.fr
5	!
6	! LICENCE : IPSL (2006). This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
7	!
8	!>\BRIEF This module manages the beginning of the growing season (leaf onset).
9	!!
10	!!\n DESCRIPTION: None
11	!!
12	!! RECENT CHANGE(S): None
13	!!
14	!! SVN :
15	!! $HeadURL$
16	!! $Date$
17	!! $Revision$
18	!! \n
19	!_ =================================================================================================================================
20
21	MODULE stomate_phenology
22
23	! modules used:
24
25	USE ioipsl_para
26	USE xios_orchidee
27	USE ioipsl_para
28	USE stomate_data
29	USE constantes
30	USE pft_parameters
31	USE sapiens_agriculture, ONLY: crop_planting
32	USE function_library, ONLY: wood_to_height_eff, calculate_c0_alloc
33
34	IMPLICIT NONE
35
36	! private & public routines
37
38	PRIVATE
39	PUBLIC phenology_diagnostic, phenology_prognostic, phenology_clear
40
41	! first call
42	LOGICAL, SAVE :: firstcall = .TRUE.
43	!$OMP THREADPRIVATE(firstcall)
44	LOGICAL, SAVE :: firstcall_hum = .TRUE.
45	!$OMP THREADPRIVATE(firstcall_hum)
46	LOGICAL, SAVE :: firstcall_moi = .TRUE.
47	!$OMP THREADPRIVATE(firstcall_moi)
48	LOGICAL, SAVE :: firstcall_humgdd = .TRUE.
49	!$OMP THREADPRIVATE(firstcall_humgdd)
50	LOGICAL, SAVE :: firstcall_moigdd = .TRUE.
51	!$OMP THREADPRIVATE(firstcall_moigdd)
52
53	CONTAINS
54
55
56	!! ================================================================================================================================
57	!! SUBROUTINE : phenology_clear
58	!!
59	!>\BRIEF Flags setting
60	!!
61	!! DESCRIPTION : This subroutine sets flags
62	!! ::firstcall, ::firstcall_hum, ::firstcall_moi, ::firstcall_humgdd,
63	!! ::firstcall_moigdd to .TRUE., and therefore activates section 1.1 of each
64	!! subroutine which writes messages to the output. \n
65	!! This subroutine is called at the beginning of ::stomateLpj_clear in the
66	!! ::stomate_lpj module.
67	!!
68	!! RECENT CHANGE(S): None
69	!!
70	!! MAIN OUTPUT VARIABLE(S): ::firstcall, ::firstcall_hum, ::firstcall_moi, ::firstcall_humgdd,
71	!! ::firstcall_moigdd
72	!!
73	!! REFERENCE(S) : None
74	!!
75	!! FLOWCHART : None
76	!! \n
77	!_ ================================================================================================================================
78
79	SUBROUTINE phenology_clear
80	firstcall=.TRUE.
81	firstcall_hum=.TRUE.
82	firstcall_moi = .TRUE.
83	firstcall_humgdd = .TRUE.
84	firstcall_moigdd = .TRUE.
85	END SUBROUTINE phenology_clear
86
87
88	!! ================================================================================================================================
89	!! SUBROUTINE : phenology_diagnostic
90	!!
91	!>\BRIEF This subroutine controls the detection of the beginning of the growing season
92	!! (if dormancy has been long enough), leaf onset, given favourable biometeorological
93	!! conditions, and leaf growth and biomass allocation when leaf biomass is low (i.e.
94	!! at the start of the growing season.
95	!!
96	!! DESCRIPTION : This subroutine is called by the module ::stomate_lpj and deals with the beginning of the
97	!! growing season. First it is established whether the beginning of the growing season is
98	!! allowed. This occurs if the dormance period has been long enough (i.e. greater
99	!! than a minimum PFT-dependent threshold, specified by ::lowgpp_time),
100	!! AND if the last beginning of the growing season was a sufficiently long time ago
101	!! (i.e. when the growing season length is greater than a minimum threshold, specified
102	!! by ::min_growthinit_time, which is defined in this module to be 300 days. \n
103	!! The dormancy time-length is represented by the variable
104	!! ::time_lowgpp, which is calculated in ::stomate_season. It is increased by
105	!! the stomate time step when the weekly GPP is lower than a threshold. Otherwise
106	!! it is set to zero. \n
107	!! ::lowgpp_time is set for each PFT in ::stomate_data from a table of all
108	!! PFT values (::lowgpp_time_tab), which is defined in ::stomate_constants. \n
109	!! The growing season length is given by ::when_growthinit, which increases
110	!! by the stomate time-step at each call to this phenology module, except for when
111	!! leaf onset is detected, when it is set to 0. \n
112	!! If these two conditions are met, leaf onset occurs if the biometeorological
113	!! conditions are also met. This is determined by the leaf onset models, which are
114	!! biome-specific. Each PFT is looped over (ignoring bare soil).
115	!! The onset phenology model is selected, (according to the parameter
116	!! ::pheno_model, which is initialised in stomate_data), and called. \n
117	!! There are six leaf onset phenology models currently being used by ORCHIDEE.
118	!! These are: 'hum' and 'moi', which are based exclusively on moisture conditions,
119	!! 'humgdd' and 'moigdd', which are based on both temperature and moisture conditions,
120	!! 'ncdgdd', which is based on a "chilling" requirement for leaf onset, and
121	!! 'ngd', which is based on the number of growing days since the temperature was
122	!! above a certain threshold, to account for the end of soil frost.
123	!! Those models which are based mostly on temperature conditions are used for
124	!! temperate and boreal biomes, and those which include a moisture condition are used
125	!! for tropical biomes. More detail on the biometeorological conditions is provided
126	!! in the sections on the individual onset models. \n
127	!! The moisture conditions are based on the concept of plant "moisture availability".
128	!! This is based on the soil humidity (relative soil moisture), but is moderated by
129	!! the root density profile, as per the equation:
130	!! \latexonly
131	!! \input{phenology_moiavail_eqn1.tex}
132	!! \endlatexonly
133	!! \n
134	!! Although some studies have shown that the length of the photoperiod is important
135	!! in determining onset (and senescence) dates, this is not considered in the current
136	!! versions of the onset models (Krinner et al., 2005). \n
137	!! If conditions are favourable, leaf onset occurs (::begin_leaves is set to TRUE),
138	!! ::when_growthinit is set to 0.0, and the growing season has begun. \n
139	!! Following the detection of leaf onset, biomass is allocated from the carbohydrate
140	!! reserves equally to the leaves and roots IF the leaf biomass is lower than a minimum
141	!! threshold, which is calculated in this subroutine from the parameter
142	!! ::lai_initmin, divided by the specific leaf area (both of which are
143	!! PFT-dependent and set in ::stomate_constants). \n
144	!! Finally, if biomass is required to be allocated from the carbohydrate reserve
145	!! because the leaf biomass is too low, the leaf age and leaf age distribution is
146	!! re-set. In this case the youngest age class fraction is set to 1 and all other
147	!! leaf age class fractions are set to 0. All leaf ages are set to 0. If there is
148	!! no biomass in the carbohydrate reserve, leaf onset will not occur and the PFT
149	!! will disappear from the grid cell (Krinner et al., 2005). \n
150	!! This subrouting is called in ::stomate_lpj.
151	!!
152	!! RECENT CHANGE(S): None
153	!!
154	!! MAIN OUTPUT VARIABLE(S): ::biomass,
155	!! ::when_growthinit,
156	!! ::leaf age distribution
157	!! ::leaf fraction
158	!!
159	!! REFERENCE(S) :
160	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
161	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
162	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
163	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
164	!!
165	!! FLOWCHART :
166	!! \latexonly
167	!! \includegraphics[scale = 1]{phenology_flowchart.png}
168	!! \endlatexonly
169	!! \n
170	!_ ================================================================================================================================
171
172	SUBROUTINE phenology_diagnostic (npts, dt, PFTpresent, veget_max, tlong_ref, &
173	t2m_month, t2m_week, gpp, maxmoiavail_lastyear, minmoiavail_lastyear, &
174	moiavail_month, moiavail_week, gdd_m5_dormance, gdd_midwinter, ncd_dormance, &
175	ngd_minus5, senescence, time_hum_min, biomass, leaf_frac, &
176	leaf_age, when_growthinit, co2_to_bm)
177
178	!! 0. Variable and parameter declaration
179
180	!! 0.1 Input variables
181
182	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid
183	!! cells (unitless)
184	REAL(r_std), INTENT(in) :: dt !! time step (dt_days)
185	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
186	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_max !! "maximal" coverage fraction of a
187	!! PFT (LAI -> infinity) on ground
188	!! (0-1, unitless)
189	REAL(r_std), DIMENSION(npts), INTENT(in) :: tlong_ref !! "long term" 2 meter reference
190	!! temperatures (K)
191	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures
192	!! (K)
193	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
194	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: gpp !! daily gross primary productivity
195	!! @tex ($gC m^{-2} of
196	!! ground/day$) @endtex
197	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: maxmoiavail_lastyear !! last year's maximum moisture
198	!! availability (0-1, unitless)
199	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: minmoiavail_lastyear !! last year's minimum moisture
200	!! availability (0-1, unitless)
201	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_month !! "monthly" moisture availability
202	!! (0-1, unitless)
203	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_week !! "weekly" moisture availability
204	!! (0-1, unitless)
205	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: gdd_m5_dormance !! growing degree days above a
206	!! threshold of -5 deg C (C)
207	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gdd_midwinter !! growing degree days, since
208	!! midwinter (C)
209	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: ncd_dormance !! number of chilling days since
210	!! leaves were lost (days)
211	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ngd_minus5 !! number of growing days above a
212	!! threshold of -5 deg C (days)
213	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: senescence !! is the plant senescent? (only
214	!! for deciduous trees -
215	!! carbohydrate reserve)
216	!! (true/false)
217	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: time_hum_min !! time elapsed since strongest
218	!! moisture availability (days)
219
220	!
221	!! 0.2 Ouput variables
222	!
223
224	!
225	!! 0.3 Modified variables
226	!
227	REAL(r_std), DIMENSION(npts,nvm,nparts,nelements), &
228	INTENT(inout) :: biomass !! biomass @tex ($gC m^{-2} of
229	!! ground$) @endtex
230	REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_frac !! fraction of leaves in leaf age
231	!! class (0-1, unitless)
232	REAL(r_std), DIMENSION(npts,nvm,nleafages), INTENT(inout) :: leaf_age !! leaf age (days)
233	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: when_growthinit !! how many days since the
234	!! beginning of the growing season
235	!! (days)
236	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: co2_to_bm !! co2 taken up by carbohydrate
237	!! reserve at the beginning of the
238	!! growing season @tex ($gC m^{-2}
239	!! of total ground/day$) @endtex
240
241	!! 0.4 Local variables
242
243	LOGICAL, DIMENSION(npts,nvm) :: allow_initpheno !! are we allowed to decalre the
244	!! beginning of the growing
245	!! season? (true/false)
246	LOGICAL(r_std), DIMENSION(npts) :: age_reset !! does the leaf age distribution
247	!! have to be reset? (true/false)
248	INTEGER(i_std) :: i,j,k,l,m !! indices (unitless)
249	REAL(r_std), DIMENSION(npts) :: bm_wanted !! biomass we would like to have
250	!! @tex ($gC m^{-2} of ground$)
251	!! @endtex
252	REAL(r_std), DIMENSION(npts) :: bm_use !! biomass we use (from
253	!! carbohydrate reserve or from
254	!! atmosphere) @tex ($gC m^{-2} of
255	!! ground$) @endtex
256	REAL(r_std), DIMENSION(npts) :: lm_min !! minimum leaf mass @tex ($gC
257	!! m^{-2} of ground$) @endtex
258	LOGICAL, DIMENSION(npts,nvm) :: begin_leaves !! signal to start putting leaves
259	!! on (true/false)
260	REAL(r_std), DIMENSION(npts,nvm) :: histvar !! controls the history output
261	!! level - 0: nothing is written;
262	!! 10: everything is written
263	!! (0-10, unitless)
264	!_ ================================================================================================================================
265
266	IF (bavard.GE.2) WRITE(numout,*) 'Entering phenology'
267
268	!! 2. first call
269
270	! output message giving the setting of the ::always_init
271	! and ::min_growthinit_time parameters.
272
273	IF ( firstcall ) THEN
274
275	WRITE(numout,*) 'phenology:'
276
277	WRITE(numout,*) ' > take carbon from atmosphere if carbohydrate' // &
278	' reserve too small (::always_init): ', always_init
279
280	WRITE(numout,*) ' > minimum time since last beginning of a growing' // &
281	' season (d) (::min_growthinit_time): ', min_growthinit_time
282
283	firstcall = .FALSE.
284
285	ENDIF
286
287
288	!! 3. Detection of the beginning of the growing season.
289
290	!! 3.1 Allow detection of the beginning of the growing season
291	! if dormance was long enough (i.e. when ::time_lowgpp, which is calculated in ::stomate_season,
292	! is above a certain PFT-dependent threshold, ::lowgpp_time, which is given in ::stomate_constants),
293	! AND the last beginning of growing season was a sufficiently long time ago
294	! (i.e. when ::when_growthinit, which is calculated in this module,
295	! is greater than ::min_growthinit_time, which is declared at the beginning of this module).
296	! If these conditions are met, allow_initpheno is set to TRUE. Each PFT is looped over
297	allow_initpheno(:,1) = .FALSE.
298	DO j = 2,nvm
299
300	WHERE ( when_growthinit(:,j) .GT. min_growthinit_time )
301	allow_initpheno(:,j) = .TRUE.
302	ELSEWHERE
303	allow_initpheno(:,j) = .FALSE.
304	ENDWHERE
305
306	ENDDO
307
308	WHERE(allow_initpheno)
309	histvar=un
310	ELSEWHERE
311	histvar=zero
312	ENDWHERE
313	CALL histwrite_p (hist_id_stomate, 'ALLOW_INITPHENO', itime, histvar, npts*nvm, horipft_index)
314
315	!! 3.2 increase the ::when_growthinit counter, which gives the number of days since the beginning of the growing season.
316	! Needed for allocation and for the detection of the beginning of the growing season.
317	when_growthinit(:,:) = when_growthinit(:,:) + dt
318
319
320	!! 4. Leaf onset
321
322	! Check biometeorological conditions using the onset phenological models,
323	! which are different for each PFT group (i.e. grass versus tropical etc.
324	! See below for more detail on the different models and which PFTs use each model).
325	! By default: phenology does not start (::begin_leaves set to FALSE).
326	begin_leaves(:,:) = .FALSE.
327
328	! - The onset phenology model is selected, (according to the parameter ::pheno_model,
329	! which is initialised in stomate_data), and called.
330	! Each PFT is looped over (ignoring bare soil).
331	! If conditions are favourable, begin_leaves is set to TRUE.
332	! parameter used in all the differents models of phenology
333	t_always = ZeroCelsius + t_always_add
334
335	DO j = 2,nvm ! Loop over # PFTs
336
337	SELECT CASE ( pheno_model(j) )
338
339	CASE ( 'hum' )
340
341	CALL pheno_hum (npts, j, PFTpresent, allow_initpheno, &
342	moiavail_month, moiavail_week, &
343	maxmoiavail_lastyear, minmoiavail_lastyear, &
344	begin_leaves)
345
346	CASE ( 'moi' )
347
348	CALL pheno_moi (npts, j, PFTpresent, allow_initpheno, &
349	time_hum_min, &
350	moiavail_month, moiavail_week, &
351	begin_leaves)
352
353
354	CASE ( 'ncdgdd' )
355
356	CALL pheno_ncdgdd (npts, j, PFTpresent, allow_initpheno, &
357	ncd_dormance, gdd_midwinter, &
358	t2m_month, t2m_week, begin_leaves)
359
360	CASE ( 'ngd' )
361
362	CALL pheno_ngd (npts, j, PFTpresent, allow_initpheno, ncd_dormance, &
363	ngd_minus5, t2m_month, t2m_week, begin_leaves)
364
365	CASE ( 'humgdd' )
366
367	CALL pheno_humgdd (npts, j, PFTpresent, allow_initpheno, gdd_m5_dormance, &
368	maxmoiavail_lastyear, minmoiavail_lastyear, &
369	tlong_ref, t2m_month, t2m_week, &
370	moiavail_week, moiavail_month, &
371	begin_leaves)
372
373	CASE ( 'moigdd' )
374
375	CALL pheno_moigdd (npts, j, PFTpresent, allow_initpheno, gdd_m5_dormance, &
376	time_hum_min, &
377	tlong_ref, t2m_month, t2m_week, &
378	moiavail_week, moiavail_month, &
379	begin_leaves)
380
381	CASE ( 'none' )
382
383	! no action
384
385	CASE default
386
387	WRITE(numout,*) 'phenology: don''t know how to treat this PFT.'
388	WRITE(numout,*) ' number: (::j)',j
389	WRITE(numout,*) ' phenology model (::pheno_model(j)) : ',pheno_model(j)
390	STOP
391
392	END SELECT
393
394	ENDDO
395
396	WHERE(begin_leaves)
397
398	histvar=un
399
400	ELSEWHERE
401
402	histvar=zero
403
404	ENDWHERE
405
406	CALL histwrite_p (hist_id_stomate, 'BEGIN_LEAVES', itime, histvar, npts*nvm, horipft_index)
407
408	!! 5. Leaf growth and biomass allocation when leaf biomass is low.
409
410	! Leaves start to grow if biometeorological conditions are favourable (::begin_leaves == TRUE) and if
411	! leaf growth is allowed (::allow_initpheno == TRUE).
412	! PFTs and then grid cells are looped over.
413	DO j = 2,nvm ! Loop over # PFTs
414
415	age_reset(:) = .FALSE.
416
417	DO i = 1, npts
418
419	IF ( begin_leaves(i,j) ) THEN
420
421	!! 5.1 First minimum biomass is calculated using the following equation:
422	!! \latexonly
423	!! \input{phenology_lm_min_eqn2.tex}
424	!! \endlatexonly
425	!! \n
426	lm_min(i) = lai_initmin(j) / sla(j)
427
428
429	! If leaf biomass is lower than the minimum biomass then biomass must
430	! be allocated from the carbohydrate reserves to leaves and roots.
431	IF ( biomass(i,j,ileaf,icarbon) .LT. lm_min(i) ) THEN
432
433	! Determine how much biomass is available to use
434	! First calculate how much biomass is wanted/required
435	! (::bm_wanted = 2 x the minimum leaf biomass).
436	bm_wanted(i) = 2. * lm_min(i)
437
438	! If the biomass in the carbohydrate reserves is less than the required biomass
439	! take the required amount of carbon from the atmosphere and put it into the
440	! carbohydrate reserve. This only occurs if the parameter ::always_init
441	! (set at beginning of this ::subroutine) is TRUE. Default is FALSE.
442	IF ( always_init .AND. ( biomass(i,j,icarbres,icarbon) .LT. bm_wanted(i) ) ) THEN
443
444	co2_to_bm(i,j) = co2_to_bm(i,j) + ( bm_wanted(i) - biomass(i,j,icarbres,icarbon) ) / dt
445	biomass(i,j,icarbres,icarbon) = bm_wanted(i)
446
447	ENDIF
448
449	! The biomass available to use is set to be the minimum of the biomass of
450	! the carbohydrate reservoir (if carbon not taken from the atmosphere), and
451	! the wanted biomass.
452	bm_use(i) = MIN( biomass(i,j,icarbres,icarbon), bm_wanted(i) )
453
454	! Divide the biomass which is available to use equally between the leaves
455	! and roots.
456	biomass(i,j,ileaf,icarbon) = biomass(i,j,ileaf,icarbon) + bm_use(i) / 2.
457	biomass(i,j,iroot,icarbon) = biomass(i,j,iroot,icarbon) + bm_use(i) / 2.
458
459	! Decrease carbohydrate reservoir biomass by the amount that's been allocated
460	! to the leaves and roots
461	biomass(i,j,icarbres,icarbon) = biomass(i,j,icarbres,icarbon) - bm_use(i)
462
463	! set reset leaf age distribution (::age_reset) flag. Default is TRUE.
464	! done later for better vectorization)
465	age_reset(i) = .TRUE.
466
467	ENDIF ! leaf mass is very low
468
469	!! 5.2 reset when_growthinit counter: start of the growing season
470	when_growthinit(i,j) = zero
471
472	ENDIF ! start of the growing season
473
474	ENDDO ! loop over grid points
475
476	!! 5.3 reset leaf age distribution where necessary (i.e. when age_reset is TRUE)
477	! simply say that everything is in the youngest age class
478	! Set the youngest age class fraction to 1 and all other leaf age class fractions to 0.
479	WHERE ( age_reset(:) )
480
481	leaf_frac(:,j,1) = un
482
483	ENDWHERE
484
485	DO m = 2, nleafages
486
487	WHERE ( age_reset(:) )
488
489	leaf_frac(:,j,m) = zero
490
491	ENDWHERE
492
493	ENDDO ! nleafages
494
495	! Ages - set all leaf ages to 0.
496	DO m = 1, nleafages
497
498	WHERE ( age_reset(:) )
499
500	leaf_age(:,j,m) = zero
501
502	ENDWHERE
503
504	ENDDO ! nleafages
505
506	ENDDO ! loop over # PFTs
507
508	IF (bavard.GE.2) WRITE(numout,*) 'Leaving phenology'
509
510	END SUBROUTINE phenology_diagnostic
511
512
513	!! ================================================================================================================================
514	!! SUBROUTINE : phenology_prognostic
515	!!
516	!>\BRIEF This subroutine controls the detection of the beginning of the growing season
517	!! (if dormancy has been long enough), leaf onset, given favourable biometeorological
518	!! conditions, and leaf growth and biomass allocation when leaf biomass is low (i.e.
519	!! at the start of the growing season.
520	!!
521	!! DESCRIPTION : This subroutine is called by the module ::stomate_lpj and deals with the beginning of the
522	!! growing season. First it is established whether the beginning of the growing season is
523	!! allowed. This occurs if the dormance period has been long enough (i.e. greater
524	!! than a minimum PFT-dependent threshold, specified by ::lowgpp_time),
525	!! AND if the last beginning of the growing season was a sufficiently long time ago
526	!! (i.e. when the growing season length is greater than a minimum threshold, specified
527	!! by ::min_growthinit_time, which is defined in this module to be 300 days. \n
528	!! The dormancy time-length is represented by the variable
529	!! ::time_lowgpp, which is calculated in ::stomate_season. It is increased by
530	!! the stomate time step when the weekly GPP is lower than a threshold. Otherwise
531	!! it is set to zero. \n
532	!! ::lowgpp_time is set for each PFT in ::stomate_data from a table of all
533	!! PFT values (::lowgpp_time_tab), which is defined in ::stomate_constants. \n
534	!! The growing season length is given by ::when_growthinit, which increases
535	!! by the stomate time-step at each call to this phenology module, except for when
536	!! leaf onset is detected, when it is set to 0. \n
537	!! If these two conditions are met, leaf onset occurs if the biometeorological
538	!! conditions are also met. This is determined by the leaf onset models, which are
539	!! biome-specific. Each PFT is looped over (ignoring bare soil).
540	!! The onset phenology model is selected, (according to the parameter
541	!! ::pheno_model, which is initialised in stomate_data), and called. \n
542	!! There are six leaf onset phenology models currently being used by ORCHIDEE.
543	!! These are: 'hum' and 'moi', which are based exclusively on moisture conditions,
544	!! 'humgdd' and 'moigdd', which are based on both temperature and moisture conditions,
545	!! 'ncdgdd', which is based on a "chilling" requirement for leaf onset, and
546	!! 'ngd', which is based on the number of growing days since the temperature was
547	!! above a certain threshold, to account for the end of soil frost.
548	!! Those models which are based mostly on temperature conditions are used for
549	!! temperate and boreal biomes, and those which include a moisture condition are used
550	!! for tropical biomes. More detail on the biometeorological conditions is provided
551	!! in the sections on the individual onset models. \n
552	!! The moisture conditions are based on the concept of plant "moisture availability".
553	!! This is based on the soil humidity (relative soil moisture), but is moderated by
554	!! the root density profile, as per the equation:
555	!! \latexonly
556	!! \input{phenology_moiavail_eqn1.tex}
557	!! \endlatexonly
558	!! \n
559	!! Although some studies have shown that the length of the photoperiod is important
560	!! in determining onset (and senescence) dates, this is not considered in the current
561	!! versions of the onset models (Krinner et al., 2005). \n
562	!! If conditions are favourable, leaf onset occurs (::begin_leaves is set to TRUE),
563	!! ::when_growthinit is set to 0.0, and the growing season has begun. \n
564	!! Following the detection of leaf onset, biomass is allocated from the carbohydrate
565	!! reserves equally to the leaves and roots IF the leaf biomass is lower than a minimum
566	!! threshold, which is calculated in this subroutine from the parameter
567	!! ::lai_initmin, divided by the specific leaf area (both of which are
568	!! PFT-dependent and set in ::stomate_constants). \n
569	!! Finally, if biomass is required to be allocated from the carbohydrate reserve
570	!! because the leaf biomass is too low, the leaf age and leaf age distribution is
571	!! re-set. In this case the youngest age class fraction is set to 1 and all other
572	!! leaf age class fractions are set to 0. All leaf ages are set to 0. If there is
573	!! no biomass in the carbohydrate reserve, leaf onset will not occur and the PFT
574	!! will disappear from the grid cell (Krinner et al., 2005). \n
575	!! This subrouting is called in ::stomate_lpj.
576	!!
577	!! RECENT CHANGE(S): None
578	!!
579	!! MAIN OUTPUT VARIABLE(S): ::biomass,
580	!! ::when_growthinit,
581	!! ::leaf age distribution
582	!! ::leaf fraction
583	!!
584	!! REFERENCE(S) :
585	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
586	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
587	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
588	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
589	!!
590	!! FLOWCHART :
591	!! \latexonly
592	!! \includegraphics[scale = 1]{phenology_flowchart.png}
593	!! \endlatexonly
594	!! \n
595	!_ ================================================================================================================================
596
597	SUBROUTINE phenology_prognostic (npts, dt, PFTpresent, veget_max, &
598	tlong_ref, t2m_month, t2m_week, gpp, &
599	maxmoiavail_lastyear, minmoiavail_lastyear, moiavail_month, moiavail_week, &
600	gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, &
601	senescence, time_hum_min, biomass, leaf_frac, &
602	leaf_age, when_growthinit, co2_to_bm, circ_class_n, &
603	circ_class_biomass, ind, KF, allow_initpheno, &
604	tau_eff_leaf, tau_eff_sap, tau_eff_root, age, &
605	everywhere, npp_longterm, lm_lastyearmax, k_latosa_adapt)
606
607
608	!! 0. Variable and parameter declaration
609
610
611	!! 0.1 Input variables
612
613	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells (unitless)
614	REAL(r_std), INTENT(in) :: dt !! time step (dt_days)
615	REAL(r_std), DIMENSION(:), INTENT(in) :: tlong_ref !! "long term" 2 meter reference temperatures (K)
616	REAL(r_std), DIMENSION(:), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures (K)
617	REAL(r_std), DIMENSION(:), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
618	REAL(r_std), DIMENSION(:,:), INTENT(in) :: veget_max !! "maximal" coverage fraction of a
619	!! PFT (LAI -> infinity) on ground
620	!! (0-1, unitless)
621	REAL(r_std), DIMENSION(:,:), INTENT(in) :: gpp !! daily gross primary productivity
622	!! @tex ($gC m^{-2} of
623	!! ground/day$) @endtex
624	REAL(r_std), DIMENSION(:,:), INTENT(in) :: maxmoiavail_lastyear !! last year's maximum moisture
625	!! availability (0-1, unitless)
626	REAL(r_std), DIMENSION(:,:), INTENT(in) :: minmoiavail_lastyear !! last year's minimum moisture
627	!! availability (0-1, unitless)
628	REAL(r_std), DIMENSION(:,:), INTENT(in) :: moiavail_month !! "monthly" moisture availability
629	!! (0-1, unitless)
630	REAL(r_std), DIMENSION(:,:), INTENT(in) :: moiavail_week !! "weekly" moisture availability
631	!! (0-1, unitless)
632	REAL(r_std), DIMENSION(:,:), INTENT(in) :: gdd_m5_dormance !! growing degree days above a
633	!! threshold of -5 deg C (C)
634	REAL(r_std), DIMENSION(:,:), INTENT(in) :: ncd_dormance !! number of chilling days since
635	!! leaves were lost (days)
636	REAL(r_std), DIMENSION(:,:), INTENT(in) :: tau_eff_root !! Effective root turnover time that accounts
637	!! waterstress (days)
638	REAL(r_std), DIMENSION(:,:), INTENT(in) :: tau_eff_sap !! Effective sapwood turnover time that accounts
639	!! waterstress (days)
640	REAL(r_std), DIMENSION(:,:), INTENT(in) :: tau_eff_leaf !! Effective leaf turnover time that accounts
641	!! waterstress (days)
642
643	!! 0.2 Ouput variables
644
645	LOGICAL, DIMENSION(:,:), INTENT(out) :: allow_initpheno !! are we allowed to declare the
646	!! beginning of the growing season?
647	!! and the end of senescence (true/false)
648
649	!! 0.3 Modified variables
650
651	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: gdd_midwinter !! growing degree days, since midwinter (C)
652	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: when_growthinit !! how many days since the
653	!! beginning of the growing season (days)
654	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: co2_to_bm !! co2 taken up by carbohydrate
655	!! reserve at the beginning of the
656	!! growing season @tex ($gC m^{-2}
657	!! of total ground/day$) @endtex
658	LOGICAL, DIMENSION(:,:), INTENT(inout) :: PFTpresent !! PFT exists (true/false)
659	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: ngd_minus5 !! number of growing days above a
660	!! threshold of -5 deg C (days)
661	LOGICAL, DIMENSION(:,:), INTENT(inout) :: senescence !! is the plant senescent? (only
662	!! for deciduous trees -
663	!! carbohydrate reserve) (true/false)
664	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: ind !! Stand level number of individuals
665	!! @tex $(ind m^{-2})$ @endtex
666	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: everywhere !! is the PFT everywhere in the grid box or
667	!! very localized (after its introduction) (?)
668	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: time_hum_min !! time elapsed since strongest
669	!! moisture availability (days)
670	REAL(r_std), DIMENSION(:,:,:), INTENT(inout) :: circ_class_n !! Number of individuals in each circ class
671	!! @tex $(ind m^{-2})$ @endtex
672	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: KF !! Scaling factor to convert sapwood mass
673	!! into leaf mass (m)
674	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: k_latosa_adapt !! Leaf to sapwood area adapted for long
675	!! term water stress (m)
676	REAL(r_std), DIMENSION(:,:,:), INTENT(inout) :: leaf_frac !! fraction of leaves in leaf age
677	!! class (0-1, unitless)
678	REAL(r_std), DIMENSION(:,:,:), INTENT(inout) :: leaf_age !! leaf age (days)
679	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: age !! mean age (years)
680	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: npp_longterm !! "long term" net primary productivity
681	!! @tex ($gC m^{-2} year^{-1}$) @endtex
682	REAL(r_std), DIMENSION(:,:), INTENT(inout) :: lm_lastyearmax !! last year's maximum leaf mass for each PFT
683	!! @tex ($gC m^{-2}$) @endtex
684	REAL(r_std), DIMENSION(:,:,:,:), INTENT(inout) :: biomass !! biomass @tex ($gC m^{-2} of
685	!! ground$) @endtex
686	REAL(r_std), DIMENSION(:,:,:,:,:), INTENT(inout) :: circ_class_biomass !! Biomass components of the model tree
687	!! within a circumference class
688	!! class @tex $(g C ind^{-1})$ @endtex
689
690	!! 0.4 Local variables
691
692	REAL(r_std), DIMENSION(npts,nvm) :: c0_alloc !! Root to sapwood tradeoff parameter
693	LOGICAL(r_std), DIMENSION(npts) :: age_reset !! does the leaf age distribution
694	!! have to be reset? (true/false)
695	LOGICAL, DIMENSION(npts,nvm) :: begin_leaves !! signal to start putting leaves
696	!! on (true/false)
697	INTEGER(i_std) :: i,j,k,l,m !! indices (unitless)
698	INTEGER(i_std) :: ipts, ivm, ipar !! indices (unitless)
699	INTEGER(i_std) :: imbc, iele !! indices (unitless)
700	REAL(r_std) :: deficit !! Carbon that needs to be taken from the
701	!! labile and or reserve pools
702	!! @tex $(gC m^{-2})$ @endtex
703	REAL(r_std), DIMENSION(npts,nvm) :: LF !! Scaling factor to convert sapwood mass
704	!! into root mass (unitless)
705	REAL(r_std), DIMENSION(npts,nvm) :: histvar !! controls the history output
706	!! level - 0: nothing is written;
707	!! 10: everything is written
708	!! (0-10, unitless)
709	!!$ REAL(r_std), DIMENSION(npts,nvm) :: lstress_fac !! Light stress factor, based on total
710	!!$ !! transmitted light (unitless, 0-1)
711	REAL(r_std) :: bm_wanted !! biomass we would like to have
712	!! @tex ($gC m^{-2} of ground$)
713	!! @endtex
714	REAL(r_std) :: bm_use !! biomass we use (from
715	!! carbohydrate reserve or from
716	!! atmosphere) @tex ($gC m^{-2} of
717	!! ground$) @endtex
718	REAL(r_std), DIMENSION(ncirc) :: Cl_tree !! Individual plant, leaf compartment
719	!! @tex $(gC tree^{-1})$ @endtex
720	REAL(r_std), DIMENSION(ncirc) :: Cr_tree !! Individual plant, root compartment
721	!! @tex $(gC tree^{-1})$ @endtex
722	REAL(r_std), DIMENSION(ncirc) :: Cs_tree !! Individual plant, sapwood compartment
723	!! @tex $(gC. tree^{-1})$ @endtex
724	REAL(r_std) :: Cs_grass !! Individual plant, sapwood compartment
725	!! @tex $(gC. ind^{-1})$ @endtex
726	REAL(r_std) :: Cl_init !! Initial leaf carbon required to start
727	!! the growing season
728	!! @tex $(gC tree^{-1})$ @endtex
729	REAL(r_std) :: Cr_init !! Initial root carbon required to start
730	!! the growing season
731	!! @tex $(gC tree^{-1})$ @endtex
732	REAL(r_std) :: lm_min !! minimum leaf mass @tex ($gC
733	!! m^{-2} of ground$) @endtex
734	REAL(r_std), DIMENSION(npts,nvm,ncirc) :: height_eff !! Effective tree height calculated from allometric
735	!! relationships (m) but total biomass (for
736	!! allocation calculations)
737	REAL(r_std), DIMENSION(npts,nvm,nmbcomp,nelements) :: check_intern !! Contains the components of the internal
738	!! mass balance chech for this routine
739	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
740	REAL(r_std), DIMENSION(npts,nvm,nelements) :: closure_intern !! Check closure of internal mass balance
741	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
742	REAL(r_std), DIMENSION(npts,nvm,nelements) :: pool_start !! Start and end pool of this routine
743	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
744	REAL(r_std), DIMENSION(npts,nvm,nelements) :: pool_end !! Start and end pool of this routine
745	!! @tex $(gC pixel^{-1} dt^{-1})$ @endtex
746	REAL(r_std) :: temp_total_biomass !! Total biomass of a PFT/grid
747	!! @tex $(gC m^{-2}$) @endtex
748	REAL(r_std) :: temp_class_biomass !! Total biomass of a circ_class
749	!! @tex $(gC m^{-2}$) @endtex
750	REAL(r_std) :: temp_share !! A weighting factor based on the biomass in a
751	!! circ class @tex - @endtex
752	REAL(r_std) :: bm_use_circ !! biomass we use for a given circ class
753	!! $) @endtex
754	!_ ================================================================================================================================
755
756	IF (bavard.GE.2) WRITE(numout,*) 'Entering phenology'
757
758	!! 1. Write current values to history file
759
760	! Current values for ::when_growthinit
761	CALL xios_orchidee_send_field("WHEN_GROWTHINIT",when_growthinit)
762	CALL histwrite_p (hist_id_stomate, 'WHEN_GROWTHINIT', itime, when_growthinit, npts*nvm, horipft_index)
763
764	! Set and write values for ::PFTpresent
765	WHERE(PFTpresent)
766	histvar=un
767	ELSEWHERE
768	histvar=zero
769	ENDWHERE
770	CALL xios_orchidee_send_field("PFTPRESENT",histvar)
771	CALL histwrite_p (hist_id_stomate, 'PFTPRESENT', itime, histvar, npts*nvm, horipft_index)
772
773	! Set and write values for gdd_midwinter
774	WHERE(gdd_midwinter.EQ.undef)
775	histvar=val_exp
776	ELSEWHERE
777	histvar=gdd_midwinter
778	ENDWHERE
779	CALL xios_orchidee_send_field("GDD_MIDWINTER",histvar)
780	CALL histwrite_p (hist_id_stomate, 'GDD_MIDWINTER', itime, histvar, npts*nvm, horipft_index)
781
782	! Set and write values for gdd_m5_dormance
783	WHERE(gdd_m5_dormance.EQ.undef)
784	histvar=val_exp
785	ELSEWHERE
786	histvar=gdd_m5_dormance
787	ENDWHERE
788	CALL xios_orchidee_send_field('GDD_M5_DORMANCE',histvar)
789	CALL histwrite_p (hist_id_stomate, 'GDD_M5_DORMANCE', itime, histvar, npts*nvm, horipft_index)
790
791	! Set and write values for ncd_dormance
792	WHERE(ncd_dormance.EQ.undef)
793	histvar=val_exp
794	ELSEWHERE
795	histvar=ncd_dormance
796	ENDWHERE
797	CALL xios_orchidee_send_field("NCD_DORMANCE",histvar)
798	CALL histwrite_p (hist_id_stomate, 'NCD_DORMANCE', itime, histvar, npts*nvm, horipft_index)
799
800
801	!! 2. first call
802
803	!! 2.1 Initialize check for mass balance closure
804	! The mass balance is calculated at the end of this routine
805	! in section 6
806	! Initial biomass pool
807	pool_start(:,:,:) = zero
808	DO ipar = 1,nparts
809	DO iele = 1,nelements
810	pool_start(:,:,iele) = pool_start(:,:,iele) + &
811	(biomass(:,:,ipar,iele) * veget_max(:,:))
812	ENDDO
813	ENDDO
814
815	! co2_to_bm is defined as intent inout, the variable accumulates
816	! carbon over the course of a day. Use the difference between
817	! start and the end of this routine
818	check_intern(:,:,iatm2land,icarbon) = - un * co2_to_bm(:,:) * veget_max(:,:) * dt
819
820	!! 2.2 Output messages
821	! giving the setting of the ::always_init
822	! and ::min_growthinit_time parameters.
823	IF ( firstcall ) THEN
824
825	WRITE(numout,*) 'phenology:'
826	WRITE(numout,*) ' > take carbon from atmosphere if carbohydrate' // &
827	' reserve too small (::always_init): ', always_init
828	WRITE(numout,*) ' > minimum time since last beginning of a growing' // &
829	' season (d) (::min_growthinit_time): ', min_growthinit_time
830	firstcall = .FALSE.
831
832	ENDIF
833
834
835	!! 3. Detection of the beginning of the growing season
836
837	!! 3.1 allow detection of the beginning of the growing season if dormance was
838	! long enough (i.e. when ::time_lowgpp, which is calculated in ::stomate_season,
839	! is above a certain PFT-dependent threshold, ::lowgpp_time,
840	! which is given in ::stomate_constants),
841	! AND the last beginning of growing season was a sufficiently long time ago
842	! (i.e. when ::when_growthinit, which is calculated in this module,
843	! is greater than ::min_growthinit_time, which is declared at the beginning of this module).
844	! If these conditions are met, allow_initpheno is set to TRUE. Each PFT is looped over.
845	allow_initpheno(:,1) = .FALSE.
846
847	DO j = 2,nvm
848
849	WHERE ( when_growthinit(:,j) .GT. min_growthinit_time )
850
851	! Change the values of ::allow_initpheno and ::senescence
852	allow_initpheno(:,j) = .TRUE.
853	senescence (:,j) = .FALSE.
854
855	ELSEWHERE
856
857	! Keep the values of ::allow_initpheno and ::senescence
858	allow_initpheno(:,j) = .FALSE.
859
860	ENDWHERE
861
862	ENDDO
863
864	!+++TEMP+++
865	! SHORT CIRCUIT phenology
866	!!$ IF ( firstcall ) THEN
867	!!$ allow_initpheno(:,j) = .TRUE.
868	!!$ ENDIF
869	!++++++++++++
870
871	WHERE(allow_initpheno)
872
873	histvar=un
874
875	ELSEWHERE
876
877	histvar=zero
878
879	ENDWHERE
880
881	CALL histwrite_p (hist_id_stomate, 'ALLOW_INITPHENO', itime, histvar, npts*nvm, horipft_index)
882
883
884	!! 3.2 increase the ::when_growthinit counter, which gives the number of days since the beginning
885	! of the growing season.
886	! Needed for allocation and for the detection of the beginning of the growing season.
887	when_growthinit(:,:) = when_growthinit(:,:) + dt
888
889	! Keep when_wrowthinit undefined for bare soil
890	when_growthinit(:,1) = undef
891
892
893	!! 4. Leaf onset
894
895	! Check biometeorological conditions using the onset phenological models,
896	! which are different for each PFT group (i.e. grass versus tropical etc.
897	! See below for more detail on the different models and which PFTs use each model).
898	! - By default: phenology does not start (::begin_leaves set to FALSE).
899	begin_leaves(:,:) = .FALSE.
900
901	! - The onset phenology model is selected, (according to the parameter ::pheno_model,
902	! which is initialised in stomate_data), and called.
903	! Each PFT is looped over (ignoring bare soil).
904	! If conditions are favourable, begin_leaves is set to TRUE.
905	! parameter used in all the differents models of phenology
906	t_always = ZeroCelsius + t_always_add
907
908	DO j = 2,nvm ! Loop over # PFTs
909
910	SELECT CASE ( pheno_model(j) )
911
912	CASE ( 'hum' )
913
914	CALL pheno_hum (npts, j, PFTpresent, allow_initpheno, &
915	moiavail_month, moiavail_week, maxmoiavail_lastyear, minmoiavail_lastyear, &
916	begin_leaves)
917
918	CASE ( 'moi' )
919
920	CALL pheno_moi (npts, j, PFTpresent, allow_initpheno, &
921	time_hum_min, moiavail_month, moiavail_week, begin_leaves)
922
923	CASE ( 'ncdgdd' )
924
925	CALL pheno_ncdgdd (npts, j, PFTpresent, allow_initpheno, &
926	ncd_dormance, gdd_midwinter, t2m_month, t2m_week, &
927	begin_leaves)
928
929	CASE ( 'ngd' )
930
931	CALL pheno_ngd (npts, j, PFTpresent, allow_initpheno, ncd_dormance, &
932	ngd_minus5, t2m_month, t2m_week, begin_leaves)
933
934	CASE ( 'humgdd' )
935
936	CALL pheno_humgdd (npts, j, PFTpresent, allow_initpheno, &
937	gdd_m5_dormance, maxmoiavail_lastyear, minmoiavail_lastyear, tlong_ref, &
938	t2m_month, t2m_week, moiavail_week, moiavail_month, &
939	begin_leaves)
940
941	CASE ( 'moigdd' )
942
943	CALL pheno_moigdd (npts, j, PFTpresent, allow_initpheno, &
944	gdd_m5_dormance, time_hum_min, tlong_ref, t2m_month, &
945	t2m_week, moiavail_week, moiavail_month, begin_leaves)
946
947	CASE ( 'none' )
948
949	! no action
950
951	CASE default
952
953	WRITE(numout,*) ' phenology: don''t know how to treat this PFT.'
954	WRITE(numout,*) ' number: (::j)',j
955	WRITE(numout,*) ' phenology model (::pheno_model(j)) : ',pheno_model(j)
956	STOP
957
958	END SELECT
959
960	ENDDO
961
962
963	WHERE(begin_leaves)
964
965	histvar=un
966
967	ELSEWHERE
968
969	histvar=zero
970
971	ENDWHERE
972
973	CALL histwrite_p (hist_id_stomate, 'BEGIN_LEAVES', itime, histvar, npts*nvm, horipft_index)
974
975	!+++TEMP+++
976	!!$ ! SHORT CIRCUIT phenology
977	!!$ IF ( firstcall ) THEN
978	!!$ begin_leaves(:,:) = .TRUE.
979	!!$ WRITE(numout,*) 'begin leaves, ', begin_leaves(:,:)
980	!!$ ! when deleting this line remember to uncomment the firstcall above
981	!!$ firstcall = .FALSE.
982	!!$ ENDIF
983	!++++++++++++
984
985
986	!! 5. Leaf growth and biomass allocation when leaf biomass is low.
987
988	! Leaves start to grow if biometeorological conditions are
989	! favourable (::begin_leaves == TRUE) and if
990	! leaf growth is allowed (::allow_initpheno == TRUE).
991	! PFTs and then grid cells are looped over.
992	DO j = 2,nvm ! Loop over # PFTs
993
994	age_reset(:) = .FALSE.
995	bm_use = zero
996	bm_wanted = zero
997
998	DO i = 1,npts
999
1000	! We might need the c0_alloc factor, so let's
1001	! calculate it now.
1002	c0_alloc(i,j)=calculate_c0_alloc(i, j, tau_eff_root(i,j), &
1003	tau_eff_sap(i,j))
1004
1005	IF ( begin_leaves(i,j) ) THEN
1006
1007	!! 4.1 First minimum biomass is calculated using the following equation:
1008	! \latexonly
1009	! \input{phenology_lm_min_eqn2.tex}
1010	! \endlatexonly
1011	! \n
1012	lm_min = lai_initmin(j) / sla(j)
1013
1014	! The minimum leaf biomass is prescribed by ::lm_min which in turn
1015	! is basically prescribed through ::lai_initmin. However, lm_min could exceed
1016	! the leaf mass that is required to respect the allometric relationships
1017	! therefore this leaf biomass should be calculated
1018
1019	! Calculate the allocation factors see stomate_prsecribe.f90 and
1020	! stomate_prescribe.f90 for more details
1021	LF(i,j) = c0_alloc(i,j) * KF(i,j)
1022
1023	IF ( is_tree(j) ) THEN
1024
1025	! Calculate the stand structure (gC tree-1)
1026	Cs_tree(:) = ( circ_class_biomass(i,j,:,isapabove,icarbon) + &
1027	circ_class_biomass(i,j,:,isapbelow,icarbon) )
1028	height_eff(i,j,:) = wood_to_height_eff(circ_class_biomass(i,j,:,:,icarbon),j)
1029
1030	! Calculate the leaves and roots that need to be grown. See
1031	! stomate_growth_fun_alloc.f90 for more details
1032	! Note that the units for Cl_init and Cr_init are gC m-2
1033	Cl_init = MAX(zero, MIN( lm_min, SUM( (KF(i,j) * Cs_tree(:) / height_eff(i,j,:) - &
1034	circ_class_biomass(i,j,:,ileaf,icarbon)) * circ_class_n(i,j,:) ) ) )
1035	Cr_init = MAX(zero, (biomass(i,j,ileaf,icarbon) + Cl_init) / LF(i,j) - &
1036	biomass(i,j,iroot,icarbon) )
1037
1038
1039	! Grasses
1040	! Only used when grasses are described as a mixed phenology
1041	! in the DOFOCO branch grasses are simulated as an evergreen
1042	! biome
1043	ELSEIF (.NOT. is_tree(j) .AND. natural(j) ) THEN
1044
1045	! Calculate the available biomass in sapwood (gC m-2)
1046	Cs_grass = biomass(i,j,isapabove,icarbon)
1047
1048	! Calculate structure of the crop/grassland i.e. the leaves and roots
1049	! that need to be grown see stomate_growth_fun_alloc.f90 for more details
1050	! Note that the units for Cl_init and Cr_init are gC m-2. For grasses
1051	! :: heigh_init(j) is the minimal height_eff and is therefore used at the start of
1052	! growing season (it could be considered the height_eff of Cs)
1053	Cl_init = MAX(zero, MAX( lm_min, Cs_grass * KF(i,j) - &
1054	biomass(i,j,ileaf,icarbon) ) )
1055	Cr_init = MAX(zero, (biomass(i,j,ileaf,icarbon) + Cl_init) / LF(i,j) - &
1056	biomass(i,j,iroot,icarbon))
1057
1058	!+++TEMP+++
1059	IF (j.EQ.test_pft .AND. ld_pheno) THEN
1060	WRITE(numout,*) 'Cs_grass, ', biomass(i,j,isapabove,icarbon)
1061	WRITE(numout,) 'Cl_init', lm_min, biomass(i,j,ileaf,icarbon), Cs_grass KF(i,j) - &
1062	biomass(i,j,ileaf,icarbon)
1063	WRITE(numout,*) 'Cr_init, ', biomass(i,j,ileaf,icarbon), Cl_init, LF(i,j), biomass(i,j,iroot,icarbon), &
1064	(biomass(i,j,ileaf,icarbon) + Cl_init) / LF(i,j) - &
1065	biomass(i,j,iroot,icarbon)
1066	ENDIF
1067	!++++++++++
1068
1069	ENDIF
1070
1071	! Finalize trees and grasslands which are both considered to be natural.
1072	! Because grasses and trees live on during dormance they will consume
1073	! reserves and turnover may continue. There is no guarantee that the
1074	! required amounts of C for phenology are present and can be used. This
1075	! is now being checked.
1076	IF ( natural(j) ) THEN
1077
1078	! If leaf biomass is lower than the minimum biomass then biomass must
1079	! be allocated from the labile pool to leaves and roots.
1080	IF ( biomass(i,j,ileaf,icarbon) .LT. Cl_init ) THEN
1081
1082	! First calculate how much biomass is wanted/required
1083	bm_wanted = Cl_init + Cr_init
1084
1085	! Specific setting for forced phenology
1086	! If the biomass in the carbohydrate reserves is less than the required biomass
1087	! take the required amount of carbon from the atmosphere and put it into the
1088	! labile pool. This only occurs if the parameter ::always_init
1089	! (set at beginning of this ::subroutine) is TRUE. Default is FALSE.
1090	IF ( always_init .AND. ( biomass(i,j,ilabile,icarbon) .LT. bm_wanted ) ) THEN
1091
1092	co2_to_bm(i,j) = co2_to_bm(i,j) + ( bm_wanted - biomass(i,j,ilabile,icarbon) ) / dt
1093	biomass(i,j,ilabile,icarbon) = biomass(i,j,ilabile,icarbon) + bm_wanted
1094
1095	ENDIF
1096
1097	! The biomass available to use is set to be the minimum of the biomass of
1098	! the labile pool (if carbon not taken from the atmosphere), and
1099	! the wanted biomass. Convert biomass to gc tree-1 so the allometric
1100	! relationships can be applied
1101	bm_use = MIN( biomass(i,j,ilabile,icarbon) + biomass(i,j,icarbres,icarbon), &
1102	bm_wanted ) / ind(i,j)
1103
1104	! This bm_use is for the whole stand. We need to determine how much
1105	! to distribute for each circ class. We do this based on the
1106	! allometry for each circ class.
1107	! Total biomass across parts and circumference classes
1108	temp_total_biomass = zero
1109
1110	DO l = 1,ncirc
1111
1112	DO k = 1,nparts
1113
1114	temp_total_biomass = temp_total_biomass + &
1115	circ_class_biomass(i,j,l,k,icarbon) * circ_class_n(i,j,l)
1116
1117	ENDDO
1118
1119	ENDDO
1120
1121	! zero the root and leaf biomass pools, since we'll sync them
1122	! to circ_class_biomass. Do the same for the labile and
1123	! the reserve pools
1124	biomass(i,j,ileaf,icarbon) = zero
1125	biomass(i,j,iroot,icarbon) = zero
1126	biomass(i,j,icarbres,icarbon) = zero
1127	biomass(i,j,ilabile,icarbon) = zero
1128
1129	DO l=1,ncirc
1130
1131	temp_class_biomass = zero
1132
1133	DO k = 1,nparts
1134
1135	temp_class_biomass = temp_class_biomass + &
1136	circ_class_biomass(i,j,l,k,icarbon) * circ_class_n(i,j,l)
1137
1138	ENDDO
1139
1140	! If there is no biomass in this circ class, we'll skip it
1141	! This should ALWAYS be the case for ncirc > 1 for grasses.
1142	IF (temp_total_biomass .EQ. zero) CYCLE
1143
1144	! Share of this circumference class to the total biomass
1145	temp_share = temp_class_biomass / temp_total_biomass
1146
1147	bm_use_circ = bm_use * temp_share
1148
1149	! bm_use_circ is g C / tree for this circumference class. These
1150	! are the exact same units as circ_class_biomass.
1151
1152	! bm_use may differ from Cr_init and Cl_init, the final allocation will depend
1153	! on bm_use. Distribute bm_use over leaves and roots following allometric
1154	! relationships. Cl_init and Cr_init take the same units as bm_use
1155	! i.e. gC tree-1
1156	! (i) bm_use = Cl_incp + Cr_incp
1157	! (ii) Cr_incp = (Cl_incp+Cl)/LF - Cr
1158	! Substitue (ii) in (i) and solve for Cl_inc
1159	! <=> Cl_incp = (LF*(b_incp+Cr)-Cl)/(1+LF)
1160	! Because this is the start of the growing season, Cr = Cl = 0
1161	Cl_init = ((LF(i,j) * bm_use_circ) - circ_class_biomass(i,j,l,ileaf,icarbon)) &
1162	/ (un + LF(i,j))
1163	Cr_init = (Cl_init + circ_class_biomass(i,j,l,ileaf,icarbon)) / LF(i,j)
1164
1165	! the sum of Cl_init and Cr_init should equal bm_use_circ
1166	IF (bm_use_circ - Cl_init - Cr_init .LT. -min_stomate .OR. &
1167	bm_use_circ - Cl_init - Cr_init .GT. min_stomate) THEN
1168	WRITE(numout,*) 'over or underspening in phenology, ', &
1169	bm_use_circ - Cl_init - Cr_init
1170	STOP
1171	ENDIF
1172
1173	! Decrease labile pool biomass by the amount that's been allocated
1174	! to the leaves and roots. If the labile pool is depleted use carbon
1175	! from the reserve pool.
1176	deficit = bm_use_circ - circ_class_biomass(i,j,l,ilabile,icarbon)
1177
1178	! There is enough carbon in the labile pool
1179	IF (deficit .LT. zero) THEN
1180
1181	circ_class_biomass(i,j,l,ilabile,icarbon) = &
1182	circ_class_biomass(i,j,l,ilabile,icarbon) - bm_use_circ
1183
1184	! Deplete the labile pool, use the reserve pool
1185	ELSE
1186
1187	circ_class_biomass(i,j,l,ilabile,icarbon) = zero
1188	circ_class_biomass(i,j,l,icarbres,icarbon) = &
1189	circ_class_biomass(i,j,l,icarbres,icarbon) - deficit
1190
1191	ENDIF
1192
1193	! Distribute the biomass over the leaves and roots (gC tree-1)
1194	! Since this whole loop is already over circ class, Cl_init and
1195	! Cr_init are exactly what we need.
1196	circ_class_biomass(i,j,l,ileaf,icarbon) = &
1197	circ_class_biomass(i,j,l,ileaf,icarbon) + &
1198	Cl_init
1199	circ_class_biomass(i,j,l,iroot,icarbon) = &
1200	circ_class_biomass(i,j,l,iroot,icarbon) + &
1201	Cr_init
1202
1203	! Here we sync the biomass
1204	biomass(i,j,ileaf,icarbon) = biomass(i,j,ileaf,icarbon) + &
1205	circ_class_biomass(i,j,l,ileaf,icarbon) * circ_class_n(i,j,l)
1206	biomass(i,j,iroot,icarbon) = biomass(i,j,iroot,icarbon) + &
1207	circ_class_biomass(i,j,l,iroot,icarbon) * circ_class_n(i,j,l)
1208	biomass(i,j,icarbres,icarbon) = biomass(i,j,icarbres,icarbon) + &
1209	circ_class_biomass(i,j,l,icarbres,icarbon) * circ_class_n(i,j,l)
1210	biomass(i,j,ilabile,icarbon) = biomass(i,j,ilabile,icarbon) + &
1211	circ_class_biomass(i,j,l,ilabile,icarbon) * circ_class_n(i,j,l)
1212
1213	ENDDO
1214
1215	! set reset leaf age distribution (::age_reset) flag. Default is TRUE.
1216	! done later for better vectorization)
1217	age_reset(i) = .TRUE.
1218
1219	! Leaf mass is high enough
1220	ELSE
1221
1222	! Nothing should be done because the minimum amount of leaves
1223	! to start the growing season is present.
1224
1225	ENDIF ! leaf mass is very low
1226
1227	!! Reset when_growthinit counter: start of the growing season
1228	when_growthinit(i,j) = zero
1229	senescence(i,j) = .FALSE.
1230
1231	ELSEIF ( .NOT. natural(j) ) THEN
1232
1233	! Crops get planted the day that begin_leaves is true. If not, too
1234	! much of there reserves are used before the start of the growing
1235	! season.
1236	CALL crop_planting(npts, dt, i, j, &
1237	veget_max, PFTpresent, c0_alloc, when_growthinit, &
1238	time_hum_min, everywhere, senescence, ind, &
1239	circ_class_n, KF, leaf_frac, age, &
1240	npp_longterm, lm_lastyearmax, biomass, circ_class_biomass, &
1241	co2_to_bm, k_latosa_adapt)
1242
1243	ELSE
1244
1245	CALL ipslerr_p (3,'stomate_phenology', 'Begin_leaves is true but case is undefined','','')
1246
1247	END IF
1248
1249	!---TEMP---
1250	IF (j.EQ.test_pft .AND. ld_pheno) THEN
1251	WRITE(numout,*) 'circ_class_biomass, ',SUM(circ_class_biomass(i,test_pft,:,:,icarbon))
1252	ENDIF
1253	!----------
1254
1255	ENDIF ! start of the growing season
1256
1257	ENDDO ! loop over grid point
1258
1259	!! 5.4 reset leaf age distribution where necessary (i.e. when age_reset is TRUE)
1260	! simply say that everything is in the youngest age class
1261	! Set the youngest age class fraction to 1 and all other leaf age class fractions to 0.
1262	WHERE ( age_reset(:) )
1263
1264	leaf_frac(:,j,1) = un
1265
1266	ENDWHERE
1267
1268	DO m = 2, nleafages
1269
1270	WHERE ( age_reset(:) )
1271
1272	leaf_frac(:,j,m) = zero
1273
1274	ENDWHERE
1275
1276	ENDDO ! nleafages
1277
1278	! Ages - set all leaf ages to 0.
1279	DO m = 1, nleafages
1280
1281	WHERE ( age_reset(:) )
1282
1283	leaf_age(:,j,m) = zero
1284
1285	ENDWHERE
1286
1287	ENDDO ! nleafages
1288
1289	ENDDO ! loop over # PFTs
1290
1291	!! 6. Check mass balance closure
1292
1293	!! 6.1 Calculate final biomass
1294	pool_end(:,:,:) = zero
1295	DO ipar = 1,nparts
1296	DO iele = 1,nelements
1297	pool_end(:,:,iele) = pool_end(:,:,iele) + &
1298	(biomass(:,:,ipar,iele) * veget_max(:,:))
1299	ENDDO
1300	ENDDO
1301
1302	!! 6.2 Calculate mass balance
1303	check_intern(:,:,iatm2land,icarbon) = check_intern(:,:,iatm2land,icarbon) + &
1304	co2_to_bm(:,:) * veget_max(:,:) * dt
1305	check_intern(:,:,iland2atm,icarbon) = -un * zero
1306	check_intern(:,:,ilat2out,icarbon) = zero
1307	check_intern(:,:,ilat2in,icarbon) = -un * zero
1308	check_intern(:,:,ipoolchange,icarbon) = -un * (pool_end(:,:,icarbon) - pool_start(:,:,icarbon))
1309	closure_intern = zero
1310	DO imbc = 1,nmbcomp
1311	closure_intern(:,:,icarbon) = closure_intern(:,:,icarbon) + check_intern(:,:,imbc,icarbon)
1312	ENDDO
1313
1314	! Write the verdict
1315	DO ipts=1,npts
1316	DO ivm=1,nvm
1317	IF(ABS(closure_intern(ipts,ivm,icarbon)) .LE. min_stomate)THEN
1318	IF (ld_massbal) WRITE(numout,*) 'Mass balance closure in phenology_prognostic'
1319	ELSE
1320	WRITE(numout,*) 'Error: mass balance is not closed in phenology_prognostic'
1321	WRITE(numout,*) ' ipts,ivm; ', ipts,ivm
1322	WRITE(numout,*) ' Difference is, ', closure_intern(ipts,ivm,icarbon)
1323	WRITE(numout,*) ' pool_end,pool_start: ', pool_end(ipts,ivm,icarbon), pool_start(ipts,ivm,icarbon)
1324	WRITE(numout,*) ' check_intern,co2_to_bm,veget_max: ', &
1325	check_intern(ipts,ivm,iatm2land,icarbon),co2_to_bm(ipts,ivm), veget_max(ipts,ivm)
1326	IF(ld_stop)THEN
1327	CALL ipslerr_p (3,'phenology_prognostic', 'Mass balance error.','','')
1328	ENDIF
1329	ENDIF
1330	ENDDO
1331	ENDDO
1332
1333	IF (bavard.GE.2) WRITE(numout,*) 'Leaving phenology'
1334
1335	END SUBROUTINE phenology_prognostic
1336
1337
1338	!! ================================================================================================================================
1339	!! SUBROUTINE : pheno_hum
1340	!!
1341	!>\BRIEF The 'hum' onset model initiate leaf onset based exclusively on moisture
1342	!! availability criteria.
1343	!! Currently no PFTs are assigned to this onset model.
1344	!!
1345	!! DESCRIPTION : This model is for tropical biomes, where temperatures are high but moisture
1346	!! might be a limiting factor on growth. It is based on leaf onset model 4a in
1347	!! Botta et al. (2000), which adopts the approach of Le Roux (1995). \n
1348	!! Leaf onset occurs if the monthly moisture availability is still quite
1349	!! low (i.e. lower than the weekly availability), but the weekly availability is
1350	!! higher than the critical threshold ::availability_crit (as it reacts faster),
1351	!! which indicates the weekly moisture availability is increasing.
1352	!! OR if the monthly moisture availability is high enough (i.e. above the
1353	!! threshold value ::moiavail_always), leaf onset is initiated if this has not
1354	!! already happened. This allows vegetation in arid areas to respond to rapidly
1355	!! changing soil moisture conditions (Krinner et al., 2005). \n
1356	!! The critical weekly moisture availability threshold (::availability_crit), is
1357	!! calculated in this subroutine, and is a function of last year's maximum and
1358	!! minimum moisture availability and the PFT-dependent parameter
1359	!! ::hum_frac, which specifies how much of last year's available
1360	!! moisture is required for leaf onset, as per the equation:
1361	!! \latexonly
1362	!! \input{phenology_moi_availcrit_eqn3.tex}
1363	!! \endlatexonly
1364	!! \n
1365	!! ::hum_frac is set for each PFT in ::stomate_data from a table
1366	!! which contains all the PFT values (::hum_frac_tab) in ::stomate_constants. \n
1367	!! Last year's maximum and minimum moisture availability and the monthly and
1368	!! weekly moisture availability are
1369	!! The ::pheno_hum subroutine is called in the subroutine ::phenology.
1370	!!
1371	!! RECENT CHANGE(S): None
1372	!!
1373	!! MAIN OUTPUT VARIABLE(S): ::begin_leaves - specifies whether leaf growth can start.
1374	!!
1375	!! REFERENCE(S) :
1376	!! - Botta, A., N. Viovy, P. Ciais, P. Friedlingstein and P. Monfray (2000),
1377	!! A global prognostic scheme of leaf onset using satellite data,
1378	!! Global Change Biology, 207, 337-347.
1379	!! - Le Roux, X. (1995), Etude et modelisation des echanges d'eau et d'energie
1380	!! sol-vegetation-atmosphere dans une savane humide, PhD Thesis, University
1381	!! Pierre et Marie Curie, Paris, France.
1382	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
1383	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
1384	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
1385	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
1386	!!
1387	!! FLOWCHART :
1388	!! \latexonly
1389	!! \includegraphics[scale = 1]{pheno_hum.png}
1390	!! \endlatexonly
1391	!! \n
1392	!_ ================================================================================================================================
1393
1394	SUBROUTINE pheno_hum (npts, j, PFTpresent, allow_initpheno, &
1395	moiavail_month, moiavail_week, &
1396	maxmoiavail_lastyear, minmoiavail_lastyear, &
1397	begin_leaves)
1398
1399	!
1400	!! 0. Variable and parameter declarations
1401	!
1402
1403	!
1404	!! 0.1 Input variables
1405	!
1406	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of
1407	!! grid cells (unitless)
1408	INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
1409	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
1410	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to
1411	!! declare the beginning of
1412	!! the growing season?
1413	!! (true/false)
1414	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_month !! "monthly" moisture
1415	!! availability (0-1, unitless)
1416	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_week !! "weekly" moisture
1417	!! availability (0-1, unitless)
1418	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: maxmoiavail_lastyear !! last year's maximum
1419	!! moisture availability
1420	!! (0-1, unitless)
1421	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: minmoiavail_lastyear !! last year's minimum
1422	!! moisture availability
1423	!! (0-1, unitless)
1424
1425	!
1426	!! 0.2 Output variables
1427	!
1428
1429	!
1430	!! 0.3 Modified variables
1431	!
1432	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting
1433	!! leaves on (true/false)
1434
1435	!
1436	!! 0.4 Local variables
1437	!
1438	REAL(r_std) :: moiavail_always !! critical monthly
1439	!! moisture availability - set
1440	!! for tree or grass
1441	!! (0-1, unitless)
1442	REAL(r_std), DIMENSION(npts) :: availability_crit !! critical weekly moisture
1443	!! availability (0-1, unitless)
1444	INTEGER(i_std) :: i !! index (unitless)
1445
1446	!_ ================================================================================================================================
1447
1448	IF (bavard.GE.2) WRITE(numout,*) 'Entering hum'
1449
1450	!
1451	!! 1. Initializations
1452	!
1453
1454	!
1455	!! 1.1 first call - outputs the name of onset model and the moisture availability
1456	!! parameters for tree and grass
1457	!
1458
1459	IF ( firstcall_hum ) THEN
1460
1461	WRITE(numout,*) 'pheno_hum:'
1462	WRITE(numout,*) ' > moisture availability above which moisture tendency doesn''t matter: '
1463	WRITE(numout,*) ' trees (::moiavail_always_tree): ', moiavail_always_tree
1464	WRITE(numout,*) ' grasses (::moiavail_always_grass):', moiavail_always_grass
1465
1466	firstcall_hum = .FALSE.
1467
1468	ENDIF
1469
1470	!
1471	!! 1.2 initialize output
1472	!
1473
1474	begin_leaves(:,j) = .FALSE.
1475
1476	!
1477	!! 1.3 check the critical value ::hum_frac is defined. If not, stop.
1478	!
1479
1480	IF ( hum_frac(j) .EQ. undef ) THEN
1481
1482	WRITE(numout,*) 'hum: hum_frac is undefined for PFT (::j)',j
1483	WRITE(numout,*) 'We stop.'
1484	STOP
1485
1486	ENDIF
1487
1488	!
1489	!! 1.4 set the critical monthly moisture availability above which we always detect the beginning of the
1490	!! growing season - set as the moisture availability for trees or grass.
1491	!
1492
1493	IF ( is_tree(j) ) THEN
1494	moiavail_always = moiavail_always_tree
1495	ELSE
1496	moiavail_always = moiavail_always_grass
1497	ENDIF
1498
1499	!
1500	!! 2. Check if biometeorological conditions are favourable for leaf growth.
1501	!! The PFT has to be there and start of growing season must be allowed
1502	!
1503
1504	DO i = 1, npts
1505
1506	IF ( PFTpresent(i,j) .AND. allow_initpheno(i,j) ) THEN
1507
1508	!! 2.1 Calculate the critical weekly moisture availability: depends linearly on the last year
1509	!! minimum and maximum moisture availabilities, and on the parameter ::hum_frac.
1510
1511	availability_crit(i) = minmoiavail_lastyear(i,j) + hum_frac(j) * &
1512	( maxmoiavail_lastyear(i,j) - minmoiavail_lastyear(i,j) )
1513
1514	!! 2.2 Determine if growing season should start (if so, ::begin_leaves set to TRUE).
1515	!! Leaf onset occurs if the monthly moisture availability is still quite
1516	!! low (i.e. lower than the weekly availability), but the weekly availability is
1517	!! already higher than the critical threshold ::availability_crit (as it reacts faster),
1518	!! which indicates the weekly moisture availability is increasing.
1519	!! OR if the monthly moisture availability is high enough (i.e. above the threshold value
1520	!! ::moiavail_always), leaf onset is initiated if this has not already happened.
1521
1522	IF ( ( ( moiavail_week(i,j) .GE. availability_crit(i) ) .AND. &
1523	( moiavail_month(i,j) .LT. moiavail_week(i,j) ) ) .OR. &
1524	( moiavail_month(i,j) .GE. moiavail_always ) ) THEN
1525	begin_leaves(i,j) = .TRUE.
1526	ENDIF
1527
1528	ENDIF ! PFT there and start of growing season allowed
1529
1530	ENDDO ! end loop over grid points
1531
1532	IF (bavard.GE.2) WRITE(numout,*) 'Leaving hum'
1533
1534	END SUBROUTINE pheno_hum
1535
1536
1537	!! ================================================================================================================================
1538	!! SUBROUTINE : pheno_moi
1539	!!
1540	!>\BRIEF The 'moi' onset model (::pheno_moi) initiates leaf onset based exclusively
1541	!! on moisture availability criteria.
1542	!! It is very similar to the 'hum' onset model but instead of the weekly moisture
1543	!! availability being higher than a constant threshold, the condition is that the
1544	!! moisture minimum happened a sufficiently long time ago.
1545	!! Currently PFT 3 (Tropical Broad-leaved Raingreen) is assigned to this model.
1546	!!
1547	!! DESCRIPTION : This model is for tropical biomes, where temperatures are high but moisture
1548	!! might be a limiting factor on growth. It is based on leaf onset model 4b in
1549	!! Botta et al. (2000).
1550	!! Leaf onset begins if the plant moisture availability minimum was a sufficiently
1551	!! time ago, as specified by the PFT-dependent parameter ::hum_min_time
1552	!! AND if the "monthly" moisture availability is lower than the "weekly"
1553	!! availability (indicating that soil moisture is increasing).
1554	!! OR if the monthly moisture availability is high enough (i.e. above the threshold
1555	!! value ::moiavail_always), leaf onset is initiated if this has not already
1556	!! happened. \n
1557	!! ::hum_min_time is set for each PFT in ::stomate_data, and is
1558	!! defined in the table ::hum_min_time_tab in ::stomate_constants. \n
1559	!! ::moiavail_always is defined for both tree and grass in this subroutine
1560	!! (set to 1. and 0.6 respectively). \n
1561	!! The ::pheno_moi subroutine is called in the subroutine ::phenology.
1562	!!
1563	!! RECENT CHANGE(S): None
1564	!!
1565	!! MAIN OUTPUT VARIABLE(S): ::begin_leaves - specifies whether leaf growth can start.
1566	!!
1567	!! REFERENCE(S) :
1568	!! - Botta, A., N. Viovy, P. Ciais, P. Friedlingstein and P. Monfray (2000),
1569	!! A global prognostic scheme of leaf onset using satellite data,
1570	!! Global Change Biology, 207, 337-347.
1571	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
1572	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
1573	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
1574	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
1575	!!
1576	!! FLOWCHART :
1577	!! \latexonly
1578	!! \includegraphics[scale = 1]{pheno_moi.png}
1579	!! \endlatexonly
1580	!! \n
1581	!_ ================================================================================================================================
1582
1583	SUBROUTINE pheno_moi (npts, j, PFTpresent, allow_initpheno, &
1584	time_hum_min, &
1585	moiavail_month, moiavail_week, &
1586	begin_leaves)
1587
1588	!
1589	!! 0. Variable and parameter declaration
1590	!
1591
1592	!
1593	!! 0.1 Input variables
1594	!
1595	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells (unitless)
1596	INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
1597	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
1598	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to declare the beginning of the
1599	!! growing season? (true/false)
1600	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: time_hum_min !! time elapsed since strongest moisture
1601	!! availability (days)
1602	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_month !! "monthly" moisture availability (0-1, unitless)
1603	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_week !! "weekly" moisture availability (0-1, unitless)
1604
1605	!
1606	!! 0.2 Output variables
1607	!
1608	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting leaves on (true/false)
1609
1610	!
1611	!! 0.3 Modified variables
1612	!
1613
1614	!
1615	!! 0.4 Local variables
1616	!
1617	REAL(r_std) :: moiavail_always !! critical moisture availability -
1618	!! set for tree or grass
1619	!! (0-1, unitless)
1620	INTEGER(i_std) :: i !! index (unitless)
1621
1622	!_ ================================================================================================================================
1623
1624	IF (bavard.GE.2) WRITE(numout,*) 'Entering moi'
1625
1626	!
1627	!! 1. Initializations
1628	!
1629
1630	!
1631	!! 1.1 first call - outputs the name of onset model and the moisture availability
1632	!! parameters for tree and grass
1633	!
1634
1635	IF ( firstcall_moi ) THEN
1636
1637	WRITE(numout,*) 'pheno_moi:'
1638	WRITE(numout,*) ' > moisture availability above which moisture tendency doesn''t matter: '
1639	WRITE(numout,*) ' trees (::moiavail_always_tree):', moiavail_always_tree
1640	WRITE(numout,*) ' grasses (::moiavail_always_grass):', moiavail_always_grass
1641
1642	firstcall_moi = .FALSE.
1643
1644	ENDIF
1645
1646	!
1647	!! 1.2 initialize output
1648	!
1649
1650	begin_leaves(:,j) = .FALSE.
1651
1652	!
1653	!! 1.3 check the critical value ::hum_min_time is definded. If not, stop
1654	!
1655
1656	IF ( hum_min_time(j) .EQ. undef ) THEN
1657
1658	WRITE(numout,*) 'moi: hum_min_time is undefined for PFT (::j) ',j
1659	WRITE(numout,*) 'We stop.'
1660	STOP
1661
1662	ENDIF
1663
1664	!
1665	!! 1.4 set the critical monthly moisture availability above which we always detect the beginning of the
1666	!! growing season - set as the moisture availability for trees or grass.
1667	!
1668
1669	IF ( is_tree(j) ) THEN
1670	moiavail_always = moiavail_always_tree
1671	ELSE
1672	moiavail_always = moiavail_always_grass
1673	ENDIF
1674
1675	!
1676	!! 2. Check if biometeorological conditions are favourable for leaf growth.
1677	!! The PFT has to be there and start of growing season must be allowed.
1678	!
1679
1680	DO i = 1, npts
1681
1682	IF ( PFTpresent(i,j) .AND. allow_initpheno(i,j) ) THEN
1683
1684	!! 2.1 Determine if growing season should start (if so, ::begin_leaves set to TRUE).
1685	!! The favorable season starts if the moisture minimum (::time_hum_min) was a sufficiently long
1686	!! time ago, i.e. greater than the threshold specified by the parameter ::hum_min_time
1687	!! and if the "monthly" moisture availability is lower than the "weekly"
1688	!! availability (indicating that soil moisture is increasing).
1689	!! OR if the monthly moisture availability is high enough (i.e. above the threshold value
1690	!! ::moiavail_always), initiate the growing season if this has not happened yet.
1691
1692	IF ( ( ( moiavail_week(i,j) .GT. moiavail_month(i,j) ) .AND. &
1693	!!$ ( biomass(i,j,ileaf,icarbon) .LT. min_stomate) .AND. &
1694	( time_hum_min(i,j) .GT. hum_min_time(j) ) ) .OR. &
1695	( moiavail_month(i,j) .GE. moiavail_always ) ) THEN
1696	begin_leaves(i,j) = .TRUE.
1697	ENDIF
1698
1699	ENDIF ! PFT there and start of growing season allowed
1700
1701	ENDDO ! end loop over grid points
1702
1703	IF (bavard.GE.2) WRITE(numout,*) 'Leaving moi'
1704
1705	END SUBROUTINE pheno_moi
1706
1707
1708	!! ================================================================================================================================
1709	!! SUBROUTINE : pheno_humgdd
1710	!!
1711	!>\BRIEF The 'humgdd' onset model initiates leaf onset based on mixed conditions of
1712	!! temperature and moisture availability criteria.
1713	!! Currently no PFTs are assigned to this onset model.
1714	!!
1715	!! DESCRIPTION : In this model the Growing Degree Day (GDD) model (Chuine, 2000) is combined
1716	!! with the 'hum' onset model (::pheno_hum), which has previously been described,
1717	!! in order to account for dependence on both temperature and moisture conditions
1718	!! in warmer climates. \n.
1719	!! The GDD model specifies that daily temperatures above a threshold of -5
1720	!! degrees C are summed, minus this threshold, giving the GDD, starting from
1721	!! the beginning of the dormancy period (::time_lowgpp>0), i.e. since the leaves
1722	!! were lost. \n.
1723	!! The dormancy time-length is represented by the variable
1724	!! ::time_lowgpp, which is calculated in ::stomate_season. It is increased by
1725	!! the stomate time step when the weekly GPP is lower than a threshold. Otherwise
1726	!! it is set to zero. \n
1727	!! Leaf onset begins when the a PFT-dependent GDD-threshold is reached.
1728	!! In addition there are temperature and moisture conditions.
1729	!! The temperature condition specifies that the monthly temperature has to be
1730	!! higher than a constant threshold (::t_always) OR
1731	!! the weekly temperature is higher than the monthly temperature.
1732	!! There has to be at least some moisture. The moisture condition
1733	!! is exactly the same as the 'hum' onset model (::pheno_hum), which has already
1734	!! been described. \n
1735	!! The GDD (::gdd_m5_dormance) is calculated in ::stomate_season. GDD is set to
1736	!! undef if beginning of the growing season detected, i.e. when there is GPP
1737	!! (::time_lowgpp>0).
1738	!! The parameter ::t_always is defined as 10 degrees C in this subroutine,
1739	!! as are the parameters ::moisture_avail_tree and ::moisture_avail_grass
1740	!! (set to 1 and 0.6 respectively), which are used in the moisture condition
1741	!! (see ::pheno_moi onset model description). \n
1742	!! The PFT-dependent GDD threshold (::gdd_crit) is calculated as in the onset
1743	!! model ::pheno_humgdd, using the equation:
1744	!! \latexonly
1745	!! \input{phenology_hummoigdd_gddcrit_eqn4.tex}
1746	!! \endlatexonly
1747	!! \n
1748	!! The three GDDcrit parameters (::gdd(j,*)) are set for each PFT in
1749	!! ::stomate_data, and three tables defining each of the three critical GDD
1750	!! parameters for each PFT is given in ::gdd_crit1_tab, ::gdd_crit2_tab and
1751	!! ::gdd_crit3_tab in ::stomate_constants. \n
1752	!! The ::pheno_humgdd subroutine is called in the subroutine ::phenology.
1753	!!
1754	!! RECENT CHANGES: None
1755	!!
1756	!! MAIN OUTPUT VARIABLES: ::begin_leaves - specifies whether leaf growth can start
1757	!!
1758	!! REFERENCE(S) :
1759	!! - Botta, A., N. Viovy, P. Ciais, P. Friedlingstein and P. Monfray (2000),
1760	!! A global prognostic scheme of leaf onset using satellite data,
1761	!! Global Change Biology, 207, 337-347.
1762	!! - Chuine, I (2000), A unified model for the budburst of trees, Journal of
1763	!! Theoretical Biology, 207, 337-347.
1764	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
1765	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
1766	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
1767	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
1768	!!
1769	!! FLOWCHART :
1770	!! \latexonly
1771	!! \includegraphics[scale = 1]{pheno_humgdd.png}
1772	!! \endlatexonly
1773	!! \n
1774	!_ ================================================================================================================================
1775
1776	SUBROUTINE pheno_humgdd (npts, j, PFTpresent, allow_initpheno, gdd, &
1777	maxmoiavail_lastyear, minmoiavail_lastyear, &
1778	tlong_ref, t2m_month, t2m_week, &
1779	moiavail_week, moiavail_month, &
1780	begin_leaves)
1781
1782	!
1783	!! 0. Variable and parameter declaration
1784	!
1785
1786	!
1787	!! 0.1 Input variables
1788	!
1789	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells
1790	!! (unitless)
1791	INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
1792	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
1793	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to declare the beginning
1794	!! of the growing season? (true/false)
1795	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: gdd !! growing degree days, calculated since
1796	!! leaves have fallen (C)
1797	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: maxmoiavail_lastyear !! last year's maximum moisture
1798	!! availability (0-1, unitless)
1799	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: minmoiavail_lastyear !! last year's minimum moisture
1800	!! availability (0-1, unitless)
1801	REAL(r_std), DIMENSION(npts), INTENT(in) :: tlong_ref !! "long term" 2 meter temperatures (K)
1802	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures (K)
1803	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
1804	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_week !! "weekly" moisture availability
1805	!! (0-1, unitless)
1806	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_month !! "monthly" moisture availability
1807	!! (0-1, unitless)
1808
1809	!
1810	!! 0.2 Output variables
1811	!
1812	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting leaves on
1813	!! (true/false)
1814
1815	!
1816	!! 0.3 Modified variables
1817	!
1818
1819	!
1820	!! 0.4 Local variables
1821	!
1822	REAL(r_std) :: moiavail_always !! critical moisture availability -
1823	!! set for tree or grass
1824	!! (0-1, unitless)
1825	REAL(r_std), DIMENSION(npts) :: moiavail_crit !! critical moisture availability
1826	!! (0-1, unitless)
1827	REAL(r_std), DIMENSION(npts) :: tl !! long term temperature (C)
1828	REAL(r_std), DIMENSION(npts) :: gdd_crit !! critical GDD (C)
1829	INTEGER(i_std) :: i !! index (unitless)
1830
1831	!_ ================================================================================================================================
1832
1833	IF (bavard.GE.2) WRITE(numout,*) 'Entering humgdd'
1834
1835	!
1836	!! 1. Initializations
1837	!
1838
1839	!
1840	!! 1.1 first call - outputs the name of the onset model, the values of the
1841	!! moisture availability parameters for tree and grass, and the value of the
1842	!! critical monthly temperature.
1843	!
1844
1845	IF ( firstcall_humgdd ) THEN
1846
1847	WRITE(numout,*) 'pheno_humgdd:'
1848	WRITE(numout,*) ' > moisture availability above which moisture tendency doesn''t matter: '
1849	WRITE(numout,*) ' trees (::moiavail_always_tree): ', moiavail_always_tree
1850	WRITE(numout,*) ' grasses (::moiavail_always_grass): ', moiavail_always_grass
1851	WRITE(numout,*) ' > monthly temp. above which temp. tendency doesn''t matter: ', &
1852	t_always
1853
1854	firstcall_humgdd = .FALSE.
1855
1856	ENDIF
1857
1858	!
1859	!! 1.2 initialize output
1860	!
1861
1862	begin_leaves(:,j) = .FALSE.
1863
1864	!
1865	!! 1.3 check the critical values ::gdd and ::pheno_crit_hum_frac are defined.
1866	!! If not, stop.
1867	!
1868
1869	IF ( ANY(pheno_gdd_crit(j,:) .EQ. undef) ) THEN
1870
1871	WRITE(numout,*) 'humgdd: pheno_gdd_crit is undefined for PFT (::j) ',j
1872	WRITE(numout,*) 'We stop.'
1873	STOP
1874
1875	ENDIF
1876
1877	IF ( hum_frac(j) .EQ. undef ) THEN
1878
1879	WRITE(numout,*) 'humgdd: hum_frac is undefined for PFT (::j) ',j
1880	WRITE(numout,*) 'We stop.'
1881	STOP
1882
1883	ENDIF
1884
1885	!
1886	!! 1.4 set the critical moisture availability above which we always detect the beginning of the
1887	!! growing season - set as the moisture availability for trees or grass.
1888	!
1889
1890	IF ( is_tree(j) ) THEN
1891	moiavail_always = moiavail_always_tree
1892	ELSE
1893	moiavail_always = moiavail_always_grass
1894	ENDIF
1895
1896	!
1897	!! 2. Check if biometeorological conditions are favourable for leaf growth.
1898	!! The PFT has to be there, start of growing season must be allowed,
1899	!! and GDD has to be defined.
1900	!
1901
1902	DO i = 1, npts
1903
1904	IF ( PFTpresent(i,j) .AND. allow_initpheno(i,j) .AND. &
1905	( gdd(i,j) .NE. undef ) ) THEN
1906
1907	!! 2.1 Calculate the critical weekly moisture availability: depends linearly on the last year
1908	!! minimum and maximum moisture availabilities, and on the parameter ::hum_frac.,
1909	!! (as in the ::pheno_hum model), as per the equation:
1910
1911	moiavail_crit(i) = minmoiavail_lastyear(i,j) + hum_frac(j) * &
1912	( maxmoiavail_lastyear(i,j) - minmoiavail_lastyear(i,j) )
1913
1914	!! 2.2 Calculate the critical GDD (::gdd_crit), which is a function of the PFT-dependent
1915	!! critical GDD and the "long term" 2 meter air temperatures.
1916
1917	tl(i) = tlong_ref(i) - ZeroCelsius
1918	gdd_crit(i) = pheno_gdd_crit(j,1) + tl(i)*pheno_gdd_crit(j,2) + &
1919	tl(i)tl(i)pheno_gdd_crit(j,3)
1920
1921	!!Added to make use of the optimised phenology parameters of Natasha MacBean
1922	!!In ORCHIS, a multiplicative factor (opti_kpheno_crit) was optimised in stead of ggd_crit
1923	gdd_crit(i) = gdd_crit(i)*opti_kpheno_crit(j)
1924
1925	!! 2.3 Determine if the growing season should start (if so, ::begin_leaves set to TRUE).
1926	!! - Has the critical gdd been reached and is the temperature increasing?
1927	!! - Is there at least some humidity/moisture availability?
1928	!! This occurs if the critical gdd (::gdd_crit) has been reached
1929	!! AND that is temperature increasing, which is true either if the monthly
1930	!! temperature being higher than the threshold ::t_always, OR if the weekly
1931	!! temperature is higher than the monthly,
1932	!! AND finally that there is sufficient moisture availability, which is
1933	!! the same condition as for the ::pheno_hum onset model.
1934
1935	IF ( ( gdd(i,j) .GE. gdd_crit(i) ) .AND. &
1936	( ( t2m_week(i) .GT. t2m_month(i) ) .OR. &
1937	( t2m_month(i) .GT. t_always ) ) .AND. &
1938	( ( ( moiavail_week(i,j) .GE. moiavail_crit(i) ) .AND. &
1939	( moiavail_month(i,j) .LT. moiavail_crit(i) ) ) .OR. &
1940	( moiavail_month(i,j) .GE. moiavail_always ) ) ) THEN
1941	begin_leaves(i,j) = .TRUE.
1942	ENDIF
1943
1944	ENDIF ! PFT there and start of growing season allowed
1945
1946	ENDDO ! End loop over grid points
1947
1948	IF (bavard.GE.2) WRITE(numout,*) 'Leaving humgdd'
1949
1950	END SUBROUTINE pheno_humgdd
1951
1952
1953	!! ================================================================================================================================
1954	!! SUBROUTINE : pheno_moigdd
1955	!!
1956	!>\BRIEF The 'moigdd' onset model initiates leaf onset based on mixed temperature
1957	!! and moisture availability criteria.
1958	!! Currently PFTs 10 - 13 (C3 and C4 grass, and C3 and C4 agriculture)
1959	!! are assigned to this model.
1960	!!
1961	!! DESCRIPTION : This onset model combines the GDD model (Chuine, 2000), as described for
1962	!! the 'humgdd' onset model (::pheno_humgdd), and the 'moi' model, in order
1963	!! to account for dependence on both temperature and moisture conditions in
1964	!! warmer climates. \n
1965	!! Leaf onset begins when the a PFT-dependent GDD threshold is reached.
1966	!! In addition there are temperature and moisture conditions.
1967	!! The temperature condition specifies that the monthly temperature has to be
1968	!! higher than a constant threshold (::t_always) OR
1969	!! the weekly temperature is higher than the monthly temperature.
1970	!! There has to be at least some moisture. The moisture condition
1971	!! is exactly the same as the 'moi' onset model (::pheno_moi), which has
1972	!! already been described. \n
1973	!! GDD is set to undef if beginning of the growing season detected.
1974	!! As in the ::pheno_humgdd model, the parameter ::t_always is defined as
1975	!! 10 degrees C in this subroutine, as are the parameters ::moisture_avail_tree
1976	!! and ::moisture_avail_grass (set to 1 and 0.6 respectively), which are used
1977	!! in the moisture condition (see ::pheno_moi onset model description). \n
1978	!! The PFT-dependent GDD threshold (::gdd_crit) is calculated as in the onset
1979	!! model ::pheno_humgdd, using the equation:
1980	!! \latexonly
1981	!! \input{phenology_hummoigdd_gddcrit_eqn4.tex}
1982	!! \endlatexonly
1983	!! \n
1984	!! where i and j are the grid cell and PFT respectively.
1985	!! The three GDDcrit parameters (::gdd(j,*)) are set for each PFT in
1986	!! ::stomate_data, and three tables defining each of the three critical GDD
1987	!! parameters for each PFT is given in ::gdd_crit1_tab, ::gdd_crit2_tab and
1988	!! ::gdd_crit3_tab in ::stomate_constants. \n
1989	!! The ::pheno_moigdd subroutine is called in the subroutine ::phenology.
1990	!!
1991	!! RECENT CHANGE(S): None
1992	!!
1993	!! MAIN OUTPUT VARIABLE(S): ::begin_leaves - specifies whether leaf growth can start
1994	!!
1995	!! REFERENCE(S) :
1996	!! - Botta, A., N. Viovy, P. Ciais, P. Friedlingstein and P. Monfray (2000),
1997	!! A global prognostic scheme of leaf onset using satellite data,
1998	!! Global Change Biology, 207, 337-347.
1999	!! - Chuine, I (2000), A unified model for the budburst of trees, Journal of
2000	!! Theoretical Biology, 207, 337-347.
2001	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
2002	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
2003	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
2004	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
2005	!!
2006	!! FLOWCHART :
2007	!! \latexonly
2008	!! \includegraphics[scale = 1]{pheno_moigdd.png}
2009	!! \endlatexonly
2010	!! \n
2011	!_ ================================================================================================================================
2012
2013	SUBROUTINE pheno_moigdd (npts, j, PFTpresent, allow_initpheno, gdd, &
2014	time_hum_min, &
2015	tlong_ref, t2m_month, t2m_week, &
2016	moiavail_week, moiavail_month, &
2017	begin_leaves)
2018
2019	!
2020	!! 0. Variable and parameter declaration
2021	!
2022
2023	!
2024	!! 0.1 Input variables
2025	!
2026	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells (unitless)
2027	INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
2028	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
2029	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to decalre the beginning of the
2030	!! growing season? (true/false)
2031	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: gdd !! growing degree days, calculated since leaves
2032	!! have fallen (C)
2033	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: time_hum_min !! time elapsed since strongest moisture
2034	!! availability (days)
2035	REAL(r_std), DIMENSION(npts), INTENT(in) :: tlong_ref !! "long term" 2 meter temperatures (K)
2036	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures (K)
2037	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
2038	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_week !! "weekly" moisture availability (0-1, unitless)
2039	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_month !! "monthly" moisture availability (0-1, unitless)
2040
2041	!
2042	!! 0.2 Output variables
2043	!
2044
2045	!
2046	!! 0.3 Modified variables
2047	!
2048	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting leaves on (true/false)
2049
2050	!
2051	!! 0.4 Local variables
2052	!
2053	REAL(r_std) :: moiavail_always !! critical moisture availability -
2054	!! set for tree or grass
2055	!! (0-1, unitless)
2056	REAL(r_std), DIMENSION(npts) :: tl !! long term temperature (C)
2057	REAL(r_std), DIMENSION(npts) :: gdd_crit !! critical GDD (C)
2058	INTEGER(i_std) :: i !! index (unitless)
2059
2060	!_ ================================================================================================================================
2061
2062	IF (bavard.GE.2) WRITE(numout,*) 'Entering moigdd'
2063
2064	!! 1. Initializations
2065
2066	!! 1.1 first call
2067	!! Outputs the name of the onset model, the values of the
2068	!! moisture availability parameters for tree and grass, and the value of the
2069	!! critical monthly temperature.
2070	IF ( firstcall_moigdd ) THEN
2071
2072	WRITE(numout,*) 'pheno_moigdd:'
2073	WRITE(numout,*) ' > moisture availability above which moisture tendency doesn''t matter: '
2074	WRITE(numout,*) ' trees (::moiavail_always_tree) :', moiavail_always_tree
2075	WRITE(numout,*) ' grasses (::moiavail_always_grass) :', moiavail_always_grass
2076	WRITE(numout,*) ' > monthly temp. above which temp. tendency doesn''t matter (::t_always): ', &
2077	t_always
2078
2079	firstcall_moigdd = .FALSE.
2080
2081	ENDIF
2082
2083	!! 1.2 initialize output
2084	begin_leaves(:,j) = .FALSE.
2085
2086	!! 1.3 check the critical values ::gdd and ::pheno_crit_hum_min_time are defined.
2087	!! If not, stop.
2088	IF ( ANY(pheno_gdd_crit(j,:) .EQ. undef) ) THEN
2089
2090	WRITE(numout,*) 'moigdd: pheno_gdd_crit is undefined for PFT',j
2091	WRITE(numout,*) 'We stop.'
2092	STOP
2093
2094	ENDIF
2095
2096	IF ( hum_min_time(j) .EQ. undef ) THEN
2097
2098	WRITE(numout,*) 'moigdd: hum_min_time is undefined for PFT',j
2099	WRITE(numout,*) 'We stop.'
2100	STOP
2101
2102	ENDIF
2103
2104	!! 1.4 set the critical moisture availability above which we always detect the beginning of the
2105	!! growing season - set as the moisture availability for trees or grass
2106	IF ( is_tree(j) ) THEN
2107
2108	moiavail_always = moiavail_always_tree
2109
2110	ELSE
2111
2112	moiavail_always = moiavail_always_grass
2113
2114	ENDIF
2115
2116	!+++TEMP+++
2117	IF (ld_pheno .AND. j.EQ.test_pft) THEN
2118	WRITE(numout,*) 'phenology PFTpresent, ', j, PFTpresent(:,j)
2119	WRITE(numout,*) 'phenology allow_initpheno, ', allow_initpheno(:,j)
2120	WRITE(numout,*) 'phenology gdd, ', gdd(:,j)
2121	ENDIF
2122	!++++++++++
2123
2124
2125	!! 2. Check if biometeorological conditions are favourable for leaf growth.
2126	!! The PFT has to be there, the start of growing season must be allowed,
2127	!! and GDD has to be defined.
2128	DO i = 1, npts
2129
2130	IF ( (PFTpresent(i,j) .AND. allow_initpheno(i,j)) .AND. &
2131	(gdd(i,j) .NE. undef) ) THEN
2132
2133	!! 2.1 Calculate the critical GDD (::gdd_crit), which is a function of the PFT-dependent
2134	!! critical GDD and the "long term" 2 meter air temperatures
2135	tl(i) = tlong_ref(i) - ZeroCelsius
2136	gdd_crit(i) = pheno_gdd_crit(j,1) + tl(i)*pheno_gdd_crit(j,2) + &
2137	tl(i)tl(i)pheno_gdd_crit(j,3)
2138
2139	!!Added to make use of the optimised phenology parameters of Natasha MacBean
2140	!!In ORCHIS, a multiplicative factor (opti_kpheno_crit) was optimised in stead of ggd_crit
2141	gdd_crit(i) = gdd_crit(i)*opti_kpheno_crit(j)
2142
2143	!+++TEMP+++
2144	IF (ld_pheno .AND. j.EQ.test_pft) THEN
2145	WRITE(numout,*) 'phenology gdd_crit, ', gdd_crit(:)
2146	WRITE(numout,*) 'phenology t2m_week, ', t2m_week(:)
2147	WRITE(numout,*) 'phenology t2m_month, ', t2m_month(:)
2148	WRITE(numout,*) 'phenology t_always, ', t_always
2149	WRITE(numout,*) 'phenology time_hum_min, ', time_hum_min(:,j)
2150	WRITE(numout,*) 'phenology hum_min_time, ', hum_min_time(j)
2151	WRITE(numout,*) 'phenology moiavail_week, ', moiavail_week(:,j)
2152	WRITE(numout,*) 'phenology moiavail_month, ', moiavail_month(:,j)
2153	WRITE(numout,*) 'phenology moiavail_always, ', moiavail_always
2154	WRITE(numout,*) 'phenology, is_temperate, ', is_temperate(j)
2155	WRITE(numout,*) 'phenology, is_boreal, ', is_boreal(j)
2156	ENDIF
2157	!++++++++++
2158
2159	!! 2.2 Determine if the growing season should start (if so, ::begin_leaves set to TRUE).
2160	!! This occurs if the critical gdd (::gdd_crit) has been reached
2161	!! AND that is temperature increasing, which is true either if the monthly
2162	!! temperature being higher than the threshold ::t_always, OR if the weekly
2163	!! temperature is higher than the monthly,
2164	!! AND finally that there is sufficient moisture availability, which is
2165	!! the same condition as for the ::pheno_moi onset model.
2166
2167	!+++CHECK+++
2168	!!$ t_always = 283
2169	IF ( ( is_temperate(j) ) .AND. &
2170	( gdd(i,j) .GE. gdd_crit(i) ) .AND. &
2171	( ( t2m_week(i) .GT. t2m_month(i) ) .OR. &
2172	( t2m_month(i) .GT. 283 ) ) .AND. &
2173	( moiavail_month(i,j) .GE. moiavail_always ) ) THEN
2174	!!$ ( time_hum_min(i,j) .GT. hum_min_time(j) .OR. &
2175	!!$ ( moiavail_week(i,j) .GT. moiavail_month(i,j) .OR. &
2176
2177	begin_leaves(i,j) = .TRUE.
2178
2179
2180	!+++TEMP+++
2181	IF (ld_pheno .AND. j.EQ.test_pft) THEN
2182	WRITE(numout,*) 'phenology updated begin_leaves, ', begin_leaves(:,j)
2183	ENDIF
2184	!++++++++++
2185
2186	ENDIF
2187	!++++++++++++
2188
2189	!+++CHECK+++
2190	!!$ t_always = 279
2191	IF ( ( is_boreal(j) ) .AND. &
2192	( gdd(i,j) .GE. gdd_crit(i) ) .AND. &
2193	( ( t2m_week(i) .GT. t2m_month(i) ) .OR. &
2194	( t2m_month(i) .GT. 279 ) ) .AND. &
2195	( moiavail_month(i,j) .GE. moiavail_always ) ) THEN
2196	!!$ ( ( time_hum_min(i,j) .GT. hum_min_time(j) ) .OR. &
2197
2198	begin_leaves(i,j) = .TRUE.
2199
2200
2201	!+++TEMP+++
2202	IF (ld_pheno .AND. j.EQ.test_pft) THEN
2203	WRITE(numout,*) 'phenology updated begin_leaves, ', begin_leaves(:,j)
2204	ENDIF
2205	!++++++++++
2206
2207	ENDIF
2208	!++++++++++++
2209
2210	ENDIF ! PFT there and start of growing season allowed
2211
2212	ENDDO
2213
2214	IF (bavard.GE.2) WRITE(numout,*) 'Leaving moigdd'
2215
2216	END SUBROUTINE pheno_moigdd
2217
2218
2219	!! ================================================================================================================================
2220	!! SUBROUTINE : pheno_ncdgdd
2221	!!
2222	!>\BRIEF The Number of Chilling Days - Growing Degree Day (NCD-GDD) model initiates
2223	!! leaf onset if a certain relationship between the number of chilling days (NCD)
2224	!! since leaves were lost, and the growing degree days (GDD) since midwinter, is
2225	!! fulfilled.
2226	!! Currently PFT 6 (Temperate Broad-leaved Summergreen) and PFT 8 (Boreal Broad-
2227	!! leaved Summergreen) are assigned to this model.
2228	!!
2229	!! DESCRIPTION : Experiments have shown that some
2230	!! species have a "chilling" requirement, i.e. their physiology needs cold
2231	!! temperatures to trigger the mechanism that will allow the following budburst
2232	!! (e.g. Orlandi et al., 2004).
2233	!! An increase in chilling days, defined as a day with a daily mean air
2234	!! temperature below a PFT-dependent threshold, reduces a plant's GDD demand
2235	!! (Cannell and Smith, 1986; Murray et al., (1989); Botta et al., 2000).
2236	!! The GDD threshold therefore decreases as NCD
2237	!! increases, using the following empirical negative explonential law:
2238	!! \latexonly
2239	!! \input{phenology_ncdgdd_gddmin_eqn5.tex}
2240	!! \endlatexonly
2241	!! \n
2242	!! The constants used have been calibrated against data CHECK FOR REFERENCE OR PERSON WHO DID UPDATE.
2243	!! Leaf onset begins if the GDD is higher than the calculated minimum GDD
2244	!! (dependent upon NCD) AND if the weekly temperature is higher than the monthly
2245	!! temperature. This is to ensure the temperature is increasing. \n
2246	!! The dormancy time-length is represented by the variable
2247	!! ::time_lowgpp, which is calculated in ::stomate_season. It is increased by
2248	!! the stomate time step when the weekly GPP is lower than a threshold. Otherwise
2249	!! it is set to zero. \n
2250	!! The NCD (::ncd_dormance) is calculated in ::stomate_season as
2251	!! the number of days with a temperature below a PFT-dependent constant threshold
2252	!! (::ncdgdd_temp), starting from the beginning of the dormancy period
2253	!! (::time_lowgpp>0), i.e. since the leaves were lost. \n
2254	!! The growing degree day sum of the temperatures higher than
2255	!! ::ncdgdd_temp (GDD) since midwinter (::gdd_midwinter)
2256	!! is also calculated in ::stomate_season.
2257	!! Midwinter is detected if the monthly temperature is lower than the weekly
2258	!! temperature AND the monthly temperature is lower than the long-term
2259	!! temperature. ::gdd_minter is therefore set to 0 at the beginning of midwinter
2260	!! and increased with each temperature greater than the PFT-dependent threshold.
2261	!! When midsummer is detected (the opposite of the above conditions),
2262	!! ::gdd_midwinter is set to undef.
2263	!! CHECK! WHEN TO START OF DORMANCY BEEN MODIFIED FROM BOTTA- ADD IN?
2264	!! The ::pheno_ncdgdd subroutine is called in the subroutine ::phenology.
2265	!!
2266	!! RECENT CHANGE(S): None
2267	!!
2268	!! MAIN OUTPUT VARIABLE(S): ::begin_leaves - specifies whether leaf growth can start
2269	!!
2270	!! REFERENCE(S) :
2271	!! - Botta, A., N. Viovy, P. Ciais, P. Friedlingstein and P. Monfray (2000),
2272	!! A global prognostic scheme of leaf onset using satellite data,
2273	!! Global Change Biology, 207, 337-347.
2274	!! - Cannell, M.J.R. and R.I. Smith (1986), Climatic warming, spring budburst and
2275	!! frost damage on trees, Journal of Applied Ecology, 23, 177-191.
2276	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
2277	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
2278	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
2279	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
2280	!! - Murray, M.B., G.R. Cannell and R.I. Smith (1989), Date of budburst of fifteen
2281	!! tree species in Britain following climatic warming, Journal of Applied Ecology,
2282	!! 26, 693-700.
2283	!! - Orlandi, F., H. Garcia-Mozo, L.V. Ezquerra, B. Romano, E. Dominquez, C. Galan,
2284	!! and M. Fornaciari (2004), Phenological olive chilling requirements in Umbria
2285	!! (Italy) and Andalusia (Spain), Plant Biosystems, 138, 111-116.
2286	!!
2287	!! FLOWCHART :
2288	!! \latexonly
2289	!! \includegraphics[scale = 1]{pheno_ncdgdd.png}
2290	!! \endlatexonly
2291	!! \n
2292	!_ ================================================================================================================================
2293
2294	SUBROUTINE pheno_ncdgdd (npts, j, PFTpresent, allow_initpheno, &
2295	ncd_dormance, gdd_midwinter, &
2296	t2m_month, t2m_week, begin_leaves)
2297
2298	!
2299	!! 0. Variable and parameter declaration
2300	!
2301
2302	!
2303	!! 0.1 Input variables
2304	!
2305	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells (unitless)
2306	INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
2307	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
2308	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to declare the beginning of the
2309	!! growing season? (true/false)
2310	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: ncd_dormance !! number of chilling days since leaves were lost
2311	!! (days)
2312	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gdd_midwinter !! growing degree days since midwinter (C)
2313	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures (K)
2314	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
2315
2316	!
2317	!! 0.2 Output variables
2318	!
2319
2320	!
2321	!! 0.3 Modified variables
2322	!
2323	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting leaves on (true/false)
2324
2325	!
2326	!! 0.4 Local variables
2327	!
2328	INTEGER(i_std) :: i !! index (unitless)
2329	REAL(r_std) :: gdd_min !! critical gdd (C)
2330
2331	!_ ================================================================================================================================
2332
2333	IF (bavard.GE.2) WRITE(numout,*) 'Entering ncdgdd'
2334
2335	!! 1. Initializations
2336
2337	!! 1.1 initialize output
2338	begin_leaves(:,j) = .FALSE.
2339
2340	!! 1.2 check the critical value ::ncdgdd_temp is defined.
2341	!! If not, stop.
2342	IF ( ncdgdd_temp(j) .EQ. undef ) THEN
2343
2344	WRITE(numout,*) 'ncdgdd: ncdgdd_temp is undefined for PFT (::j) ',j
2345	WRITE(numout,*) 'We stop.'
2346	STOP
2347
2348	ENDIF
2349
2350	!! 2. Check if biometeorological conditions are favourable for leaf growth.
2351	!! PFT has to be there and start of growing season must be allowed.
2352	DO i = 1, npts ! loop over grid points
2353
2354	IF ( PFTpresent(i,j) .AND. allow_initpheno(i,j) .AND. &
2355	( gdd_midwinter(i,j) .NE. undef ) .AND. &
2356	( ncd_dormance(i,j) .NE. undef )) THEN
2357
2358	!! 2.1 Calculate the critical gdd, which is related to ::ncd_dormance
2359	!! using an empirical negative exponential law as described above.
2360	!! kpheno_crit was added to make use of the optimised phenology parameters of Natasha MacBean
2361	!! In ORCHIS, a multiplicative factor (opti_kpheno_crit) was optimised in stead of ggd_crit
2362	gdd_min = ( gddncd_ref / exp(gddncd_curvencd_dormance(i,j)) - gddncd_offset) opti_kpheno_crit(j)
2363
2364	!! 2.2 Determine if the growing season should start (if so, ::begin_leaves set to TRUE).
2365	!! This occurs if the critical GDD been reached AND the temperatures are increasing.
2366	!! If the growing season has started, ::gdd_midwinter is set to "undef".
2367	IF ( ( gdd_midwinter(i,j) .GE. gdd_min ) .AND. &
2368	( t2m_week(i) .GT. t2m_month(i) ) ) THEN
2369
2370	begin_leaves(i,j) = .TRUE.
2371	gdd_midwinter(i,j) = undef
2372
2373	ENDIF
2374
2375	ENDIF ! PFT there and start of growing season allowed
2376
2377	ENDDO ! end loop over grid points
2378
2379	IF (bavard.GE.2) WRITE(numout,*) 'Leaving ncdgdd'
2380
2381	END SUBROUTINE pheno_ncdgdd
2382
2383
2384	!! ================================================================================================================================
2385	!! SUBROUTINE : pheno_ngd
2386	!!
2387	!>\BRIEF The Number of Growing Days (NGD) leaf onset model initiates leaf onset if the NGD,
2388	!! defined as the number of days with temperature above a constant threshold,
2389	!! exceeds a critical value.
2390	!! Currently PFT 9 (Boreal Leedleleaf Summergreen) is assigned to this model.
2391	!!
2392	!! DESCRIPTION +++CHECK+++
2393	!! Description without dormance period include drormance (see ncdgdd)
2394	!! The NGD model is a variant of the GDD model. The model was proposed by Botta et
2395	!! al. (2000) for boreal and arctic biomes, and is designed to estimate
2396	!! leaf onset after the end of soil frost.
2397	!! The NDG (::ngd_minus5) is the number of days with a daily mean air
2398	!! temperature of greater than -5 degrees C,
2399	!! starting from the beginning of the dormancy period (i.e. time since the leaves
2400	!! were lost/GPP below a certain threshold).
2401	!! Leaf onset begins if the NGD is higher than the PFT-dependent constant threshold,
2402	!! ::ngd, AND if the weekly temperature is higher than the monthly
2403	!! temperature. \n
2404	!! The dormancy time-length is represented by the variable
2405	!! ::time_lowgpp, which is calculated in ::stomate_season. It is increased by
2406	!! the stomate time step when the weekly GPP is lower than a threshold. Otherwise
2407	!! it is set to zero. \n
2408	!! ::ngd_minus5 is also calculated in ::stomate_season. It is initialised at the
2409	!! beginning of the dormancy period (::time_lowgpp>0), and increased by the
2410	!! stomate time step when the temperature > -5 degrees C. \n
2411	!! ::ngd is set for each PFT in ::stomate_data, and a
2412	!! table defining the minimum NGD for each PFT is given in ::ngd_crit_tab
2413	!! in ::stomate_constants. \n
2414	!! The ::pheno_ngd subroutine is called in the subroutine ::phenology.
2415	!!
2416	!! RECENT CHANGE(S): None
2417	!!
2418	!! MAIN OUTPUT VARIABLE(S): ::begin_leaves - specifies whether leaf growth can start
2419	!!
2420	!! REFERENCE(S) :
2421	!! - Botta, A., N. Viovy, P. Ciais, P. Friedlingstein and P. Monfray (2000),
2422	!! A global prognostic scheme of leaf onset using satellite data,
2423	!! Global Change Biology, 207, 337-347.
2424	!! - Krinner, G., N. Viovy, N. de Noblet-Ducoudre, J. Ogee, J. Polcher, P.
2425	!! Friedlingstein, P. Ciais, S. Sitch and I.C. Prentice (2005), A dynamic global
2426	!! vegetation model for studies of the coupled atmosphere-biosphere system, Global
2427	!! Biogeochemical Cycles, 19, doi:10.1029/2003GB002199.
2428	!!
2429	!! FLOWCHART :
2430	!! \latexonly
2431	!! \includegraphics[scale = 1]{pheno_ngd.png}
2432	!! \endlatexonly
2433	!! \n
2434	!_ ================================================================================================================================
2435
2436	SUBROUTINE pheno_ngd (npts, j, PFTpresent, allow_initpheno, ncd_dormance, &
2437	ngd_minus5, t2m_month, t2m_week, begin_leaves)
2438
2439
2440	!! 0. Variable and parameter declaration
2441
2442	!! 0.1 Input variables
2443
2444	INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells (unitless)
2445	INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
2446	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
2447	LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to declare the beginning of the
2448	!! growing season? (true/false)
2449	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: ncd_dormance !! number of chilling days since leaves were lost
2450	!! (days)
2451	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ngd_minus5 !! growing degree days since midwinter (C)
2452	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures (K)
2453	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
2454
2455
2456	!! 0.2 Output variables
2457
2458	!! 0.3 Modified variables
2459
2460	LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting leaves on (true/false)
2461
2462
2463	!! 0.4 Local variables
2464
2465	INTEGER(i_std) :: i !! index (unitless)
2466	REAL(r_std) :: gdd_min !! critical gdd (C)
2467
2468	!_ ================================================================================================================================
2469
2470	IF (bavard.GE.2) WRITE(numout,*) 'Entering ngd'
2471
2472	!! 1. Initializations
2473
2474	!! 1.1 initialize output
2475	begin_leaves(:,j) = .FALSE.
2476
2477	!! 1.2 check the critical value ::ngd_crit is defined.
2478	! If not, stop.
2479	IF ( ngd_crit(j) .EQ. undef ) THEN
2480
2481	WRITE(numout,*) 'ngd: ngd_crit is undefined for PFT (::j) ',j
2482	WRITE(numout,*) 'We stop.'
2483	STOP
2484
2485	ENDIF
2486
2487	!+++TEMP+++
2488	IF (ld_pheno .AND. j.EQ.test_pft) THEN
2489	WRITE(numout,*) 'phenology PFTpresent, ', j, PFTpresent(:,j)
2490	WRITE(numout,*) 'phenology allow_initpheno, ', allow_initpheno(:,j)
2491	WRITE(numout,*) 'phenology ngd, ', ngd_minus5(:,j)
2492	WRITE(numout,*) 'phenology ngd_crit, ',ngd_crit(j)
2493	WRITE(numout,*) 'phenology ncd_dormance, ', ncd_dormance(:,j)
2494	WRITE(numout,*) 'phenology ngd_min_dormance, ',ngd_min_dormance
2495	WRITE(numout,*) 'phenology begin_leaves, ', begin_leaves(:,j)
2496	WRITE(numout,*) 'phenology t2m_week, ', t2m_week(:)
2497	WRITE(numout,*) 'phenology ncdgdd_temp, ', ncdgdd_temp(j)+ZeroCelsius
2498	WRITE(numout,*) 'phenology t2m_month, ', t2m_month(:)
2499	ENDIF
2500	!++++++++++
2501
2502
2503	!! 2. Check if biometeorological conditions are favourable for leaf growth.
2504	!! PFT has to be there and start of growing season must be allowed.
2505
2506
2507	DO i = 1, npts ! loop over grid points
2508
2509	IF ( PFTpresent(i,j) .AND. allow_initpheno(i,j) .AND. &
2510	( ngd_minus5(i,j) .NE. undef ) .AND. &
2511	( ncd_dormance(i,j) .NE. undef )) THEN
2512
2513	!! 2.1 Determine if the growing season should start
2514	! (if so, ::begin_leaves set to TRUE). This occurs if the critical GDD
2515	! been reached AND the temperatures are increasing. If the growing season
2516	! has started, ::gdd_midwinter is set to "undef".
2517
2518	IF ( ( ngd_minus5(i,j) .GE. ngd_crit(j)*opti_kpheno_crit(j) ) .AND. &
2519	( ( ncd_dormance(i,j) .GE. ngd_min_dormance) .OR. &
2520	( t2m_week(i) .GT. ncdgdd_temp(j)+ZeroCelsius) ) .AND. &
2521	( t2m_week(i) .GT. t2m_month(i) ) ) THEN
2522
2523	begin_leaves(i,j) = .TRUE.
2524	ngd_minus5(i,j) = undef
2525
2526	!+++TEMP+++
2527	IF (ld_pheno .AND. j.EQ.test_pft) THEN
2528	WRITE(numout,*) 'phenology updated begin_leaves, ', begin_leaves(i,j)
2529	ENDIF
2530	!++++++++++
2531
2532	ENDIF
2533
2534	ENDIF ! PFT there and start of growing season allowed
2535
2536	ENDDO ! end loop over grid points
2537
2538	IF (bavard.GE.2) WRITE(numout,*) 'Leaving ngd'
2539
2540	END SUBROUTINE pheno_ngd
2541
2542	!+++THE ORGINAL PHENOLOGY CODE++++
2543	! Note that ngd was not reset when begin_leaves was true. The begin leaves conditions
2544	! could be triggered several times per year. Furthermore, their is no forced dormance
2545	! period. This approach may have been OK in Siberia but it seems suboptimal for the DGVM
2546	! or even when prescribing larch (the only PFT with this phenology model) outside of
2547	! Siberia but still witin its growing range i.e. the mountains of Europe and N-America.
2548
2549
2550	!!$ SUBROUTINE pheno_ngd (npts, j, PFTpresent, allow_initpheno, ngd, &
2551	!!$ t2m_month, t2m_week, begin_leaves)
2552	!!$
2553	!!$ !
2554	!!$ !! 0. Variable and parameter declaration
2555	!!$ !
2556	!!$
2557	!!$ !
2558	!!$ !! 0.1 Input variables
2559	!!$ !
2560	!!$ INTEGER(i_std), INTENT(in) :: npts !! Domain size - number of grid cells (unitless)
2561	!!$ INTEGER(i_std), INTENT(in) :: j !! PFT index (unitless)
2562	!!$ LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: PFTpresent !! PFT exists (true/false)
2563	!!$ LOGICAL, DIMENSION(npts,nvm), INTENT(in) :: allow_initpheno !! are we allowed to declare the beginning of the
2564	!!$ !! growing season? (true/false)
2565	!!$ REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: ngd !! growing degree days (C)
2566	!!$ REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_month !! "monthly" 2-meter temperatures (K)
2567	!!$ REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_week !! "weekly" 2-meter temperatures (K)
2568	!!$
2569	!!$ !
2570	!!$ !! 0.2 Output variables
2571	!!$ !
2572	!!$
2573	!!$ !
2574	!!$ !! 0.3 Modified variables
2575	!!$ !
2576	!!$ LOGICAL, DIMENSION(npts,nvm), INTENT(inout) :: begin_leaves !! signal to start putting leaves on (true/false)
2577	!!$
2578	!!$ !
2579	!!$ !! 0.4 Local variables
2580	!!$ !
2581	!!$ INTEGER(i_std) :: i !! index (unitless)
2582	!!$
2583	!!$ !! =========================================================================
2584	!!$
2585	!!$ IF (bavard.GE.2) WRITE(numout,*) 'Entering ngd'
2586	!!$
2587	!!$ !
2588	!!$ !! 1. Initializations
2589	!!$ !
2590	!!$
2591	!!$ !
2592	!!$ !! 1.1 initialize output
2593	!!$ !
2594	!!$
2595	!!$ begin_leaves(:,j) = .FALSE.
2596	!!$
2597	!!$ !
2598	!!$ !! 1.2 check the critical value ::ngd_crit is defined.
2599	!!$ !! If not, stop.
2600	!!$ !
2601	!!$
2602	!!$ IF ( ngd_crit(j) .EQ. undef ) THEN
2603	!!$
2604	!!$ WRITE(numout,*) 'ngd: ngd_crit is undefined for PFT (::j) ',j
2605	!!$ WRITE(numout,*) 'We stop.'
2606	!!$ STOP
2607	!!$
2608	!!$ ENDIF
2609	!!$
2610	!!$ !
2611	!!$ !! 2. Check if biometeorological conditions are favourable for leaf growth.
2612	!!$ !! PFT has to be there and start of growing season must be allowed.
2613	!!$ !
2614	!!$
2615	!!$ DO i = 1, npts
2616	!!$
2617	!!$ IF ( PFTpresent(i,j) .AND. allow_initpheno(i,j) ) THEN
2618	!!$
2619	!!$ !! 2.1 Determine if the growing season should start (if so, ::begin_leaves set to TRUE).
2620	!!$ !! This occurs if the critical NGD has been reached AND are temperatures increasing.
2621	!!$
2622	!!$ IF ( ( ngd(i,j) .GE. ngd_crit(j) ) .AND. &
2623	!!$ ( t2m_week(i) .GT. t2m_month(i) ) ) THEN
2624	!!$ begin_leaves(i,j) = .TRUE.
2625	!!$ ENDIF
2626	!!$
2627	!!$ ENDIF ! PFT there and start of growing season allowed
2628	!!$
2629	!!$ ENDDO ! end loop over grid points
2630	!!$
2631	!!$ IF (bavard.GE.2) WRITE(numout,*) 'Leaving ngd'
2632	!!$
2633	!!$ END SUBROUTINE pheno_ngd
2634
2635	END MODULE stomate_phenology

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

Download in other formats: