Context Navigation

mod_1D_amma_read.F90 @ 3066

Last change on this file since 3066 was 2373, checked in by jyg, 9 years ago
From MPL: Correction of some missing initializations in dyn1D and of two 1D case bugs.
File size: 15.0 KB

Line
1	MODULE mod_1D_amma_read
2
3	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
4	!Declarations specifiques au cas AMMA
5	character*80 :: fich_amma
6	! Option du cas AMMA ou on impose la discretisation verticale (Ap,Bp)
7	integer nlev_amma, nt_amma
8
9
10	integer year_ini_amma, day_ini_amma, mth_ini_amma
11	real heure_ini_amma
12	real day_ju_ini_amma ! Julian day of amma first day
13	parameter (year_ini_amma=2006)
14	parameter (mth_ini_amma=7)
15	parameter (day_ini_amma=10) ! 10 = 10Juil2006
16	parameter (heure_ini_amma=0.) !0h en secondes
17	real dt_amma
18	parameter (dt_amma=1800.)
19
20	!profils initiaux:
21	real, allocatable:: plev_amma(:)
22
23	real, allocatable:: z_amma(:)
24	real, allocatable:: th_amma(:),q_amma(:)
25	real, allocatable:: u_amma(:)
26	real, allocatable:: v_amma(:)
27
28	real, allocatable:: th_ammai(:),q_ammai(:)
29	real, allocatable:: u_ammai(:)
30	real, allocatable:: v_ammai(:)
31	real, allocatable:: vitw_ammai(:)
32	real, allocatable:: ht_ammai(:)
33	real, allocatable:: hq_ammai(:)
34	real, allocatable:: vt_ammai(:)
35	real, allocatable:: vq_ammai(:)
36
37	!forcings
38	real, allocatable:: ht_amma(:,:)
39	real, allocatable:: hq_amma(:,:)
40	real, allocatable:: vitw_amma(:,:)
41	real, allocatable:: lat_amma(:),sens_amma(:)
42
43	!champs interpoles
44	real, allocatable:: vitw_profamma(:)
45	real, allocatable:: ht_profamma(:)
46	real, allocatable:: hq_profamma(:)
47	real lat_profamma,sens_profamma
48	real, allocatable:: vt_profamma(:)
49	real, allocatable:: vq_profamma(:)
50	real, allocatable:: th_profamma(:)
51	real, allocatable:: q_profamma(:)
52	real, allocatable:: u_profamma(:)
53	real, allocatable:: v_profamma(:)
54
55
56	CONTAINS
57
58	SUBROUTINE read_1D_cases
59	implicit none
60
61	#include "netcdf.inc"
62
63	INTEGER nid,rid,ierr
64
65	fich_amma='amma.nc'
66	print*,'fich_amma ',fich_amma
67	ierr = NF_OPEN(fich_amma,NF_NOWRITE,nid)
68	print*,'fich_amma,NF_NOWRITE,nid ',fich_amma,NF_NOWRITE,nid
69	if (ierr.NE.NF_NOERR) then
70	write(,) 'ERROR: GROS Pb opening forcings nc file '
71	write(,) NF_STRERROR(ierr)
72	stop ""
73	endif
74	!.......................................................................
75	ierr=NF_INQ_DIMID(nid,'lev',rid)
76	IF (ierr.NE.NF_NOERR) THEN
77	print*, 'Oh probleme lecture dimension zz'
78	ENDIF
79	ierr=NF_INQ_DIMLEN(nid,rid,nlev_amma)
80	print*,'OK nid,rid,nlev_amma',nid,rid,nlev_amma
81	!.......................................................................
82	ierr=NF_INQ_DIMID(nid,'time',rid)
83	print*,'nid,rid',nid,rid
84	nt_amma=0
85	IF (ierr.NE.NF_NOERR) THEN
86	stop 'probleme lecture dimension sens'
87	ENDIF
88	ierr=NF_INQ_DIMLEN(nid,rid,nt_amma)
89	print*,'nid,rid,nlev_amma',nid,rid,nt_amma
90
91	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
92	!profils initiaux:
93	allocate(plev_amma(nlev_amma))
94
95	allocate(z_amma(nlev_amma))
96	allocate(th_amma(nlev_amma),q_amma(nlev_amma))
97	allocate(u_amma(nlev_amma))
98	allocate(v_amma(nlev_amma))
99
100	!forcings
101	allocate(ht_amma(nlev_amma,nt_amma))
102	allocate(hq_amma(nlev_amma,nt_amma))
103	allocate(vitw_amma(nlev_amma,nt_amma))
104	allocate(lat_amma(nt_amma),sens_amma(nt_amma))
105
106	!profils initiaux:
107	allocate(th_ammai(nlev_amma),q_ammai(nlev_amma))
108	allocate(u_ammai(nlev_amma))
109	allocate(v_ammai(nlev_amma))
110	allocate(vitw_ammai(nlev_amma) )
111	allocate(ht_ammai(nlev_amma))
112	allocate(hq_ammai(nlev_amma))
113	allocate(vt_ammai(nlev_amma))
114	allocate(vq_ammai(nlev_amma))
115
116	!champs interpoles
117	allocate(vitw_profamma(nlev_amma))
118	allocate(ht_profamma(nlev_amma))
119	allocate(hq_profamma(nlev_amma))
120	allocate(vt_profamma(nlev_amma))
121	allocate(vq_profamma(nlev_amma))
122	allocate(th_profamma(nlev_amma))
123	allocate(q_profamma(nlev_amma))
124	allocate(u_profamma(nlev_amma))
125	allocate(v_profamma(nlev_amma))
126
127	print*,'Allocations OK'
128	call read_amma(nid,nlev_amma,nt_amma &
129	& ,z_amma,plev_amma,th_amma,q_amma,u_amma,v_amma,vitw_amma &
130	& ,ht_amma,hq_amma,sens_amma,lat_amma)
131
132	END SUBROUTINE read_1D_cases
133
134
135
136	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
137	SUBROUTINE deallocate_1D_cases
138	!profils initiaux:
139	deallocate(plev_amma)
140
141	deallocate(z_amma)
142	deallocate(th_amma,q_amma)
143	deallocate(u_amma)
144	deallocate(v_amma)
145
146	deallocate(th_ammai,q_ammai)
147	deallocate(u_ammai)
148	deallocate(v_ammai)
149	deallocate(vitw_ammai )
150	deallocate(ht_ammai)
151	deallocate(hq_ammai)
152	deallocate(vt_ammai)
153	deallocate(vq_ammai)
154
155	!forcings
156	deallocate(ht_amma)
157	deallocate(hq_amma)
158	deallocate(vitw_amma)
159	deallocate(lat_amma,sens_amma)
160
161	!champs interpoles
162	deallocate(vitw_profamma)
163	deallocate(ht_profamma)
164	deallocate(hq_profamma)
165	deallocate(vt_profamma)
166	deallocate(vq_profamma)
167	deallocate(th_profamma)
168	deallocate(q_profamma)
169	deallocate(u_profamma)
170	deallocate(v_profamma)
171	END SUBROUTINE deallocate_1D_cases
172
173
174	END MODULE mod_1D_amma_read
175	!=====================================================================
176	subroutine read_amma(nid,nlevel,ntime &
177	& ,zz,pp,temp,qv,u,v,dw &
178	& ,dt,dq,sens,flat)
179
180	!program reading forcings of the AMMA case study
181	implicit none
182	#include "netcdf.inc"
183
184	integer ntime,nlevel
185
186	real zz(nlevel)
187	real temp(nlevel),pp(nlevel)
188	real qv(nlevel),u(nlevel)
189	real v(nlevel)
190	real dw(nlevel,ntime)
191	real dt(nlevel,ntime)
192	real dq(nlevel,ntime)
193	real flat(ntime),sens(ntime)
194
195
196	integer nid, ierr,rid
197	integer nbvar3d
198	parameter(nbvar3d=30)
199	integer var3didin(nbvar3d)
200
201	ierr=NF_INQ_VARID(nid,"zz",var3didin(1))
202	if(ierr/=NF_NOERR) then
203	write(,) NF_STRERROR(ierr)
204	stop 'lev'
205	endif
206
207
208	ierr=NF_INQ_VARID(nid,"temp",var3didin(2))
209	if(ierr/=NF_NOERR) then
210	write(,) NF_STRERROR(ierr)
211	stop 'temp'
212	endif
213
214	ierr=NF_INQ_VARID(nid,"qv",var3didin(3))
215	if(ierr/=NF_NOERR) then
216	write(,) NF_STRERROR(ierr)
217	stop 'qv'
218	endif
219
220	ierr=NF_INQ_VARID(nid,"u",var3didin(4))
221	if(ierr/=NF_NOERR) then
222	write(,) NF_STRERROR(ierr)
223	stop 'u'
224	endif
225
226	ierr=NF_INQ_VARID(nid,"v",var3didin(5))
227	if(ierr/=NF_NOERR) then
228	write(,) NF_STRERROR(ierr)
229	stop 'v'
230	endif
231
232	ierr=NF_INQ_VARID(nid,"dw",var3didin(6))
233	if(ierr/=NF_NOERR) then
234	write(,) NF_STRERROR(ierr)
235	stop 'dw'
236	endif
237
238	ierr=NF_INQ_VARID(nid,"dt",var3didin(7))
239	if(ierr/=NF_NOERR) then
240	write(,) NF_STRERROR(ierr)
241	stop 'dt'
242	endif
243
244	ierr=NF_INQ_VARID(nid,"dq",var3didin(8))
245	if(ierr/=NF_NOERR) then
246	write(,) NF_STRERROR(ierr)
247	stop 'dq'
248	endif
249
250	ierr=NF_INQ_VARID(nid,"sens",var3didin(9))
251	if(ierr/=NF_NOERR) then
252	write(,) NF_STRERROR(ierr)
253	stop 'sens'
254	endif
255
256	ierr=NF_INQ_VARID(nid,"flat",var3didin(10))
257	if(ierr/=NF_NOERR) then
258	write(,) NF_STRERROR(ierr)
259	stop 'flat'
260	endif
261
262	ierr=NF_INQ_VARID(nid,"pp",var3didin(11))
263	if(ierr/=NF_NOERR) then
264	write(,) NF_STRERROR(ierr)
265	endif
266
267	!dimensions lecture
268	! call catchaxis(nid,ntime,nlevel,time,z,ierr)
269
270	#ifdef NC_DOUBLE
271	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(1),zz)
272	#else
273	ierr = NF_GET_VAR_REAL(nid,var3didin(1),zz)
274	#endif
275	if(ierr/=NF_NOERR) then
276	write(,) NF_STRERROR(ierr)
277	stop "getvarup"
278	endif
279	! write(,)'lecture z ok',zz
280
281	#ifdef NC_DOUBLE
282	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(2),temp)
283	#else
284	ierr = NF_GET_VAR_REAL(nid,var3didin(2),temp)
285	#endif
286	if(ierr/=NF_NOERR) then
287	write(,) NF_STRERROR(ierr)
288	stop "getvarup"
289	endif
290	! write(,)'lecture th ok',temp
291
292	#ifdef NC_DOUBLE
293	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(3),qv)
294	#else
295	ierr = NF_GET_VAR_REAL(nid,var3didin(3),qv)
296	#endif
297	if(ierr/=NF_NOERR) then
298	write(,) NF_STRERROR(ierr)
299	stop "getvarup"
300	endif
301	! write(,)'lecture qv ok',qv
302
303	#ifdef NC_DOUBLE
304	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(4),u)
305	#else
306	ierr = NF_GET_VAR_REAL(nid,var3didin(4),u)
307	#endif
308	if(ierr/=NF_NOERR) then
309	write(,) NF_STRERROR(ierr)
310	stop "getvarup"
311	endif
312	! write(,)'lecture u ok',u
313
314	#ifdef NC_DOUBLE
315	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(5),v)
316	#else
317	ierr = NF_GET_VAR_REAL(nid,var3didin(5),v)
318	#endif
319	if(ierr/=NF_NOERR) then
320	write(,) NF_STRERROR(ierr)
321	stop "getvarup"
322	endif
323	! write(,)'lecture v ok',v
324
325	#ifdef NC_DOUBLE
326	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(6),dw)
327	#else
328	ierr = NF_GET_VAR_REAL(nid,var3didin(6),dw)
329	#endif
330	if(ierr/=NF_NOERR) then
331	write(,) NF_STRERROR(ierr)
332	stop "getvarup"
333	endif
334	! write(,)'lecture w ok',dw
335
336	#ifdef NC_DOUBLE
337	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(7),dt)
338	#else
339	ierr = NF_GET_VAR_REAL(nid,var3didin(7),dt)
340	#endif
341	if(ierr/=NF_NOERR) then
342	write(,) NF_STRERROR(ierr)
343	stop "getvarup"
344	endif
345	! write(,)'lecture dt ok',dt
346
347	#ifdef NC_DOUBLE
348	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(8),dq)
349	#else
350	ierr = NF_GET_VAR_REAL(nid,var3didin(8),dq)
351	#endif
352	if(ierr/=NF_NOERR) then
353	write(,) NF_STRERROR(ierr)
354	stop "getvarup"
355	endif
356	! write(,)'lecture dq ok',dq
357
358	#ifdef NC_DOUBLE
359	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(9),sens)
360	#else
361	ierr = NF_GET_VAR_REAL(nid,var3didin(9),sens)
362	#endif
363	if(ierr/=NF_NOERR) then
364	write(,) NF_STRERROR(ierr)
365	stop "getvarup"
366	endif
367	! write(,)'lecture sens ok',sens
368
369	#ifdef NC_DOUBLE
370	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(10),flat)
371	#else
372	ierr = NF_GET_VAR_REAL(nid,var3didin(10),flat)
373	#endif
374	if(ierr/=NF_NOERR) then
375	write(,) NF_STRERROR(ierr)
376	stop "getvarup"
377	endif
378	! write(,)'lecture flat ok',flat
379
380	#ifdef NC_DOUBLE
381	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(11),pp)
382	#else
383	ierr = NF_GET_VAR_REAL(nid,var3didin(11),pp)
384	#endif
385	if(ierr/=NF_NOERR) then
386	write(,) NF_STRERROR(ierr)
387	stop "getvarup"
388	endif
389	! write(,)'lecture pp ok',pp
390
391	return
392	end subroutine read_amma
393	!======================================================================
394	SUBROUTINE interp_amma_time(day,day1,annee_ref &
395	& ,year_ini_amma,day_ini_amma,nt_amma,dt_amma,nlev_amma &
396	& ,vitw_amma,ht_amma,hq_amma,lat_amma,sens_amma &
397	& ,vitw_prof,ht_prof,hq_prof,lat_prof,sens_prof)
398	implicit none
399
400	!---------------------------------------------------------------------------------------
401	! Time interpolation of a 2D field to the timestep corresponding to day
402	!
403	! day: current julian day (e.g. 717538.2)
404	! day1: first day of the simulation
405	! nt_amma: total nb of data in the forcing (e.g. 48 for AMMA)
406	! dt_amma: total time interval (in sec) between 2 forcing data (e.g. 30min for AMMA)
407	!---------------------------------------------------------------------------------------
408
409	#include "compar1d.h"
410
411	! inputs:
412	integer annee_ref
413	integer nt_amma,nlev_amma
414	integer year_ini_amma
415	real day, day1,day_ini_amma,dt_amma
416	real vitw_amma(nlev_amma,nt_amma)
417	real ht_amma(nlev_amma,nt_amma)
418	real hq_amma(nlev_amma,nt_amma)
419	real lat_amma(nt_amma)
420	real sens_amma(nt_amma)
421	! outputs:
422	real vitw_prof(nlev_amma)
423	real ht_prof(nlev_amma)
424	real hq_prof(nlev_amma)
425	real lat_prof,sens_prof
426	! local:
427	integer it_amma1, it_amma2,k
428	real timeit,time_amma1,time_amma2,frac
429
430
431	if (forcing_type.eq.6) then
432	! Check that initial day of the simulation consistent with AMMA case:
433	if (annee_ref.ne.2006) then
434	print*,'Pour AMMA, annee_ref doit etre 2006'
435	print*,'Changer annee_ref dans run.def'
436	stop
437	endif
438	if (annee_ref.eq.2006 .and. day1.lt.day_ini_amma) then
439	print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_amma
440	print*,'Changer dayref dans run.def'
441	stop
442	endif
443	if (annee_ref.eq.2006 .and. day1.gt.day_ini_amma+1) then
444	print*,'AMMA a fini le 11 juillet'
445	print*,'Changer dayref ou nday dans run.def'
446	stop
447	endif
448	endif
449
450	! Determine timestep relative to the 1st day of AMMA:
451	! timeit=(day-day1)*86400.
452	! if (annee_ref.eq.1992) then
453	! timeit=(day-day_ini_toga)*86400.
454	! else
455	! timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
456	! endif
457	timeit=(day-day_ini_amma)*86400
458
459	! Determine the closest observation times:
460	! it_amma1=INT(timeit/dt_amma)+1
461	! it_amma2=it_amma1 + 1
462	! time_amma1=(it_amma1-1)*dt_amma
463	! time_amma2=(it_amma2-1)*dt_amma
464
465	it_amma1=INT(timeit/dt_amma)+1
466	IF (it_amma1 .EQ. nt_amma) THEN
467	it_amma2=it_amma1
468	ELSE
469	it_amma2=it_amma1 + 1
470	ENDIF
471	time_amma1=(it_amma1-1)*dt_amma
472	time_amma2=(it_amma2-1)*dt_amma
473
474	if (it_amma1 .gt. nt_amma) then
475	write(,) 'PB-stop: day, it_amma1, it_amma2, timeit: ' &
476	& ,day,day_ini_amma,it_amma1,it_amma2,timeit/86400.
477	stop
478	endif
479
480	! time interpolation:
481	IF (it_amma1 .EQ. it_amma2) THEN
482	frac=0.
483	ELSE
484	frac=(time_amma2-timeit)/(time_amma2-time_amma1)
485	frac=max(frac,0.0)
486	ENDIF
487
488	lat_prof = lat_amma(it_amma2) &
489	& -frac*(lat_amma(it_amma2)-lat_amma(it_amma1))
490	sens_prof = sens_amma(it_amma2) &
491	& -frac*(sens_amma(it_amma2)-sens_amma(it_amma1))
492
493	do k=1,nlev_amma
494	vitw_prof(k) = vitw_amma(k,it_amma2) &
495	& -frac*(vitw_amma(k,it_amma2)-vitw_amma(k,it_amma1))
496	ht_prof(k) = ht_amma(k,it_amma2) &
497	& -frac*(ht_amma(k,it_amma2)-ht_amma(k,it_amma1))
498	hq_prof(k) = hq_amma(k,it_amma2) &
499	& -frac*(hq_amma(k,it_amma2)-hq_amma(k,it_amma1))
500	enddo
501
502	return
503	END
504

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

Download in other formats:

Original Format