Context Navigation

source: tags/ORCHIDEE_1_9_5_1/ORCHIDEE/src_stomate/stomate_season.f90 @ 2049

Last change on this file since 2049 was 42, checked in by mmaipsl, 14 years ago

MM: Replace all 0.0 by 'zero' and 1.0 by 'un',

and all 86400. by 'one_day' in the code to reduce explicit constants.

File size: 30.6 KB

Line
1	! Calculate long-term meteorological parameters from daily temperatures
2	! and precipitations (essentially for phenology)
3	!
4	! $Header: /home/ssipsl/CVSREP/ORCHIDEE/src_stomate/stomate_season.f90,v 1.15 2010/04/06 16:06:34 ssipsl Exp $
5	! IPSL (2006)
6	! This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
7	!
8	MODULE stomate_season
9
10	! modules used:
11
12	USE ioipsl
13	USE stomate_constants
14	USE constantes_veg
15
16	IMPLICIT NONE
17
18	! private & public routines
19
20	PRIVATE
21	PUBLIC season,season_clear
22
23	! first call
24	LOGICAL, SAVE :: firstcall = .TRUE.
25
26	CONTAINS
27
28
29	SUBROUTINE season_clear
30	firstcall=.TRUE.
31	END SUBROUTINE season_clear
32
33	SUBROUTINE season (npts, dt, EndOfYear, veget, veget_max, &
34	moiavail_daily, t2m_daily, tsoil_daily, soilhum_daily, &
35	precip_daily, npp_daily, biomass, turnover_daily, gpp_daily, when_growthinit, &
36	maxmoiavail_lastyear, maxmoiavail_thisyear, &
37	minmoiavail_lastyear, minmoiavail_thisyear, &
38	maxgppweek_lastyear, maxgppweek_thisyear, &
39	gdd0_lastyear, gdd0_thisyear, &
40	precip_lastyear, precip_thisyear, &
41	lm_lastyearmax, lm_thisyearmax, &
42	maxfpc_lastyear, maxfpc_thisyear, &
43	moiavail_month, moiavail_week, t2m_longterm, tlong_ref, t2m_month, t2m_week, &
44	tsoil_month, soilhum_month, &
45	npp_longterm, turnover_longterm, gpp_week, &
46	gdd_m5_dormance, gdd_midwinter, ncd_dormance, ngd_minus5, time_lowgpp, &
47	time_hum_min, hum_min_dormance, herbivores)
48
49	!
50	! 0 declarations
51	!
52
53	! 0.1 input
54
55	! Domain size
56	INTEGER(i_std), INTENT(in) :: npts
57	! time step in days
58	REAL(r_std), INTENT(in) :: dt
59	! update yearly variables?
60	LOGICAL, INTENT(in) :: EndOfYear
61	! coverage fraction of a PFT. Here: fraction of total ground.
62	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget
63	! "maximal" coverage fraction of a PFT (for LAI -> infinity)
64	! Here: fraction of total ground.
65	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: veget_max
66	! Daily moisture availability
67	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: moiavail_daily
68	! Daily 2 meter temperature (K)
69	REAL(r_std), DIMENSION(npts), INTENT(in) :: t2m_daily
70	! Daily soil temperature (K)
71	REAL(r_std), DIMENSION(npts,nbdl), INTENT(in) :: tsoil_daily
72	! Daily soil humidity
73	REAL(r_std), DIMENSION(npts,nbdl), INTENT(in) :: soilhum_daily
74	! Daily mean precipitation (mm/day)
75	REAL(r_std), DIMENSION(npts), INTENT(in) :: precip_daily
76	! daily net primary productivity (gC/m**2/day)
77	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: npp_daily
78	! biomass (gC/(m**2 of ground))
79	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in) :: biomass
80	! Turnover rates (gC/m**2/day)
81	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(in) :: turnover_daily
82	! daily gross primary productivity (Here: gC/(m**2 of total ground)/day)
83	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: gpp_daily
84	! how many days ago was the beginning of the growing season
85	REAL(r_std), DIMENSION(npts,nvm), INTENT(in) :: when_growthinit
86
87	! 0.2 modified fields
88
89	! last year's maximum moisture availability
90	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxmoiavail_lastyear
91	! this year's maximum moisture availability
92	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxmoiavail_thisyear
93	! last year's minimum moisture availability
94	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: minmoiavail_lastyear
95	! this year's minimum moisture availability
96	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: minmoiavail_thisyear
97	! last year's maximum weekly GPP
98	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxgppweek_lastyear
99	! this year's maximum weekly GPP
100	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxgppweek_thisyear
101	! last year's annual GDD0
102	REAL(r_std), DIMENSION(npts), INTENT(inout) :: gdd0_lastyear
103	! this year's annual GDD0
104	REAL(r_std), DIMENSION(npts), INTENT(inout) :: gdd0_thisyear
105	! last year's annual precipitation (mm/year)
106	REAL(r_std), DIMENSION(npts), INTENT(inout) :: precip_lastyear
107	! this year's annual precipitation (mm/year)
108	REAL(r_std), DIMENSION(npts), INTENT(inout) :: precip_thisyear
109	! last year's maximum leaf mass, for each PFT (gC/m**2)
110	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: lm_lastyearmax
111	! this year's maximum leaf mass, for each PFT (gC/m**2)
112	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: lm_thisyearmax
113	! last year's maximum fpc for each PFT, on ground
114	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxfpc_lastyear
115	! this year's maximum fpc for each PFT, on ground
116	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: maxfpc_thisyear
117	! "monthly" moisture availability
118	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: moiavail_month
119	! "weekly" moisture availability
120	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: moiavail_week
121	! "long term" 2-meter temperatures (K)
122	REAL(r_std), DIMENSION(npts), INTENT(inout) :: t2m_longterm
123	! "long term" refernce 2-meter temperatures (K)
124	REAL(r_std), DIMENSION(npts), INTENT(inout) :: tlong_ref
125	! "monthly" 2-meter temperatures (K)
126	REAL(r_std), DIMENSION(npts), INTENT(inout) :: t2m_month
127	! "weekly" 2-meter temperatures (K)
128	REAL(r_std), DIMENSION(npts), INTENT(inout) :: t2m_week
129	! "monthly" soil temperatures (K)
130	REAL(r_std), DIMENSION(npts,nbdl), INTENT(inout) :: tsoil_month
131	! "monthly" soil humidity
132	REAL(r_std), DIMENSION(npts,nbdl), INTENT(inout) :: soilhum_month
133	! "long term" net primary productivity (gC/m**2/year)
134	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: npp_longterm
135	! "long term" turnover rate (gC/m**2/year)
136	REAL(r_std), DIMENSION(npts,nvm,nparts), INTENT(inout) :: turnover_longterm
137	! "weekly" GPP (gC/day/(m**2 covered)
138	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gpp_week
139	! growing degree days, threshold -5 deg. C
140	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gdd_m5_dormance
141	! growing degree days since midwinter
142	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: gdd_midwinter
143	! number of chilling days since leaves were lost
144	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ncd_dormance
145	! number of growing days
146	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: ngd_minus5
147	! duration of dormance (d)
148	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: time_lowgpp
149	! time elapsed since strongest moisture availability (d)
150	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: time_hum_min
151	! minimum moisture during dormance
152	REAL(r_std), DIMENSION(npts,nvm), INTENT(inout) :: hum_min_dormance
153
154	! 0.3 output (diagnostic)
155
156	! time constant of probability of a leaf to be eaten by a herbivore (days)
157	REAL(r_std), DIMENSION(npts,nvm), INTENT(out) :: herbivores
158
159	! 0.4 local
160
161	! indices
162	INTEGER(i_std) :: j
163	! rapport maximal GPP/GGP_max pour dormance
164	REAL(r_std), PARAMETER :: gppfrac_dormance = 0.2
165	!
166	!NVADD
167	! minimum gpp considered as not "lowgpp"
168	REAL(r_std), PARAMETER :: min_gpp_allowed = 0.3
169	! tau (year) for "climatologic variables
170	REAL(r_std), PARAMETER :: tau_climatology = 20
171	!ENDNVADD
172	! maximum ncd (d) (to avoid floating point underflows)
173	REAL(r_std) :: ncd_max
174	! parameters for herbivore activity
175	REAL(r_std), PARAMETER :: hvc1 = 0.019
176	REAL(r_std), PARAMETER :: hvc2 = 1.38
177	REAL(r_std), PARAMETER :: leaf_frac=.33
178	! sum of natural fpcs
179	REAL(r_std), DIMENSION(npts) :: sumfpc_nat
180	! weights
181	REAL(r_std), DIMENSION(npts) :: weighttot
182	! natural long-term leaf NPP ( gC/m**2/year)
183	REAL(r_std), DIMENSION(npts) :: nlflong_nat
184	! residence time of green tissue (years)
185	REAL(r_std), DIMENSION(npts) :: green_age
186	! herbivore consumption (gC/m**2/day)
187	REAL(r_std), DIMENSION(npts) :: consumption
188
189	! =========================================================================
190
191	IF (bavard.GE.3) WRITE(numout,*) 'Entering season'
192
193	!
194	! 1 Initializations
195	!
196	ncd_max = 3. * one_year
197
198	IF ( firstcall ) THEN
199
200	!
201	! 1.1 messages
202	!
203
204	IF ( bavard .GE. 1 ) THEN
205
206	WRITE(numout,*) 'season: '
207
208	WRITE(numout,*) ' > rapport maximal GPP/GGP_max pour dormance: ',gppfrac_dormance
209
210	WRITE(numout,*) ' > maximum possible ncd (d): ',ncd_max
211
212	WRITE(numout,*) ' > herbivore consumption C (gC/m2/day) as a function of NPP (gC/m2/d):'
213	WRITE(numout,) ' C=',hvc1,' NPP^',hvc2
214	WRITE(numout,) ' > for herbivores, suppose that ',leaf_frac100., &
215	'% of NPP is allocated to leaves'
216
217
218	ENDIF
219
220	!
221	! 1.2 Check whether longer-term meteorological parameters are initialized
222	! to zero
223	!
224
225	! 1.2.1 moisture availabilities
226
227	! 1.2.1.1 "monthly"
228	!MM PAS PARALLELISE!!
229	IF ( ABS( SUM( moiavail_month(:,2:nvm) ) ) .LT. min_stomate ) THEN
230
231	! in this case, set them it today's moisture availability
232	WRITE(numout,*) 'Warning! We have to initialize the ''monthly'' moisture availabilities.'
233	moiavail_month(:,:) = moiavail_daily(:,:)
234
235	ENDIF
236
237	! 1.2.1.2 "weekly"
238
239	IF ( ABS( SUM( moiavail_week(:,2:nvm) ) ) .LT. min_stomate ) THEN
240
241	! in this case, set them it today's moisture availability
242	WRITE(numout,*) 'Warning! We have to initialize the ''weekly'' moisture availabilities.'
243	moiavail_week(:,:) = moiavail_daily(:,:)
244
245	ENDIF
246
247	! 1.2.2 2-meter temperatures
248
249	! 1.2.2.1 "long term"
250
251	IF ( ABS( SUM( t2m_longterm(:) ) ) .LT. min_stomate ) THEN
252
253	! in this case, set them to today's temperature
254	WRITE(numout,*) 'Warning! We have to initialize the ''long term'' 2m temperatures.'
255	t2m_longterm(:) = t2m_daily(:)
256
257	ENDIF
258
259	! 1.2.2.2 "monthly"
260
261	IF ( ABS( SUM( t2m_month(:) ) ) .LT. min_stomate ) THEN
262
263	! in this case, set them to today's temperature
264	WRITE(numout,*) 'Warning! We have to initialize the ''monthly'' 2m temperatures.'
265	t2m_month(:) = t2m_daily(:)
266
267	ENDIF
268
269	! 1.2.2.3 "weekly"
270
271	IF ( ABS( SUM( t2m_week(:) ) ) .LT. min_stomate ) THEN
272
273	! in this case, set them to today's temperature
274	WRITE(numout,*) 'Warning! We have to initialize the ''weekly'' 2m temperatures.'
275	t2m_week(:) = t2m_daily(:)
276
277	ENDIF
278
279	! 1.2.3 "monthly" soil temperatures
280	!MM PAS PARALLELISE!!
281	IF ( ABS( SUM( tsoil_month(:,:) ) ) .LT. min_stomate ) THEN
282
283	! in this case, set them to today's temperature
284	WRITE(numout,*) 'Warning!'// &
285	' We have to initialize the ''monthly'' soil temperatures.'
286	tsoil_month(:,:) = tsoil_daily(:,:)
287
288	ENDIF
289
290	! 1.2.4 "monthly" soil humidity
291	IF ( ABS( SUM( soilhum_month(:,:) ) ) .LT. min_stomate ) THEN
292
293	! in this case, set them to today's humidity
294	WRITE(numout,*) 'Warning!'// &
295	' We have to initialize the ''monthly'' soil humidity.'
296	soilhum_month(:,:) = soilhum_daily(:,:)
297
298	ENDIF
299
300	! 1.2.5 growing degree days, threshold -5 deg C
301	IF ( ABS( SUM( gdd_m5_dormance(:,2:nvm) ) ) .LT. min_stomate ) THEN
302	WRITE(numout,*) 'Warning! Growing degree days (-5 deg) are initialized to ''undef''.'
303	gdd_m5_dormance(:,:) = undef
304	ENDIF
305
306	! 1.2.6 growing degree days since midwinter
307	IF ( ABS( SUM( gdd_midwinter(:,2:nvm) ) ) .LT. min_stomate ) THEN
308	WRITE(numout,*) 'Warning! Growing degree days since midwinter' // &
309	' are initialized to ''undef''.'
310	gdd_midwinter(:,:) = undef
311	ENDIF
312
313	! 1.2.7 number of chilling days since leaves were lost
314	IF ( ABS( SUM( ncd_dormance(:,2:nvm) ) ) .LT. min_stomate ) THEN
315	WRITE(numout,*) 'Warning! Number of chilling days is initialized to ''undef''.'
316	ncd_dormance(:,:) = undef
317	ENDIF
318
319	! 1.2.8 number of growing days, threshold -5 deg C
320	IF ( ABS( SUM( ngd_minus5(:,2:nvm) ) ) .LT. min_stomate ) THEN
321	WRITE(numout,*) 'Warning! Number of growing days (-5 deg) is initialized to 0.'
322	ENDIF
323
324	! 1.2.9 "long term" npp
325	IF ( ABS( SUM( npp_longterm(:,2:nvm) ) ) .LT. min_stomate ) THEN
326	WRITE(numout,*) 'Warning! Long term NPP is initialized to 0.'
327	ENDIF
328
329	! 1.2.10 "long term" turnover
330	IF ( ABS( SUM( turnover_longterm(:,2:nvm,:) ) ) .LT. min_stomate ) THEN
331	WRITE(numout,*) 'Warning! Long term turnover is initialized to 0.'
332	ENDIF
333
334	! 1.2.11 "weekly" GPP
335	IF ( ABS( SUM( gpp_week(:,2:nvm) ) ) .LT. min_stomate ) THEN
336	WRITE(numout,*) 'Warning! Weekly GPP is initialized to 0.'
337	ENDIF
338
339	! 1.2.12 minimum moisture availabilities
340	IF ( ABS( SUM( minmoiavail_thisyear(:,2:nvm) ) ) .LT. min_stomate ) THEN
341
342	! in this case, set them to a very high value
343	WRITE(numout,*) 'Warning! We have to initialize this year''s minimum '// &
344	'moisture availabilities.'
345	minmoiavail_thisyear(:,:) = large_value
346
347	ENDIF
348
349	!
350	! 1.3 reset flag
351	!
352
353	firstcall = .FALSE.
354
355	ENDIF
356
357	!
358	! 2 moisture availabilities
359	!
360
361	!
362	! 2.1 "monthly"
363	!
364
365	moiavail_month = ( moiavail_month * ( pheno_crit%tau_hum_month - dt ) + &
366	moiavail_daily * dt ) / pheno_crit%tau_hum_month
367
368	DO j = 2,nvm
369	WHERE ( ABS(moiavail_month(:,j)) .LT. EPSILON(zero) )
370	moiavail_month(:,j) = zero
371	ENDWHERE
372	ENDDO
373
374	!
375	! 2.2 "weekly"
376	!
377
378	moiavail_week = ( moiavail_week * ( pheno_crit%tau_hum_week - dt ) + &
379	moiavail_daily * dt ) / pheno_crit%tau_hum_week
380
381	DO j = 2,nvm
382	WHERE ( ABS(moiavail_week(:,j)) .LT. EPSILON(zero) )
383	moiavail_week(:,j) = zero
384	ENDWHERE
385	ENDDO
386
387	!
388	! 3 2-meter temperatures
389	!
390
391	!
392	! 3.1 "long term"
393	!
394
395	t2m_longterm = ( t2m_longterm * ( pheno_crit%tau_longterm - dt ) + &
396	t2m_daily * dt ) / pheno_crit%tau_longterm
397
398	WHERE ( ABS(t2m_longterm(:)) .LT. EPSILON(zero) )
399	t2m_longterm(:) = zero
400	ENDWHERE
401
402	CALL histwrite (hist_id_stomate, 'T2M_LONGTERM', itime, &
403	t2m_longterm, npts, hori_index)
404	!
405	! 3.2 "long term reference"
406	! This temperature is used for recalculating PFT-specific parameters such as
407	! critical photosynthesis temperatures of critical GDDs for phenology. This
408	! means that if the reference temperature varies, the PFTs adapt to them.
409	! Therefore the reference temperature can vary only if the vegetation is not
410	! static.
411	!
412
413	!!$ IF (control%ok_dgvm) THEN
414	tlong_ref(:) = MAX( tlong_ref_min, MIN( tlong_ref_max, t2m_longterm(:) ) )
415	!!$ ENDIF
416	!
417	! 3.3 "monthly"
418	!
419
420	t2m_month = ( t2m_month * ( pheno_crit%tau_t2m_month - dt ) + &
421	t2m_daily * dt ) / pheno_crit%tau_t2m_month
422
423	WHERE ( ABS(t2m_month(:)) .LT. EPSILON(zero) )
424	t2m_month(:) = zero
425	ENDWHERE
426
427	!
428	! 3.4 "weekly"
429	!
430
431	t2m_week = ( t2m_week * ( pheno_crit%tau_t2m_week - dt ) + &
432	t2m_daily * dt ) / pheno_crit%tau_t2m_week
433
434	WHERE ( ABS(t2m_week(:)) .LT. EPSILON(zero) )
435	t2m_week(:) = zero
436	ENDWHERE
437
438	!
439	! 4 ''monthly'' soil temperatures
440	!
441
442	tsoil_month = ( tsoil_month * ( pheno_crit%tau_tsoil_month - dt ) + &
443	tsoil_daily(:,:) * dt ) / pheno_crit%tau_tsoil_month
444
445	WHERE ( ABS(tsoil_month(:,:)) .LT. EPSILON(zero) )
446	tsoil_month(:,:) = zero
447	ENDWHERE
448
449	!
450	! 5 ''monthly'' soil humidity
451	!
452
453	soilhum_month = ( soilhum_month * ( pheno_crit%tau_soilhum_month - dt ) + &
454	soilhum_daily * dt ) / pheno_crit%tau_soilhum_month
455
456	WHERE ( ABS(soilhum_month(:,:)) .LT. EPSILON(zero) )
457	soilhum_month(:,:) = zero
458	ENDWHERE
459
460	!
461	! 6 dormance (d)
462	! when gpp is low, increase dormance time. Otherwise, set it to zero.
463	! NV: special case (3rd condition): plant is accumulating carbohydrates
464	! and does never use them. In this case, we allow the plant to
465	! detect a beginning of the growing season by declaring it dormant
466	!
467	!NVMODIF
468	DO j = 2,nvm
469	WHERE ( ( gpp_week(:,j) .LT. min_gpp_allowed ) .OR. &
470	( gpp_week(:,j) .LT. gppfrac_dormance * maxgppweek_lastyear(:,j) ) .OR. &
471	( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. &
472	( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) )
473	! WHERE ( ( gpp_week(:,j) .EQ. zero ) .OR. &
474	! ( gpp_week(:,j) .LT. gppfrac_dormance * maxgppweek_lastyear(:,j) ) .OR. &
475	! ( ( when_growthinit(:,j) .GT. 2.*one_year ) .AND. &
476	! ( biomass(:,j,icarbres) .GT. biomass(:,j,ileaf)*4. ) ) )
477
478	time_lowgpp(:,j) = time_lowgpp(:,j) + dt
479
480	ELSEWHERE
481
482	time_lowgpp(:,j) = zero
483
484	ENDWHERE
485	ENDDO
486
487	!
488	! 7 growing degree days, threshold -5 deg C
489	!
490
491	DO j = 2,nvm
492
493	! only for PFTs for which critical gdd is defined
494	! gdd_m5_dormance is set to 0 at the end of the growing season. It is set to undef
495	! at the beginning of the growing season.
496
497	IF ( ALL(pheno_crit%gdd(j,:) .NE. undef) ) THEN
498
499	!
500	! 7.1 set to zero if undef and no gpp
501	!
502
503	WHERE ( ( time_lowgpp(:,j) .GT. zero ) .AND. ( gdd_m5_dormance(:,j) .EQ. undef ) )
504
505	gdd_m5_dormance(:,j) = zero
506
507	ENDWHERE
508
509	!
510	! 7.2 set to undef if there is gpp
511	!
512
513	WHERE ( time_lowgpp(:,j) .EQ. zero )
514
515	gdd_m5_dormance(:,j) = undef
516
517	ENDWHERE
518
519	!
520	! 7.3 normal update where gdd_m5_dormance is defined
521	!
522
523	WHERE ( ( t2m_daily(:) .GT. (ZeroCelsius-5.) ) .AND. &
524	( gdd_m5_dormance(:,j) .NE. undef ) )
525	gdd_m5_dormance(:,j) = gdd_m5_dormance(:,j) + &
526	dt * ( t2m_daily(:) - (ZeroCelsius-5.) )
527	ENDWHERE
528
529	WHERE ( gdd_m5_dormance(:,j) .NE. undef )
530	gdd_m5_dormance(:,j) = gdd_m5_dormance(:,j) * &
531	( pheno_crit%tau_gdd - dt ) / pheno_crit%tau_gdd
532	ENDWHERE
533
534	ENDIF
535
536	ENDDO
537
538	DO j = 2,nvm
539	WHERE ( ABS(gdd_m5_dormance(:,j)) .LT. EPSILON(zero) )
540	gdd_m5_dormance(:,j) = zero
541	ENDWHERE
542	ENDDO
543
544	!
545	! 8 growing degree days since midwinter
546	!
547
548	DO j = 2,nvm
549
550	! only for PFTs for which ncdgdd_crittemp is defined
551
552	IF ( pheno_crit%ncdgdd_temp(j) .NE. undef ) THEN
553
554	!
555	! 8.1 set to 0 if undef and if we detect "midwinter"
556	!
557
558	WHERE ( ( gdd_midwinter(:,j) .EQ. undef ) .AND. &
559	( t2m_month(:) .LT. t2m_week(:) ) .AND. &
560	( t2m_month(:) .LT. t2m_longterm(:) ) )
561
562	gdd_midwinter(:,j) = zero
563
564	ENDWHERE
565
566	!
567	! 8.2 set to undef if we detect "midsummer"
568	!
569
570	WHERE ( ( t2m_month(:) .GT. t2m_week(:) ) .AND. &
571	( t2m_month(:) .GT. t2m_longterm(:) ) )
572
573	gdd_midwinter(:,j) = undef
574
575	ENDWHERE
576
577	!
578	! 8.3 normal update
579	!
580
581	WHERE ( ( gdd_midwinter(:,j) .NE. undef ) .AND. &
582	( t2m_daily(:) .GT. pheno_crit%ncdgdd_temp(j)+ZeroCelsius ) )
583
584	gdd_midwinter(:,j) = &
585	gdd_midwinter(:,j) + &
586	dt * ( t2m_daily(:) - ( pheno_crit%ncdgdd_temp(j)+ZeroCelsius ) )
587
588	ENDWHERE
589
590	ENDIF
591
592	ENDDO
593
594	!
595	! 9 number of chilling days since leaves were lost
596	!
597
598	DO j = 2,nvm
599
600	IF ( pheno_crit%ncdgdd_temp(j) .NE. undef ) THEN
601
602	!
603	! 9.1 set to zero if undef and no gpp
604	!
605
606	WHERE ( ( time_lowgpp(:,j) .GT. zero ) .AND. ( ncd_dormance(:,j) .EQ. undef ) )
607
608	ncd_dormance(:,j) = zero
609
610	ENDWHERE
611
612	!
613	! 9.2 set to undef if there is gpp
614	!
615
616	WHERE ( time_lowgpp(:,j) .EQ. zero )
617
618	ncd_dormance(:,j) = undef
619
620	ENDWHERE
621
622	!
623	! 9.3 normal update where ncd_dormance is defined
624	!
625
626	WHERE ( ( ncd_dormance(:,j) .NE. undef ) .AND. &
627	( t2m_daily(:) .LE. pheno_crit%ncdgdd_temp(j)+ZeroCelsius ) )
628
629	ncd_dormance(:,j) = MIN( ncd_dormance(:,j) + dt, ncd_max )
630
631	ENDWHERE
632
633	ENDIF
634
635	ENDDO
636
637	!
638	! 10 number of growing days, threshold -5 deg C
639	!
640
641	DO j = 2,nvm
642
643	!
644	! 10.1 Where there is GPP, set ngd to 0
645	! This means that we only take into account ngds when the leaves are off
646	!
647
648	WHERE ( time_lowgpp(:,j) .LT. min_stomate )
649	ngd_minus5(:,j) = zero
650	ENDWHERE
651
652	!
653	! 10.2 normal update
654	!
655
656	WHERE ( t2m_daily(:) .GT. (ZeroCelsius-5.) )
657	ngd_minus5(:,j) = ngd_minus5(:,j) + dt
658	ENDWHERE
659
660	ngd_minus5(:,j) = ngd_minus5(:,j) * ( pheno_crit%tau_ngd - dt ) / pheno_crit%tau_ngd
661
662	ENDDO
663
664	DO j = 2,nvm
665	WHERE ( ABS(ngd_minus5(:,j)) .LT. EPSILON(zero) )
666	ngd_minus5(:,j) = zero
667	ENDWHERE
668	ENDDO
669	!
670	! 11 minimum humidity since dormance began and time elapsed since this minimum
671	!
672
673	DO j = 2,nvm
674
675	IF ( pheno_crit%hum_min_time(j) .NE. undef ) THEN
676
677	!
678	! 11.1 initialize if undef and no gpp
679	!
680
681	WHERE ( ( time_lowgpp(:,j) .GT. zero ) .AND. &
682	( ( time_hum_min(:,j) .EQ. undef ) .OR. ( hum_min_dormance(:,j) .EQ. undef ) ) )
683
684	time_hum_min(:,j) = zero
685	hum_min_dormance(:,j) = moiavail_month(:,j)
686
687	ENDWHERE
688
689	!
690	! 11.2 set to undef where there is gpp
691	!
692
693	WHERE ( time_lowgpp(:,j) .EQ. zero )
694
695	time_hum_min(:,j) = undef
696	hum_min_dormance(:,j) = undef
697
698	ENDWHERE
699
700	!
701	! 11.3 normal update where time_hum_min and hum_min_dormance are defined
702	!
703
704	! 11.3.1 increase time counter
705
706	WHERE ( ( time_hum_min(:,j) .NE. undef ) .AND. &
707	( hum_min_dormance(:,j) .NE. undef ) )
708
709	time_hum_min(:,j) = time_hum_min(:,j) + dt
710
711	ENDWHERE
712
713	! 11.3.2 set time to zero if minimum is reached
714
715	WHERE ( ( time_hum_min(:,j) .NE. undef ) .AND. &
716	( hum_min_dormance(:,j) .NE. undef ) .AND. &
717	( moiavail_month(:,j) .LE. hum_min_dormance(:,j) ) )
718
719	hum_min_dormance(:,j) = moiavail_month(:,j)
720	time_hum_min(:,j) = zero
721
722	ENDWHERE
723
724	ENDIF
725
726	ENDDO
727
728	!
729	! 12 "long term" NPP. npp_daily in gC/m2/day, npp_longterm in gC/m2/year.
730	!
731
732	npp_longterm = ( npp_longterm * ( pheno_crit%tau_longterm - dt ) + &
733	(npp_dailyone_year) dt ) / &
734	pheno_crit%tau_longterm
735
736	DO j = 2,nvm
737	WHERE ( ABS(npp_longterm(:,j)) .LT. EPSILON(zero) )
738	npp_longterm(:,j) = zero
739	ENDWHERE
740	ENDDO
741
742	!
743	! 13 "long term" turnover rates, in gC/m**2/year.
744	!
745
746	turnover_longterm = ( turnover_longterm * ( pheno_crit%tau_longterm - dt ) + &
747	(turnover_dailyone_year) dt ) / &
748	pheno_crit%tau_longterm
749
750	DO j = 2,nvm
751	WHERE ( ABS(turnover_longterm(:,j,:)) .LT. EPSILON(zero) )
752	turnover_longterm(:,j,:) = zero
753	ENDWHERE
754	ENDDO
755
756	!
757	! 14 "weekly" GPP, in gC/(m**2 covered)/day (!)
758	! i.e. divide daily gpp (in gC/m**2 of total ground/day) by vegetation fraction
759	! (m2 covered/m2 of total ground)
760	!
761
762	WHERE ( veget_max .GT. zero )
763
764	gpp_week = ( gpp_week * ( pheno_crit%tau_gpp_week - dt ) + &
765	gpp_daily * dt ) / pheno_crit%tau_gpp_week
766
767	ELSEWHERE
768
769	gpp_week = zero
770
771	ENDWHERE
772
773	DO j = 2,nvm
774	WHERE ( ABS(gpp_week(:,j)) .LT. EPSILON(zero) )
775	gpp_week(:,j) = zero
776	ENDWHERE
777	ENDDO
778
779	!
780	! 15 maximum and minimum moisture availabilities
781	!
782
783	WHERE ( moiavail_daily .GT. maxmoiavail_thisyear )
784	maxmoiavail_thisyear = moiavail_daily
785	ENDWHERE
786
787	WHERE ( moiavail_daily .LT. minmoiavail_thisyear )
788	minmoiavail_thisyear = moiavail_daily
789	ENDWHERE
790
791	!
792	! 16 annual maximum weekly GPP
793	!
794
795	WHERE ( gpp_week .GT. maxgppweek_thisyear )
796	maxgppweek_thisyear = gpp_week
797	ENDWHERE
798
799	!
800	! 17 annual GDD0
801	!
802
803	WHERE ( t2m_daily .GT. ZeroCelsius )
804	gdd0_thisyear = gdd0_thisyear + dt * ( t2m_daily - ZeroCelsius )
805	ENDWHERE
806
807	!
808	! 18 annual precipitation
809	!
810
811	precip_thisyear = precip_thisyear + dt * precip_daily
812
813	!
814	! 19 annual maximum leaf mass
815	! If STOMATE is not activated, this corresponds to the maximum possible
816	! LAI of the PFT
817	!
818
819	IF ( control%ok_stomate ) THEN
820
821	DO j = 2,nvm
822	WHERE ( biomass(:,j,ileaf) .GT. lm_thisyearmax(:,j) )
823	lm_thisyearmax(:,j) = biomass(:,j,ileaf)
824	ENDWHERE
825	ENDDO
826
827	ELSE
828
829	DO j = 2,nvm
830	lm_thisyearmax(:,j) = lai_max(j) / sla(j)
831	ENDDO
832
833	ENDIF
834
835	!
836	! 20 annual maximum fpc for each PFT
837	! "veget" is defined as fraction of total ground. Therefore, maxfpc_thisyear has
838	! the same unit.
839	!
840
841	WHERE ( veget(:,:) .GT. maxfpc_thisyear(:,:) )
842	maxfpc_thisyear(:,:) = veget(:,:)
843	ENDWHERE
844
845	!
846	! 21 Every year, replace last year's maximum and minimum moisture availability,
847	! annual GDD0, annual precipitation, annual max weekly GPP, and maximum leaf mass
848
849	IF ( EndOfYear ) THEN
850
851	!
852	! 21.1 replace old values
853	!
854	!NVMODIF
855	maxmoiavail_lastyear(:,:) = (maxmoiavail_lastyear(:,:)*(tau_climatology-1)+ maxmoiavail_thisyear(:,:))/tau_climatology
856	minmoiavail_lastyear(:,:) = (minmoiavail_lastyear(:,:)*(tau_climatology-1)+ minmoiavail_thisyear(:,:))/tau_climatology
857	maxgppweek_lastyear(:,:) =( maxgppweek_lastyear(:,:)*(tau_climatology-1)+ maxgppweek_thisyear(:,:))/tau_climatology
858	! maxmoiavail_lastyear(:,:) = maxmoiavail_thisyear(:,:)
859	! minmoiavail_lastyear(:,:) = minmoiavail_thisyear(:,:)
860	! maxgppweek_lastyear(:,:) = maxgppweek_thisyear(:,:)
861
862	gdd0_lastyear(:) = gdd0_thisyear(:)
863
864	precip_lastyear(:) = precip_thisyear(:)
865
866	lm_lastyearmax(:,:) = lm_thisyearmax(:,:)
867
868	maxfpc_lastyear(:,:) = maxfpc_thisyear(:,:)
869
870	!
871	! 21.2 reset new values
872	!
873
874	maxmoiavail_thisyear(:,:) = zero
875	minmoiavail_thisyear(:,:) = large_value
876
877	maxgppweek_thisyear(:,:) = zero
878
879	gdd0_thisyear(:) = zero
880
881	precip_thisyear(:) = zero
882
883	lm_thisyearmax(:,:) = zero
884
885	maxfpc_thisyear(:,:) = zero
886
887	!
888	! 21.3 Special treatment for maxfpc.
889	!
890
891	!
892	! 21.3.1 Only take into account natural PFTs
893	!
894
895	DO j = 2,nvm
896	IF ( .NOT. natural(j) ) THEN
897	maxfpc_lastyear(:,j) = zero
898	ENDIF
899	ENDDO
900
901	! 21.3.2 In Stomate, veget is defined as a fraction of ground, not as a fraction
902	! of total ground. maxfpc_lastyear will be compared to veget in lpj_light.
903	! Therefore, we have to transform maxfpc_lastyear.
904
905
906	! 21.3.3 The sum of the maxfpc_lastyear for natural PFT must not exceed fpc_crit (=.95).
907	! However, it can slightly exceed this value as not all PFTs reach their maximum
908	! fpc at the same time. Therefore, if sum(maxfpc_lastyear) for the natural PFTs
909	! exceeds fpc_crit, we scale the values of maxfpc_lastyear so that the sum is
910	! fpc_crit.
911
912	! calculate the sum of maxfpc_lastyear
913	sumfpc_nat(:) = zero
914	DO j = 2,nvm
915	sumfpc_nat(:) = sumfpc_nat(:) + maxfpc_lastyear(:,j)
916	ENDDO
917
918	! scale so that the new sum is fpc_crit
919	DO j = 2,nvm
920	WHERE ( sumfpc_nat(:) .GT. fpc_crit )
921	maxfpc_lastyear(:,j) = maxfpc_lastyear(:,j) * (fpc_crit/sumfpc_nat(:))
922	ENDWHERE
923	ENDDO
924
925	ENDIF ! EndOfYear
926
927	!
928	! 22 diagnose herbivore activity, determined through as probability for a leaf to be
929	! eaten in a day
930	! Follows McNaughton et al., Nat 341, 142-144, 1989.
931	!
932
933	!
934	! 22.1 first calculate mean long-term leaf NPP in grid box, mean residence
935	! time (years) of green tissue (i.e. tissue that will be eaten by
936	! herbivores) (crudely approximated: 6 months for seasonal and 2 years
937	! for evergreen) and mean length of growing season (6 months for
938	! seasonal and 1 year for evergreen).
939	!
940
941	nlflong_nat(:) = zero
942	weighttot(:) = zero
943	green_age(:) = zero
944	!
945	DO j = 2,nvm
946	!
947	IF ( natural(j) ) THEN
948	!
949	weighttot(:) = weighttot(:) + lm_lastyearmax(:,j)
950	nlflong_nat(:) = nlflong_nat(:) + npp_longterm(:,j) * leaf_frac
951	!
952	IF ( pheno_crit%pheno_model(j) .EQ. 'none' ) THEN
953	green_age(:) = green_age(:) + 2. * lm_lastyearmax(:,j)
954	ELSE
955	green_age(:) = green_age(:) + .5 * lm_lastyearmax(:,j)
956	ENDIF
957	!
958	ENDIF
959	!
960	ENDDO
961	!
962	WHERE ( weighttot(:) .GT. zero )
963	green_age(:) = green_age(:) / weighttot(:)
964	ELSEWHERE
965	green_age(:) = un
966	ENDWHERE
967
968	!
969	! 22.2 McNaughton et al. give herbivore consumption as a function of annual leaf NPP.
970	! The annual leaf NPP can give us an idea about the edible biomass:
971	!
972
973	DO j = 2,nvm
974	!
975	IF ( natural(j) ) THEN
976	!
977	WHERE ( nlflong_nat(:) .GT. zero )
978	consumption(:) = hvc1 * nlflong_nat(:) ** hvc2
979	herbivores(:,j) = one_year * green_age(:) * nlflong_nat(:) / consumption(:)
980	ELSEWHERE
981	herbivores(:,j) = 100000.
982	ENDWHERE
983	!
984	ELSE
985	!
986	herbivores(:,j) = 100000.
987	!
988	ENDIF
989	!
990	ENDDO
991	herbivores(:,ibare_sechiba) = zero
992
993	IF (bavard.GE.4) WRITE(numout,*) 'Leaving season'
994
995	END SUBROUTINE season
996
997	END MODULE stomate_season

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

Download in other formats: