Context Navigation

source: branches/ORCHIDEE_EXT/ORCHIDEE_OL/teststomate.f90 @ 258

Last change on this file since 258 was 258, checked in by didier.solyga, 13 years ago
Externalized version of ORCHIDEE_OL files merged with the trunk
File size: 35.0 KB

Line
1	PROGRAM teststomate
2	!---------------------------------------------------------------------
3	! This program reads a forcing file for STOMATE
4	! which was created with STOMATE coupled to SECHIBA.
5	! It then integrates STOMATE for a
6	! given number of years using this forcing.
7	! The GPP is scaled to the new calculated LAI...
8	! Nothing is done for soil moisture which should
9	! actually evolve with the vegetation.
10	!---------------------------------------------------------------------
11	!- IPSL (2006)
12	!- This software is governed by the CeCILL licence see ORCHIDEE/ORCHIDEE_CeCILL.LIC
13	USE netcdf
14	!-
15	USE defprec
16	USE constantes
17	USE pft_parameters
18	USE stomate_data
19	USE ioipsl
20	USE grid
21	USE slowproc
22	USE stomate
23	USE intersurf, ONLY: stom_define_history , intsurf_time
24	USE parallel
25	!-
26	IMPLICIT NONE
27	!-
28	! Declarations
29	!-
30	INTEGER(i_std) :: vegax_id
31	INTEGER(i_std) :: kjpij,kjpindex
32	REAL(r_std) :: dtradia,dt_force
33	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indices
34	INTEGER(i_std),DIMENSION(:),ALLOCATABLE :: indexveg
35	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: soiltype, veget_x
36	REAL(r_std),DIMENSION(:),ALLOCATABLE :: totfrac_nobio
37	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: frac_nobio
38	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: veget_max_x
39	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: lai_x
40	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: veget_force_x
41	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: veget_max_force_x
42	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: lai_force_x
43	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: lon,lat
44	REAL(r_std),DIMENSION(:),ALLOCATABLE :: t2m,t2m_min,temp_sol
45	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: soiltemp,soilhum
46	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: humrel_x
47	REAL(r_std),DIMENSION(:),ALLOCATABLE :: litterhum
48	REAL(r_std),DIMENSION(:),ALLOCATABLE :: precip_rain,precip_snow
49	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: gpp_x
50	REAL(r_std),DIMENSION(:),ALLOCATABLE :: deadleaf_cover
51	REAL(r_std),DIMENSION(:,:,:),ALLOCATABLE :: assim_param_x
52	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: height_x
53	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: qsintmax_x
54	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: co2_flux
55	!!$ DS add for externalised version
56	REAL(r_std),DIMENSION(:),ALLOCATABLE :: hist_PFTaxis
57	!!$ DS
58
59
60	INTEGER :: ier,iret
61	INTEGER :: ncfid
62	LOGICAL :: a_er
63	CHARACTER(LEN=80) :: &
64	& dri_restname_in,dri_restname_out, &
65	& sec_restname_in,sec_restname_out, &
66	& sto_restname_in,sto_restname_out, stom_histname
67	INTEGER(i_std) :: iim,jjm
68	INTEGER(i_std),PARAMETER :: llm = 1
69	REAL(r_std),DIMENSION(llm) :: lev
70	REAL(r_std) :: dt_files
71	INTEGER(i_std) :: itau_dep,itau,itau_len,itau_step
72	REAL(r_std) :: date0
73	INTEGER(i_std) :: rest_id_sec,rest_id_sto
74	INTEGER(i_std) :: hist_id_sec,hist_id_sec2,hist_id_sto,hist_id_stom_IPCC
75	CHARACTER(LEN=30) :: time_str
76	REAL :: hist_days_stom,hist_dt_stom
77	!!$ REAL,DIMENSION(nvm) :: hist_PFTaxis
78	REAL(r_std),DIMENSION(10) :: hist_pool_10axis
79	REAL(r_std),DIMENSION(100) :: hist_pool_100axis
80	REAL(r_std),DIMENSION(11) :: hist_pool_11axis
81	REAL(r_std),DIMENSION(101) :: hist_pool_101axis
82	INTEGER :: hist_PFTaxis_id,hori_id
83	INTEGER :: hist_pool_10axis_id
84	INTEGER :: hist_pool_100axis_id
85	INTEGER :: hist_pool_11axis_id
86	INTEGER :: hist_pool_101axis_id
87	INTEGER :: hist_level
88	INTEGER,PARAMETER :: max_hist_level = 10
89	INTEGER(i_std) :: i,j,iv,id
90	CHARACTER*80 :: var_name
91	CHARACTER(LEN=40),DIMENSION(10) :: fluxop
92	LOGICAL :: ldrestart_read,ldrestart_write
93	LOGICAL :: l1d
94	INTEGER(i_std),PARAMETER :: nbvarmax=200
95	INTEGER(i_std) :: nbvar
96	CHARACTER*50,DIMENSION(nbvarmax) :: varnames
97	INTEGER(i_std) :: varnbdim
98	INTEGER(i_std),PARAMETER :: varnbdim_max=20
99	INTEGER,DIMENSION(varnbdim_max) :: vardims
100	CHARACTER*200 :: taboo_vars
101	REAL(r_std),DIMENSION(1) :: xtmp
102	INTEGER :: nsfm,nsft
103	INTEGER :: iisf
104	INTEGER(i_std) :: max_totsize,totsize_1step
105	INTEGER(i_std),DIMENSION(0:2) :: ifirst,ilast
106	INTEGER(i_std) :: iblocks,nblocks
107	INTEGER,PARAMETER :: ndm = 10
108	INTEGER,DIMENSION(ndm) :: start,count
109	INTEGER :: ndim,v_id
110	INTEGER :: force_id
111	CHARACTER(LEN=100) :: forcing_name
112	REAL :: x
113	REAL(r_std),DIMENSION(:),ALLOCATABLE :: x_indices
114	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: x_neighbours
115	REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: var_3d
116	!-
117	REAL(r_std) :: time_sec,time_step_sec
118	REAL(r_std) :: time_dri,time_step_dri
119	REAL(r_std),DIMENSION(1) :: r1d
120	REAL(r_std) :: julian,djulian
121	! REAL(r_std),DIMENSION(:,:),ALLOCATABLE :: soiltype
122	!-
123	! the following variables contain the forcing data
124	!-
125	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: clay_fm
126	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: humrel_x_fm
127	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: litterhum_fm
128	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: t2m_fm
129	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: t2m_min_fm
130	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: temp_sol_fm
131	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: soiltemp_fm
132	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: soilhum_fm
133	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:) :: precip_fm
134	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: gpp_x_fm
135	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: veget_force_x_fm
136	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: veget_max_force_x_fm
137	REAL(r_std),ALLOCATABLE,SAVE,DIMENSION(:,:,:) :: lai_force_x_fm
138	INTEGER(i_std),ALLOCATABLE,SAVE,DIMENSION(:) :: isf
139	LOGICAL,ALLOCATABLE,SAVE,DIMENSION(:) :: nf_written
140	INTEGER(i_std),ALLOCATABLE,SAVE,DIMENSION(:) :: nf_cumul
141	INTEGER(i_std) :: ji,jv
142	!---------------------------------------------------------------------
143	!- No parallelisation yet in teststomate !
144	#ifdef CPP_PARA
145	CALL ipslerr (3,'teststomate', &
146	& 'You try to run testsomate compiled with parallelisation. (CPP_PARA key)', &
147	& 'But it wasn''t programmed yet and teststomate will stop.','You must compiled it without CPP_PARA key.')
148	#endif
149
150	!-
151	! calendar
152	!-
153	CALL ioconf_calendar ('noleap')
154	CALL ioget_calendar (one_year,one_day)
155	!-
156	! open STOMATE's forcing file to read some basic info
157	!-
158	forcing_name = 'stomate_forcing.nc'
159	CALL getin ('STOMATE_FORCING_NAME',forcing_name)
160	iret = NF90_OPEN (TRIM(forcing_name),NF90_NOWRITE,force_id)
161	IF (iret /= NF90_NOERR) THEN
162	CALL ipslerr (3,'teststomate', &
163	& 'Could not open file : ', &
164	& forcing_name,'(Do you have forget it ?)')
165	ENDIF
166	ier = NF90_GET_ATT (force_id,NF90_GLOBAL,'dtradia',dtradia)
167	ier = NF90_GET_ATT (force_id,NF90_GLOBAL,'dt_slow',dt_force)
168	ier = NF90_GET_ATT (force_id,NF90_GLOBAL,'nsft',x)
169	nsft = NINT(x)
170	ier = NF90_GET_ATT (force_id,NF90_GLOBAL,'kjpij',x)
171	kjpij = NINT(x)
172	ier = NF90_GET_ATT (force_id,NF90_GLOBAL,'kjpindex',x)
173	kjpindex = NINT(x)
174	CALL init_grid( kjpindex )
175	!-
176	write(,) 'ATTENTION',dtradia,dt_force
177
178	!!$ DS added for the externalised version
179	!!$we need to know the pft parameters values used for stomate
180	CALL pft_parameters_main
181	CALL getin_stomate_pft_parameters
182	!!$ DS
183
184	!-
185	! allocate arrays
186	!-
187	a_er = .FALSE.
188	!!$ DS added for the externalised version
189	!!$ ALLOCATE (pft_to_mtc(nvm),stat=ier)
190	!!$ a_er = a_er .OR. (ier.NE.0)
191	ALLOCATE (hist_PFTaxis(nvm),stat=ier)
192	a_er = a_er .OR. (ier.NE.0)
193	!!$ end add DS
194	ALLOCATE (indices(kjpindex),stat=ier)
195	a_er = a_er .OR. (ier.NE.0)
196	ALLOCATE (indexveg(kjpindex*nvm), stat=ier)
197	a_er = a_er .OR. (ier.NE.0)
198	ALLOCATE (soiltype(kjpindex,nstm),stat=ier)
199	a_er = a_er .OR. (ier.NE.0)
200	ALLOCATE (veget_x(kjpindex,nvm),stat=ier)
201	a_er = a_er .OR. (ier.NE.0)
202	ALLOCATE (totfrac_nobio(kjpindex),stat=ier)
203	a_er = a_er .OR. (ier.NE.0)
204	ALLOCATE (frac_nobio(kjpindex,nnobio),stat=ier)
205	a_er = a_er .OR. (ier.NE.0)
206	ALLOCATE (veget_max_x(kjpindex,nvm),stat=ier)
207	a_er = a_er .OR. (ier.NE.0)
208	ALLOCATE (lai_x(kjpindex,nvm),stat=ier)
209	a_er = a_er .OR. (ier.NE.0)
210	ALLOCATE (veget_force_x(kjpindex,nvm),stat=ier)
211	a_er = a_er .OR. (ier.NE.0)
212	ALLOCATE (veget_max_force_x(kjpindex,nvm),stat=ier)
213	a_er = a_er .OR. (ier.NE.0)
214	ALLOCATE (lai_force_x(kjpindex,nvm),stat=ier)
215	a_er = a_er .OR. (ier.NE.0)
216	ALLOCATE (t2m(kjpindex),stat=ier)
217	a_er = a_er .OR. (ier.NE.0)
218	ALLOCATE (t2m_min(kjpindex),stat=ier)
219	a_er = a_er .OR. (ier.NE.0)
220	ALLOCATE (temp_sol(kjpindex),stat=ier)
221	a_er = a_er .OR. (ier.NE.0)
222	ALLOCATE (soiltemp(kjpindex,nbdl),stat=ier)
223	a_er = a_er .OR. (ier.NE.0)
224	ALLOCATE (soilhum(kjpindex,nbdl),stat=ier)
225	a_er = a_er .OR. (ier.NE.0)
226	ALLOCATE (humrel_x(kjpindex,nvm),stat=ier)
227	a_er = a_er .OR. (ier.NE.0)
228	ALLOCATE (litterhum(kjpindex),stat=ier)
229	a_er = a_er .OR. (ier.NE.0)
230	ALLOCATE (precip_rain(kjpindex),stat=ier)
231	a_er = a_er .OR. (ier.NE.0)
232	ALLOCATE (precip_snow(kjpindex),stat=ier)
233	a_er = a_er .OR. (ier.NE.0)
234	ALLOCATE (gpp_x(kjpindex,nvm),stat=ier)
235	a_er = a_er .OR. (ier.NE.0)
236	ALLOCATE (deadleaf_cover(kjpindex),stat=ier)
237	a_er = a_er .OR. (ier.NE.0)
238	ALLOCATE (assim_param_x(kjpindex,nvm,npco2),stat=ier)
239	a_er = a_er .OR. (ier.NE.0)
240	ALLOCATE (height_x(kjpindex,nvm),stat=ier)
241	a_er = a_er .OR. (ier.NE.0)
242	ALLOCATE (qsintmax_x(kjpindex,nvm),stat=ier)
243	a_er = a_er .OR. (ier.NE.0)
244	ALLOCATE (co2_flux(kjpindex,nvm),stat=ier)
245	a_er = a_er .OR. (ier.NE.0)
246	IF ( a_er ) STOP 'PROBLEM WITH ALLOCATION'
247	!
248	! prepare forcing
249	!
250	max_totsize = 50
251	CALL getin ('STOMATE_FORCING_MEMSIZE',max_totsize)
252	max_totsize = max_totsize * 1000000
253	totsize_1step = SIZE(soiltype(:,3))*KIND(soiltype(:,3)) + &
254	SIZE(humrel_x)*KIND(humrel_x) + &
255	SIZE(litterhum)*KIND(litterhum) + &
256	SIZE(t2m)*KIND(t2m) + &
257	SIZE(t2m_min)*KIND(t2m_min) + &
258	SIZE(temp_sol)*KIND(temp_sol) + &
259	SIZE(soiltemp)*KIND(soiltemp) + &
260	SIZE(soilhum)*KIND(soilhum) + &
261	SIZE(precip_rain)*KIND(precip_rain) + &
262	SIZE(precip_snow)*KIND(precip_snow) + &
263	SIZE(gpp_x)*KIND(gpp_x) + &
264	SIZE(veget_force_x)*KIND(veget_force_x) + &
265	SIZE(veget_max_force_x)*KIND(veget_max_force_x) + &
266	SIZE(lai_force_x)*KIND(lai_force_x)
267	! total number of forcing steps
268	nsft = INT(one_year/(dt_force/one_day))
269	! number of forcing steps in memory
270	nsfm = MIN(nsft,MAX(1,NINT(FLOAT(max_totsize)/FLOAT(totsize_1step))))
271	!-
272	WRITE(numout,*) 'Offline forcing of Stomate:'
273	WRITE(numout,*) ' Total number of forcing steps:',nsft
274	WRITE(numout,*) ' Number of forcing steps in memory:',nsfm
275	!-
276	a_er = .FALSE.
277	ALLOCATE (clay_fm(kjpindex,nsfm),stat=ier)
278	a_er = a_er.OR.(ier.NE.0)
279	ALLOCATE (humrel_x_fm(kjpindex,nvm,nsfm),stat=ier)
280	a_er = a_er.OR.(ier.NE.0)
281	ALLOCATE (litterhum_fm(kjpindex,nsfm),stat=ier)
282	a_er = a_er.OR.(ier.NE.0)
283	ALLOCATE (t2m_fm(kjpindex,nsfm),stat=ier)
284	a_er = a_er.OR.(ier.NE.0)
285	ALLOCATE (t2m_min_fm(kjpindex,nsfm),stat=ier)
286	a_er = a_er.OR.(ier.NE.0)
287	ALLOCATE (temp_sol_fm(kjpindex,nsfm),stat=ier)
288	a_er = a_er.OR.(ier.NE.0)
289	ALLOCATE (soiltemp_fm(kjpindex,nbdl,nsfm),stat=ier)
290	a_er = a_er.OR.(ier.NE.0)
291	ALLOCATE (soilhum_fm(kjpindex,nbdl,nsfm),stat=ier)
292	a_er = a_er.OR.(ier.NE.0)
293	ALLOCATE (precip_fm(kjpindex,nsfm),stat=ier)
294	a_er = a_er.OR.(ier.NE.0)
295	ALLOCATE (gpp_x_fm(kjpindex,nvm,nsfm),stat=ier)
296	a_er = a_er.OR.(ier.NE.0)
297	ALLOCATE (veget_force_x_fm(kjpindex,nvm,nsfm),stat=ier)
298	a_er = a_er.OR.(ier.NE.0)
299	ALLOCATE (veget_max_force_x_fm(kjpindex,nvm,nsfm),stat=ier)
300	a_er = a_er.OR.(ier.NE.0)
301	ALLOCATE (lai_force_x_fm(kjpindex,nvm,nsfm),stat=ier)
302	a_er = a_er.OR.(ier.NE.0)
303	ALLOCATE (isf(nsfm),stat=ier)
304	a_er = a_er.OR.(ier.NE.0)
305	ALLOCATE (nf_written(nsft),stat=ier)
306	a_er = a_er.OR.(ier.NE.0)
307	ALLOCATE (nf_cumul(nsft),stat=ier)
308	a_er = a_er.OR.(ier.NE.0)
309	IF (a_er) THEN
310	STOP 'stomate: error in memory allocation for forcing data'
311	ENDIF
312	! ensure that we read all new forcing states
313	iisf = nsfm
314	! initialize that table that contains the indices
315	! of the forcing states that will be in memory
316	isf(:) = (/ (i,i=1,nsfm) /)
317	!-
318	! read info about grids
319	!-
320	contfrac(:) = 1.0
321	!-
322	ALLOCATE (x_indices(kjpindex),stat=ier)
323	ier = NF90_INQ_VARID (force_id,'index',v_id)
324	ier = NF90_GET_VAR (force_id,v_id,x_indices)
325	indices(:) = NINT(x_indices(:))
326	DEALLOCATE (x_indices)
327	!-
328	DO ji=1,kjpindex
329	DO jv=1,nvm
330	indexveg((jv-1)kjpindex+ji) = indices(ji)+(jv-1)kjpij
331	ENDDO
332	ENDDO
333	!-
334	ier = NF90_INQ_VARID (force_id,'lalo',v_id)
335	ier = NF90_GET_VAR (force_id,v_id,lalo)
336	!-
337	ALLOCATE (x_neighbours(kjpindex,8),stat=ier)
338	ier = NF90_INQ_VARID (force_id,'neighbours',v_id)
339	ier = NF90_GET_VAR (force_id,v_id,x_neighbours)
340	neighbours(:,:) = NINT(x_neighbours(:,:))
341	DEALLOCATE (x_neighbours)
342	!-
343	ier = NF90_INQ_VARID (force_id,'resolution',v_id)
344	ier = NF90_GET_VAR (force_id,v_id,resolution)
345	!-
346	ier = NF90_INQ_VARID (force_id,'contfrac',v_id)
347	ier = NF90_GET_VAR (force_id,v_id,contfrac)
348	!-
349	! activate CO2, STOMATE, but not sechiba
350	!-
351	control%river_routing = .FALSE.
352	control%hydrol_cwrr = .FALSE.
353	control%ok_sechiba = .FALSE.
354	!
355	control%stomate_watchout = .TRUE.
356	control%ok_co2 = .TRUE.
357	control%ok_stomate = .TRUE.
358	!!$DS now we search the values for the scalar parameters in orchidee.def
359	CALL getin_co2_parameters
360	CALL getin_stomate_parameters
361	!!$ CALL getin('CLAYFRACTION_DEFAULT',clayfraction_default)
362	!-
363	! is DGVM activated?
364	!-
365	control%ok_dgvm = .FALSE.
366	CALL getin('STOMATE_OK_DGVM',control%ok_dgvm)
367	WRITE(,) 'LPJ is activated: ',control%ok_dgvm
368	!!$DS for the externalisation reading the scalar parameters when ok_dgvm is set to true
369	IF (control%ok_dgvm) THEN
370	CALL getin_dgvm_parameters
371	ENDIF
372	!-
373	! restart files
374	!-
375	! Sechiba's restart files
376	sec_restname_in = 'sechiba_start.nc'
377	CALL getin('SECHIBA_restart_in',sec_restname_in)
378	WRITE(,) 'SECHIBA INPUT RESTART_FILE: ',TRIM(sec_restname_in)
379	IF ( TRIM(sec_restname_in) .EQ. 'NONE' ) THEN
380	STOP 'Need a restart file for Sechiba'
381	ENDIF
382	sec_restname_out = 'sechiba_restart.nc'
383	CALL getin('SECHIBA_rest_out',sec_restname_out)
384	WRITE(,) 'SECHIBA OUTPUT RESTART_FILE: ',TRIM(sec_restname_out)
385	! Stomate's restart files
386	sto_restname_in = 'stomate_start.nc'
387	CALL getin ('STOMATE_RESTART_FILEIN',sto_restname_in)
388	WRITE(,) 'STOMATE INPUT RESTART_FILE: ',TRIM(sto_restname_in)
389	sto_restname_out = 'stomate_restart.nc'
390	CALL getin ('STOMATE_RESTART_FILEOUT',sto_restname_out)
391	WRITE(,) 'STOMATE OUTPUT RESTART_FILE: ',TRIM(sto_restname_out)
392	!-
393	! We need to know iim, jjm.
394	! Get them from the restart files themselves.
395	!-
396	iret = NF90_OPEN (sec_restname_in,NF90_NOWRITE,ncfid)
397	IF (iret /= NF90_NOERR) THEN
398	CALL ipslerr (3,'teststomate', &
399	& 'Could not open file : ', &
400	& sec_restname_in,'(Do you have forget it ?)')
401	ENDIF
402	iret = NF90_INQUIRE_DIMENSION (ncfid,1,len=iim)
403	iret = NF90_INQUIRE_DIMENSION (ncfid,2,len=jjm)
404	iret = NF90_CLOSE (ncfid)
405	! Allocate longitudes and latitudes
406	ALLOCATE (lon(iim,jjm),stat=ier)
407	a_er = a_er.OR.(ier.NE.0)
408	ALLOCATE (lat(iim,jjm),stat=ier)
409	a_er = a_er.OR.(ier.NE.0)
410	lon(:,:) = 0.0
411	lat(:,:) = 0.0
412	lev(1) = 0.0
413	!-
414	CALL restini &
415	& (sec_restname_in, iim, jjm, lon, lat, llm, lev, &
416	& sec_restname_out, itau_dep, date0, dt_files, rest_id_sec)
417	!-
418	CALL restini &
419	& (sto_restname_in, iim, jjm, lon, lat, llm, lev, &
420	& sto_restname_out, itau_dep, date0, dt_files, rest_id_sto)
421	!-
422	IF ( dt_files .NE. dtradia ) THEN
423	WRITE(,) 'dt_files',dt_files
424	WRITE(,) 'dtradia',dtradia
425	STOP 'PROBLEM with time steps.'
426	ENDIF
427	!-
428	! a time step for STOMATE corresponds to itau_step timesteps in SECHIBA
429	!-
430	itau_step = NINT(dt_force/dtradia)
431	!
432	CALL ioconf_startdate(date0)
433	CALL intsurf_time( itau_dep, date0, dtradia )
434	!-
435	! Length of integration
436	!-
437	WRITE(time_str,'(a)') '1Y'
438	CALL getin ('TIME_LENGTH', time_str)
439	! transform into itau
440	CALL tlen2itau(time_str, dt_files, date0, itau_len)
441	! itau_len*dtradia must be a multiple of dt_force
442	itau_len = NINT( MAX(1.,FLOAT(NINT(itau_len*dtradia/dt_force))) &
443	& *dt_force/dtradia)
444	!-
445	! set up STOMATE history file
446	!-
447	!Config Key = STOMATE_OUTPUT_FILE
448	!Config Desc = Name of file in which STOMATE's output is going
449	!Config to be written
450	!Config Def = stomate_history.nc
451	!Config Help = This file is going to be created by the model
452	!Config and will contain the output from the model.
453	!Config This file is a truly COADS compliant netCDF file.
454	!Config It will be generated by the hist software from
455	!Config the IOIPSL package.
456	!-
457	stom_histname='stomate_history.nc'
458	CALL getin ('STOMATE_OUTPUT_FILE', stom_histname)
459	WRITE(,) 'STOMATE_OUTPUT_FILE', TRIM(stom_histname)
460	!-
461	!Config Key = STOMATE_HIST_DT
462	!Config Desc = STOMATE history time step (d)
463	!Config Def = 10.
464	!Config Help = Time step of the STOMATE history file
465	!-
466	hist_days_stom = 10.
467	CALL getin ('STOMATE_HIST_DT', hist_days_stom)
468	IF ( hist_days_stom == -1. ) THEN
469	hist_dt_stom = -1.
470	WRITE(numout,*) 'output frequency for STOMATE history file (d): one month.'
471	ELSE
472	hist_dt_stom = NINT( hist_days_stom ) * one_day
473	WRITE(numout,*) 'output frequency for STOMATE history file (d): ', &
474	hist_dt_stom/one_day
475	ENDIF
476	!-
477	! initialize
478	CALL histbeg(stom_histname, iim, lon, jjm, lat, 1, iim, 1, jjm, &
479	& itau_dep, date0, dt_files, hori_id, hist_id_sto)
480	! define PFT axis
481	hist_PFTaxis = (/ ( REAL(i,r_std), i=1,nvm ) /)
482	! declare this axis
483	CALL histvert (hist_id_sto, 'PFT', 'Plant functional type', &
484	& '-', nvm, hist_PFTaxis, hist_PFTaxis_id)
485	!- define Pool_10 axis
486	hist_pool_10axis = (/ ( REAL(i,r_std), i=1,10 ) /)
487	!- declare this axis
488	CALL histvert (hist_id_sto, 'P10', 'Pool 10 years', &
489	& '-', 10, hist_pool_10axis, hist_pool_10axis_id)
490	!- define Pool_100 axis
491	hist_pool_100axis = (/ ( REAL(i,r_std), i=1,100 ) /)
492	!- declare this axis
493	CALL histvert (hist_id_sto, 'P100', 'Pool 100 years', &
494	& '-', 100, hist_pool_100axis, hist_pool_100axis_id)
495	!- define Pool_11 axis
496	hist_pool_11axis = (/ ( REAL(i,r_std), i=1,11 ) /)
497	!- declare this axis
498	CALL histvert (hist_id_sto, 'P11', 'Pool 10 years + 1', &
499	& '-', 11, hist_pool_11axis, hist_pool_11axis_id)
500	!- define Pool_101 axis
501	hist_pool_101axis = (/ ( REAL(i,r_std), i=1,101 ) /)
502	!- declare this axis
503	CALL histvert (hist_id_sto, 'P101', 'Pool 100 years + 1', &
504	& '-', 101, hist_pool_101axis, hist_pool_101axis_id)
505	! define STOMATE history file
506	CALL stom_define_history (hist_id_sto, nvm, iim, jjm, &
507	& dt_files, hist_dt_stom, hori_id, hist_PFTaxis_id, &
508	& hist_pool_10axis_id, hist_pool_100axis_id, &
509	& hist_pool_11axis_id, hist_pool_101axis_id)
510	! end definition
511	CALL histend(hist_id_sto)
512	hist_id_sec = -1
513	hist_id_sec2 = -1
514	hist_id_stom_IPCC = -1
515	!-
516	! read some variables we need from SECHIBA's restart file
517	!-
518	CALL ioget_expval(val_exp)
519	!-
520	! first call of slowproc to initialize variables
521	!-
522	itau = itau_dep
523	ldrestart_read = .TRUE.
524	ldrestart_write = .FALSE.
525	!-
526	!MM Problem here with dpu which depends on soil type
527	DO jv = 1, nbdl-1
528	! first 2.0 is dpu
529	! second 2.0 is average
530	diaglev(jv) = 2.0/(2*(nbdl-1) -1) ( ( 2(jv-1) -1) + ( 2(jv) -1) ) / 2.0
531	ENDDO
532	diaglev(nbdl) = 2.0
533	!-
534	! For sequential use only, we must initialize data_para :
535	!
536	CALL init_para(.FALSE.)
537	!
538	CALL init_data_para(iim,jjm,kjpindex,indices)
539	!
540	!- global index index_g is the index_l of root proc
541	IF (is_root_prc) index_g(:)=indices(1:kjpindex)
542	!-
543
544	! DS : getin the value of slowproc parameters
545	CALL getin('CLAYFRACTION_DEFAULT',clayfraction_default)
546
547	CALL slowproc_main &
548	& (itau, kjpij, kjpindex, dt_force, date0, &
549	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
550	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
551	& t2m, t2m_min, temp_sol, soiltemp, &
552	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
553	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
554	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
555	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_sto, hist_id_stom_IPCC, co2_flux)
556	! correct date
557	day_counter = one_day - dt_force
558	date=1
559	!-
560	! time loop
561	!-
562	DO itau = itau_dep+itau_step,itau_dep+itau_len,itau_step
563	!-- next forcing state
564	iisf = iisf+1
565	!---
566	IF (iisf .GT. nsfm) THEN
567	!-----
568	!---- we have to read new forcing states from the file
569	!-----
570	!---- determine blocks of forcing states that are contiguous in memory
571	!-----
572	nblocks = 0
573	ifirst(:) = 1
574	ilast(:) = 1
575	!-----
576	DO iisf=1,nsfm
577	IF ( (nblocks .NE. 0) ) THEN
578	IF ( (isf(iisf) .EQ. isf(ilast(nblocks))+1) ) THEN
579	!-------- element is contiguous with last element found
580	ilast(nblocks) = iisf
581	ELSE
582	!-------- found first element of new block
583	nblocks = nblocks+1
584	IF (nblocks .GT. nsfm) THEN
585	! IF (nblocks .GT. 2) THEN
586	STOP 'Problem in teststomate'
587	ENDIF
588	ifirst(nblocks) = iisf
589	ilast(nblocks) = iisf
590	ENDIF
591	ELSE
592	!-------- found first element of new block
593	nblocks = nblocks+1
594	IF (nblocks .GT. nsfm) THEN
595	! IF (nblocks .GT. 2) THEN
596	STOP 'Problem in teststomate'
597	ENDIF
598	ifirst(nblocks) = iisf
599	ilast(nblocks) = iisf
600	ENDIF
601	ENDDO
602	!-----
603	DO iblocks=1,nblocks
604	IF (ifirst(iblocks) .NE. ilast(iblocks)) THEN
605	ndim = 2;
606	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
607	count(1:ndim) = SHAPE(clay_fm)
608	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
609	ier = NF90_INQ_VARID (force_id,'clay',v_id)
610	a_er = a_er.OR.(ier.NE.0)
611	ier = NF90_GET_VAR (force_id,v_id, &
612	& clay_fm(:,ifirst(iblocks):ilast(iblocks)), &
613	& start=start(1:ndim), count=count(1:ndim))
614	a_er = a_er.OR.(ier.NE.0)
615	!---------
616	ndim = 3;
617	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
618	count(1:ndim) = SHAPE(humrel_x_fm)
619	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
620	ier = NF90_INQ_VARID (force_id,'humrel',v_id)
621	a_er = a_er.OR.(ier.NE.0)
622	ier = NF90_GET_VAR (force_id,v_id, &
623	& humrel_x_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
624	& start=start(1:ndim), count=count(1:ndim))
625	a_er = a_er.OR.(ier.NE.0)
626	!---------
627	ndim = 2;
628	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
629	count(1:ndim) = SHAPE(litterhum_fm)
630	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
631	ier = NF90_INQ_VARID (force_id,'litterhum',v_id)
632	a_er = a_er.OR.(ier.NE.0)
633	ier = NF90_GET_VAR (force_id,v_id, &
634	& litterhum_fm(:,ifirst(iblocks):ilast(iblocks)), &
635	& start=start(1:ndim), count=count(1:ndim))
636	a_er = a_er.OR.(ier.NE.0)
637	!---------
638	ndim = 2;
639	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
640	count(1:ndim) = SHAPE(t2m_fm)
641	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
642	ier = NF90_INQ_VARID (force_id,'t2m',v_id)
643	a_er = a_er.OR.(ier.NE.0)
644	ier = NF90_GET_VAR (force_id,v_id, &
645	& t2m_fm(:,ifirst(iblocks):ilast(iblocks)), &
646	& start=start(1:ndim), count=count(1:ndim))
647	a_er = a_er.OR.(ier.NE.0)
648	!---------
649	ndim = 2;
650	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
651	count(1:ndim) = SHAPE(t2m_min_fm)
652	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
653	ier = NF90_INQ_VARID (force_id,'t2m_min',v_id)
654	a_er = a_er.OR.(ier.NE.0)
655	ier = NF90_GET_VAR (force_id,v_id, &
656	& t2m_min_fm(:,ifirst(iblocks):ilast(iblocks)), &
657	& start=start(1:ndim), count=count(1:ndim))
658	a_er = a_er.OR.(ier.NE.0)
659	!---------
660	ndim = 2;
661	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
662	count(1:ndim) = SHAPE(temp_sol_fm)
663	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
664	ier = NF90_INQ_VARID (force_id,'tsurf',v_id)
665	a_er = a_er.OR.(ier.NE.0)
666	ier = NF90_GET_VAR (force_id,v_id, &
667	& temp_sol_fm(:,ifirst(iblocks):ilast(iblocks)), &
668	& start=start(1:ndim), count=count(1:ndim))
669	a_er = a_er.OR.(ier.NE.0)
670	!---------
671	ndim = 3;
672	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
673	count(1:ndim) = SHAPE(soiltemp_fm)
674	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
675	ier = NF90_INQ_VARID (force_id,'tsoil',v_id)
676	a_er = a_er.OR.(ier.NE.0)
677	ier = NF90_GET_VAR (force_id,v_id, &
678	& soiltemp_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
679	& start=start(1:ndim), count=count(1:ndim))
680	a_er = a_er.OR.(ier.NE.0)
681	!---------
682	ndim = 3;
683	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
684	count(1:ndim) = SHAPE(soilhum_fm)
685	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
686	ier = NF90_INQ_VARID (force_id,'soilhum',v_id)
687	a_er = a_er.OR.(ier.NE.0)
688	ier = NF90_GET_VAR (force_id,v_id, &
689	& soilhum_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
690	& start=start(1:ndim), count=count(1:ndim))
691	a_er = a_er.OR.(ier.NE.0)
692	!---------
693	ndim = 2;
694	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
695	count(1:ndim) = SHAPE(precip_fm)
696	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
697	ier = NF90_INQ_VARID (force_id,'precip',v_id)
698	a_er = a_er.OR.(ier.NE.0)
699	ier = NF90_GET_VAR (force_id,v_id, &
700	& precip_fm(:,ifirst(iblocks):ilast(iblocks)), &
701	& start=start(1:ndim), count=count(1:ndim))
702	a_er = a_er.OR.(ier.NE.0)
703	!---------
704	ndim = 3;
705	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
706	count(1:ndim) = SHAPE(gpp_x_fm)
707	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
708	ier = NF90_INQ_VARID (force_id,'gpp',v_id)
709	a_er = a_er.OR.(ier.NE.0)
710	ier = NF90_GET_VAR (force_id,v_id, &
711	& gpp_x_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
712	& start=start(1:ndim), count=count(1:ndim))
713	a_er = a_er.OR.(ier.NE.0)
714	!---------
715	ndim = 3;
716	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
717	count(1:ndim) = SHAPE(veget_force_x_fm)
718	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
719	ier = NF90_INQ_VARID (force_id,'veget',v_id)
720	a_er = a_er.OR.(ier.NE.0)
721	ier = NF90_GET_VAR (force_id,v_id, &
722	& veget_force_x_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
723	& start=start(1:ndim), count=count(1:ndim))
724	a_er = a_er.OR.(ier.NE.0)
725	!---------
726	ndim = 3;
727	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
728	count(1:ndim) = SHAPE(veget_max_force_x_fm)
729	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
730	ier = NF90_INQ_VARID (force_id,'veget_max',v_id)
731	a_er = a_er.OR.(ier.NE.0)
732	ier = NF90_GET_VAR (force_id,v_id, &
733	& veget_max_force_x_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
734	& start=start(1:ndim), count=count(1:ndim))
735	a_er = a_er.OR.(ier.NE.0)
736	!---------
737	ndim = 3;
738	start(:) = 1; start(ndim) = isf(ifirst(iblocks));
739	count(1:ndim) = SHAPE(lai_force_x_fm)
740	count(ndim) = isf(ilast(iblocks))-isf(ifirst(iblocks))+1
741	ier = NF90_INQ_VARID (force_id,'lai',v_id)
742	a_er = a_er.OR.(ier.NE.0)
743	ier = NF90_GET_VAR (force_id,v_id, &
744	& lai_force_x_fm(:,:,ifirst(iblocks):ilast(iblocks)), &
745	& start=start(1:ndim), count=count(1:ndim))
746	a_er = a_er.OR.(ier.NE.0)
747	ENDIF
748	ENDDO
749	!-----
750	!---- determine which forcing states must be read next time
751	!-----
752	isf(1) = isf(nsfm)+1
753	IF ( isf(1) .GT. nsft ) isf(1) = 1
754	DO iisf = 2, nsfm
755	isf(iisf) = isf(iisf-1)+1
756	IF ( isf(iisf) .GT. nsft ) isf(iisf) = 1
757	ENDDO
758	!---- start again at first forcing state
759	iisf = 1
760	ENDIF
761	soiltype(:,3) = clay_fm(:,iisf)
762	humrel_x(:,:) = humrel_x_fm(:,:,iisf)
763	litterhum(:) = litterhum_fm(:,iisf)
764	t2m(:) = t2m_fm(:,iisf)
765	t2m_min(:) = t2m_min_fm(:,iisf)
766	temp_sol(:) = temp_sol_fm(:,iisf)
767	soiltemp(:,:) = soiltemp_fm(:,:,iisf)
768	soilhum(:,:) = soilhum_fm(:,:,iisf)
769	precip_rain(:) = precip_fm(:,iisf)
770	gpp_x(:,:) = gpp_x_fm(:,:,iisf)
771	veget_force_x(:,:) = veget_force_x_fm(:,:,iisf)
772	veget_max_force_x(:,:) = veget_max_force_x_fm(:,:,iisf)
773	lai_force_x(:,:) = lai_force_x_fm(:,:,iisf)
774	WHERE ( t2m(:) .LT. ZeroCelsius )
775	precip_snow(:) = precip_rain(:)
776	precip_rain(:) = 0.
777	ELSEWHERE
778	precip_snow(:) = 0.
779	ENDWHERE
780	!---
781	!-- scale GPP to new lai and veget_max
782	!---
783	WHERE ( lai_x(:,:) .EQ. 0.0 ) gpp_x(:,:) = 0.0
784	!-- scale GPP to new LAI
785	WHERE (lai_force_x(:,:) .GT. 0.0 )
786	gpp_x(:,:) = gpp_x(:,:)atan(2.lai_x(:,:)) &
787	& /atan( 2.*MAX(lai_force_x(:,:),0.01))
788	ENDWHERE
789	!-- scale GPP to new veget_max
790	WHERE (veget_max_force_x(:,:) .GT. 0.0 )
791	gpp_x(:,:) = gpp_x(:,:)*veget_max_x(:,:)/veget_max_force_x(:,:)
792	ENDWHERE
793	!---
794	!-- number crunching
795	!---
796	CALL intsurf_time( itau, date0, dtradia )
797	ldrestart_read = .FALSE.
798	ldrestart_write = .FALSE.
799	CALL slowproc_main &
800	& (itau, kjpij, kjpindex, dt_force, date0, &
801	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
802	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
803	& t2m, t2m_min, temp_sol, soiltemp, &
804	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
805	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
806	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
807	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_sto, hist_id_stom_IPCC, co2_flux)
808	day_counter = one_day - dt_force
809	ENDDO ! end of the time loop
810	!-
811	itau = itau -itau_step
812	!-
813	! write restart files
814	!-
815	! first, read and write variables that are not managed otherwise
816	taboo_vars = &
817	& '$lat$ $lon$ $lev$ $veget_year$ '// &
818	& '$height$ $day_counter$ $veget$ $veget_max$ $frac_nobio$ '// &
819	& '$lai$ $soiltype_frac$ $clay_frac$ '// &
820	& '$nav_lon$ $nav_lat$ $nav_lev$ $time$ $time_steps$'
821	!-
822	CALL ioget_vname(rest_id_sec, nbvar, varnames)
823	!!$!-
824	!!$! read and write some special variables (1D or variables that we need)
825	!!$!-
826	!!$ var_name = 'day_counter'
827	!!$ CALL restget (rest_id_sto, var_name, 1, 1, 1, itau_dep, .TRUE., xtmp)
828	!!$ CALL restput (rest_id_sto, var_name, 1, 1, 1, itau_dep, xtmp)
829	!!$!-
830	!!$ var_name = 'dt_days'
831	!!$ CALL restget (rest_id_sto, var_name, 1, 1, 1, itau_dep, .TRUE., xtmp)
832	!!$ CALL restput (rest_id_sto, var_name, 1, 1, 1, itau_dep, xtmp)
833	!!$!-
834	!!$ var_name = 'date'
835	!!$ CALL restget (rest_id_sto, var_name, 1, 1, 1, itau_dep, .TRUE., xtmp)
836	!!$ CALL restput (rest_id_sto, var_name, 1, 1, 1, itau_dep, xtmp)
837	!-
838	DO iv = 1, nbvar
839	!-- check if the variable is to be written here
840	IF (INDEX( taboo_vars,'$'//TRIM(varnames(iv))//'$') .EQ. 0 ) THEN
841	!---- get variable dimensions, especially 3rd dimension
842	CALL ioget_vdim &
843	& (rest_id_sec, varnames(iv), varnbdim_max, varnbdim, vardims)
844	l1d = ALL(vardims(1:varnbdim) .EQ. 1)
845	ALLOCATE(var_3d(kjpindex,vardims(3)),stat=ier)
846	IF (ier .NE. 0) STOP 'ALLOCATION PROBLEM'
847	!---- read it
848	IF (l1d) THEN
849	CALL restget (rest_id_sec,TRIM(varnames(iv)), &
850	1,vardims(3),1,itau_dep,.TRUE.,var_3d)
851	ELSE
852	CALL restget (rest_id_sec,TRIM(varnames(iv)), &
853	kjpindex,vardims(3),1,itau_dep,.TRUE.,var_3d, &
854	"gather",kjpindex,indices)
855	ENDIF
856	!---- write it
857	IF (l1d) THEN
858	CALL restput (rest_id_sec,TRIM(varnames(iv)), &
859	1,vardims(3),1,itau,var_3d)
860	ELSE
861	CALL restput (rest_id_sec,TRIM(varnames(iv)), &
862	kjpindex,vardims(3),1,itau,var_3d, &
863	'scatter',kjpindex,indices)
864	ENDIF
865	DEALLOCATE (var_3d)
866	ENDIF
867	ENDDO
868	! call slowproc to write restart files
869	ldrestart_read = .FALSE.
870	ldrestart_write = .TRUE.
871	!-
872	CALL slowproc_main &
873	& (itau, kjpij, kjpindex, dt_force, date0, &
874	& ldrestart_read, ldrestart_write, .FALSE., .TRUE., &
875	& indices, indexveg, lalo, neighbours, resolution, contfrac, soiltype, &
876	& t2m, t2m_min, temp_sol, soiltemp, &
877	& humrel_x, soilhum, litterhum, precip_rain, precip_snow, gpp_x, &
878	& deadleaf_cover, assim_param_x, lai_x, height_x, veget_x, &
879	& frac_nobio, veget_max_x, totfrac_nobio, qsintmax_x, &
880	& rest_id_sec, hist_id_sec, hist_id_sec2, rest_id_sto, hist_id_sto, hist_id_stom_IPCC, co2_flux)
881	!-
882	! close files
883	!-
884	CALL restclo
885	CALL histclo
886	ier = NF90_CLOSE (force_id)
887	!-
888	! write a new driver restart file with correct time step
889	!-
890	write(,) 'teststomate: writing driver restart file with correct time step.'
891	dri_restname_in = 'driver_start.nc'
892	CALL getin ('RESTART_FILEIN',dri_restname_in)
893	dri_restname_out = 'driver_restart.nc'
894	CALL getin ('RESTART_FILEOUT',dri_restname_out)
895	CALL SYSTEM &
896	& ('cp '//TRIM(dri_restname_in)//' '//TRIM(dri_restname_out))
897	!-
898	iret = NF90_OPEN (TRIM(sec_restname_out),NF90_NOWRITE,ncfid)
899	iret = NF90_INQ_VARID (ncfid,'time',v_id)
900	iret = NF90_GET_VAR (ncfid,v_id,r1d)
901	time_sec = r1d(1)
902	iret = NF90_INQ_VARID (ncfid,'time_steps',v_id)
903	iret = NF90_GET_VAR (ncfid,v_id,time_step_sec)
904	iret = NF90_CLOSE (ncfid)
905	!-
906	iret = NF90_OPEN (TRIM(dri_restname_out),NF90_WRITE,ncfid)
907	iret = NF90_INQ_VARID (ncfid,'time',v_id)
908	iret = NF90_GET_VAR (ncfid,v_id,r1d)
909	time_dri = r1d(1)
910	r1d(1) = time_sec
911	iret = NF90_PUT_VAR (ncfid,v_id,r1d)
912	iret = NF90_INQ_VARID (ncfid,'time_steps',v_id)
913	iret = NF90_GET_VAR (ncfid,v_id,time_step_dri)
914	iret = NF90_PUT_VAR (ncfid,v_id,time_step_sec)
915	iret = NF90_INQ_VARID (ncfid,'julian',v_id)
916	iret = NF90_GET_VAR (ncfid,v_id,r1d)
917	julian = r1d(1)
918	djulian = (time_step_sec-time_step_dri)*dtradia/one_day
919	julian = julian &
920	& +djulian-FLOAT(INT((julian+djulian)/one_year))*one_year
921	r1d(1) = julian
922	iret = NF90_PUT_VAR (ncfid,v_id,r1d)
923	iret = NF90_CLOSE (ncfid)
924
925	CALL getin_dump
926	!---------------
927	END PROGRAM teststomate
928	!
929	!===
930	!

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

Download in other formats: