Context Navigation

mod_1D_cases_read_std.F90 @ 3680

Last change on this file since 3680 was 3680, checked in by fhourdin, 4 years ago
Modifications pour le format standard. Frédéric pour Marie-Pierre Frédéric
File size: 55.2 KB

Line
1	!
2	! $Id: mod_1D_cases_read.F90 2373 2015-10-13 17:28:01Z jyg $
3	!
4	MODULE mod_1D_cases_read_std
5
6	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
7	!Declarations specifiques au cas standard
8	character*80 :: fich_cas
9	! Discr?tisation
10	integer nlev_cas, nt_cas
11	real zzs_cas,pp_cas
12
13
14	!profils environnementaux
15	real, allocatable:: ppforc_cas(:,:),plev_cas(:,:)
16
17	!profils initiaux
18	real, allocatable:: zzforc_cas(:,:)
19	real, allocatable:: qt0_cas(:),qv0_cas(:),ql0_cas(:),qi0_cas(:),tke_cas(:)
20	real, allocatable:: rt0_cas(:),rv0_cas(:),rl0_cas(:),ri0_cas(:),rh0_cas(:)
21	real, allocatable:: temp0_cas(:),theta0_cas(:), thetal0_cas(:)
22	real, allocatable:: u0_cas(:),v0_cas(:),w_cas(:,:),omega_cas(:,:),ug_cas(:,:), vg_cas(:,:)
23	real, allocatable:: t_cas(:),theta_cas(:), thl_cas(:),u_cas(:),v_cas(:)
24	!advections et nudging
25	real, allocatable:: uadv_cas(:,:),vadv_cas(:,:)
26	real, allocatable:: tadv_cas(:,:),thadv_cas(:,:),thladv_cas(:,:)
27	real, allocatable:: qtadv_cas(:,:),qvadv_cas(:,:)
28	real, allocatable:: rtadv_cas(:,:),rvadv_cas(:,:)
29	real, allocatable:: trad_cas(:,:),thrad_cas(:,:),thlrad_cas(:,:)
30	real, allocatable:: temp_nudg_cas(:,:),th_nudg_cas(:,:),thl_nudg_cas(:,:)
31	real, allocatable:: qv_nudg_cas(:,:),qt_nudg_cas(:,:)
32	real, allocatable:: rv_nudg_cas(:,:),rt_nudg_cas(:,:)
33	real, allocatable:: u_nudg_cas(:,:),v_nudg_cas(:,:)
34	! flux
35	real, allocatable:: lat_cas(:),sens_cas(:),ustar_cas(:)
36	real, allocatable:: ts_cas(:),ps_cas(:),ps_forc_cas(:)
37	real, allocatable:: wpthetap_cas(:),wpqvp_cas(:),wpqtp_cas(:),wprvp_cas(:),wprtp_cas(:)
38
39	!champs interpoles
40	real, allocatable:: plev_prof_cas(:)
41	real, allocatable:: plev_forc_prof_cas(:)
42	real, allocatable:: pforc_prof_cas(:)
43	real, allocatable:: t_prof_cas(:),th_prof_cas(:),thl_prof_cas(:)
44	real, allocatable:: qt_prof_cas(:),qv_prof_cas(:),ql_prof_cas(:),qi_prof_cas(:)
45	real, allocatable:: rh_prof_cas(:)
46	real, allocatable:: rt_prof_cas(:),rv_prof_cas(:),rl_prof_cas(:),ri_prof_cas(:)
47	real, allocatable:: u_prof_cas(:),v_prof_cas(:),w_prof_cas(:),omega_prof_cas(:)
48	real, allocatable:: ug_prof_cas(:),vg_prof_cas(:)
49	real, allocatable:: uadv_prof_cas(:),vadv_prof_cas(:),tadv_prof_cas(:),thadv_prof_cas(:),thladv_prof_cas(:)
50	real, allocatable:: qtadv_prof_cas(:),qvadv_prof_cas(:),rtadv_prof_cas(:),rvadv_prof_cas(:)
51	real, allocatable:: temp_nudg_prof_cas(:), th_nudg_prof_cas(:), thl_nudg_prof_cas(:)
52	real, allocatable:: qv_nudg_prof_cas(:), qt_nudg_prof_cas(:), rv_nudg_prof_cas(:), rt_nudg_prof_cas(:)
53	real, allocatable:: u_nudg_prof_cas(:),v_nudg_prof_cas(:)
54	real, allocatable:: trad_prof_cas(:),thrad_prof_cas(:),thlrad_prof_cas(:)
55
56	real lat_prof_cas,sens_prof_cas,ts_prof_cas,ps_prof_cas,ps_forc_prof_cas,ustar_prof_cas,tke_prof_cas
57	real wpqtp_prof_cas,wpqvp_prof_cas,wprtp_prof_cas,wprvp_prof_cas,wpthetap_prof_cas
58	! real o3_cas,orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,heat_rough,rugos_cas,sand_cas,clay_cas
59
60
61
62	CONTAINS
63
64
65	!**********************************************************************************************
66	SUBROUTINE read_SCM_cas
67	implicit none
68
69	#include "netcdf.inc"
70	#include "date_cas.h"
71
72	INTEGER nid,rid,ierr
73	INTEGER ii,jj,timeid
74	REAL, ALLOCATABLE :: time_val(:)
75
76	print*,'ON EST VRAIMENT LA'
77	fich_cas='cas.nc'
78	print*,'fich_cas ',fich_cas
79	ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid)
80	print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid
81	if (ierr.NE.NF_NOERR) then
82	write(,) 'ERROR: GROS Pb opening forcings nc file '
83	write(,) NF_STRERROR(ierr)
84	stop ""
85	endif
86	!.......................................................................
87	ierr=NF_INQ_DIMID(nid,'lat',rid)
88	IF (ierr.NE.NF_NOERR) THEN
89	print*, 'Oh probleme lecture dimension lat'
90	ENDIF
91	ierr=NF_INQ_DIMLEN(nid,rid,ii)
92	print*,'OK1 read2: nid,rid,lat',nid,rid,ii
93	!.......................................................................
94	ierr=NF_INQ_DIMID(nid,'lon',rid)
95	IF (ierr.NE.NF_NOERR) THEN
96	print*, 'Oh probleme lecture dimension lon'
97	ENDIF
98	ierr=NF_INQ_DIMLEN(nid,rid,jj)
99	print*,'OK2 read2: nid,rid,lat',nid,rid,jj
100	!.......................................................................
101	ierr=NF_INQ_DIMID(nid,'lev',rid)
102	IF (ierr.NE.NF_NOERR) THEN
103	print*, 'Oh probleme lecture dimension nlev'
104	ENDIF
105	ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas)
106	print*,'OK3 read2: nid,rid,nlev_cas',nid,rid,nlev_cas
107	IF ( .NOT. ( nlev_cas > 10 .AND. nlev_cas < 200000 )) THEN
108	print*,'Valeur de nlev_cas peu probable'
109	STOP
110	ENDIF
111	!.......................................................................
112	ierr=NF_INQ_DIMID(nid,'time',rid)
113	nt_cas=0
114	IF (ierr.NE.NF_NOERR) THEN
115	stop 'Oh probleme lecture dimension time'
116	ENDIF
117	ierr=NF_INQ_DIMLEN(nid,rid,nt_cas)
118	print*,'OK4 read2: nid,rid,nt_cas',nid,rid,nt_cas
119	! Lecture de l'axe des temps
120	print*,'LECTURE DU TEMPS'
121	ierr=NF_INQ_VARID(nid,'time',timeid)
122	if(ierr/=NF_NOERR) then
123	print *,'Variable time manquante dans cas.nc:'
124	ierr=NF_NOERR
125	else
126	allocate(time_val(nt_cas))
127	#ifdef NC_DOUBLE
128	ierr = NF_GET_VAR_DOUBLE(nid,timeid,time_val)
129	#else
130	ierr = NF_GET_VAR_REAL(nid,timeid,time_val)
131	#endif
132	if(ierr/=NF_NOERR) then
133	print *,'Pb a la lecture de time cas.nc: '
134	endif
135	endif
136	IF (nt_cas>1) THEN
137	pdt_cas=time_val(2)-time_val(1)
138	ELSE
139	pdt_cas=0.
140	ENDIF
141
142	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
143	allocate(zzforc_cas(nlev_cas,nt_cas))
144	allocate(ppforc_cas(nlev_cas,nt_cas))
145	!profils initiaux
146	allocate(temp0_cas(nlev_cas),theta0_cas(nlev_cas),thetal0_cas(nlev_cas),tke_cas(nlev_cas))
147	allocate(qt0_cas(nlev_cas),qv0_cas(nlev_cas),ql0_cas(nlev_cas),qi0_cas(nlev_cas),u0_cas(nlev_cas),v0_cas(nlev_cas))
148	allocate(rt0_cas(nlev_cas),rv0_cas(nlev_cas),rl0_cas(nlev_cas),ri0_cas(nlev_cas),rh0_cas(nlev_cas))
149	allocate(t_cas(nlev_cas),theta_cas(nlev_cas),thl_cas(nlev_cas),u_cas(nlev_cas),v_cas(nlev_cas))
150	allocate(w_cas(nlev_cas,nt_cas),omega_cas(nlev_cas,nt_cas))
151	allocate(ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas))
152	!advections et nudging
153	allocate(uadv_cas(nlev_cas,nt_cas),vadv_cas(nlev_cas,nt_cas))
154	allocate(tadv_cas(nlev_cas,nt_cas),thadv_cas(nlev_cas,nt_cas),thladv_cas(nlev_cas,nt_cas))
155	allocate(qtadv_cas(nlev_cas,nt_cas),qvadv_cas(nlev_cas,nt_cas))
156	allocate(rtadv_cas(nlev_cas,nt_cas),rvadv_cas(nlev_cas,nt_cas))
157	allocate(trad_cas(nlev_cas,nt_cas),thrad_cas(nlev_cas,nt_cas),thlrad_cas(nlev_cas,nt_cas))
158	allocate(temp_nudg_cas(nlev_cas,nt_cas),th_nudg_cas(nlev_cas,nt_cas),thl_nudg_cas(nlev_cas,nt_cas))
159	allocate(qv_nudg_cas(nlev_cas,nt_cas),qt_nudg_cas(nlev_cas,nt_cas))
160	allocate(rv_nudg_cas(nlev_cas,nt_cas),rt_nudg_cas(nlev_cas,nt_cas))
161	allocate(u_nudg_cas(nlev_cas,nt_cas),v_nudg_cas(nlev_cas,nt_cas))
162	! flux
163	allocate(lat_cas(nt_cas),sens_cas(nt_cas),ustar_cas(nt_cas))
164	allocate(ts_cas(nt_cas),ps_cas(nt_cas),ps_forc_cas(nt_cas))
165	allocate(wpthetap_cas(nt_cas),wpqvp_cas(nt_cas),wpqtp_cas(nt_cas),wprvp_cas(nt_cas),wprtp_cas(nt_cas))
166
167	!champs interpoles
168	allocate(plev_prof_cas(nlev_cas))
169	allocate(t_prof_cas(nlev_cas))
170	allocate(th_prof_cas(nlev_cas))
171	allocate(thl_prof_cas(nlev_cas))
172	allocate(qt_prof_cas(nlev_cas))
173	allocate(qv_prof_cas(nlev_cas))
174	allocate(ql_prof_cas(nlev_cas))
175	allocate(qi_prof_cas(nlev_cas))
176	allocate(rh_prof_cas(nlev_cas))
177	allocate(rt_prof_cas(nlev_cas))
178	allocate(rv_prof_cas(nlev_cas))
179	allocate(rl_prof_cas(nlev_cas))
180	allocate(ri_prof_cas(nlev_cas))
181	allocate(u_prof_cas(nlev_cas))
182	allocate(v_prof_cas(nlev_cas))
183	allocate(w_prof_cas(nlev_cas))
184	allocate(omega_prof_cas(nlev_cas))
185	allocate(ug_prof_cas(nlev_cas))
186	allocate(vg_prof_cas(nlev_cas))
187	allocate(temp_nudg_prof_cas(nlev_cas),th_nudg_prof_cas(nlev_cas),thl_nudg_prof_cas(nlev_cas))
188	allocate(qt_nudg_prof_cas(nlev_cas),qv_nudg_prof_cas(nlev_cas),rt_nudg_prof_cas(nlev_cas),rv_nudg_prof_cas(nlev_cas))
189	allocate(u_nudg_prof_cas(nlev_cas),v_nudg_prof_cas(nlev_cas))
190
191	print*,'Allocations OK'
192
193	CALL read_SCM(nid,nlev_cas,nt_cas, &
194	& zzs_cas,pp_cas,zzforc_cas,ppforc_cas,temp0_cas,theta0_cas,thetal0_cas,qt0_cas,qv0_cas,ql0_cas,qi0_cas, &
195	& rh0_cas,rt0_cas,rv0_cas,rl0_cas,ri0_cas, &
196	& u0_cas,v0_cas,w_cas,omega_cas,ug_cas,vg_cas,uadv_cas,vadv_cas,tadv_cas,thadv_cas,thladv_cas, &
197	& qvadv_cas,qtadv_cas,rvadv_cas,rtadv_cas, &
198	& temp_nudg_cas,th_nudg_cas,thl_nudg_cas,qv_nudg_cas,qt_nudg_cas,rv_nudg_cas,rt_nudg_cas,u_nudg_cas,v_nudg_cas, &
199	& trad_cas,thrad_cas,thlrad_cas,tke_cas,sens_cas,lat_cas,ts_cas,ps_cas,ps_forc_cas,ustar_cas, &
200	& wpthetap_cas,wpqvp_cas,wpqtp_cas,wprvp_cas,wprtp_cas)
201
202	print*,'read_SCM cas OK'
203	do ii=1,nlev_cas
204	print*,'apres read2_SCM, plev_cas=',ii,ppforc_cas(ii,1)
205	!print*,'apres read_SCM, plev_cas=',ii,omega_cas(ii,nt_cas/2+1)
206	enddo
207
208
209	END SUBROUTINE read_SCM_cas
210
211
212	!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
213	SUBROUTINE deallocate2_1D_cases
214
215	deallocate(zzforc_cas)
216	deallocate(ppforc_cas)
217	!profils initiaux
218	deallocate(temp0_cas,theta0_cas,thetal0_cas)
219	deallocate(qt0_cas,qv0_cas,ql0_cas,qi0_cas,u0_cas,v0_cas)
220	deallocate(rt0_cas,rv0_cas,rl0_cas,ri0_cas,rh0_cas,tke_cas)
221	deallocate(t_cas,theta_cas,thl_cas,u_cas,v_cas)
222	deallocate(w_cas,omega_cas)
223	deallocate(ug_cas,vg_cas)
224	!advections et nudging
225	deallocate(uadv_cas,vadv_cas)
226	deallocate(tadv_cas,thadv_cas,thladv_cas)
227	deallocate(qtadv_cas,qvadv_cas)
228	deallocate(rtadv_cas,rvadv_cas)
229	deallocate(trad_cas,thrad_cas,thlrad_cas)
230	deallocate(temp_nudg_cas,th_nudg_cas,thl_nudg_cas)
231	deallocate(qv_nudg_cas,qt_nudg_cas)
232	deallocate(rv_nudg_cas,rt_nudg_cas)
233	deallocate(u_nudg_cas,v_nudg_cas)
234	! flux
235	deallocate(lat_cas,sens_cas,ustar_cas)
236	deallocate(ts_cas,ps_cas,ps_forc_cas)
237	deallocate(wpthetap_cas,wpqvp_cas,wpqtp_cas,wprvp_cas,wprtp_cas)
238
239	!champs interpoles
240	deallocate (plev_prof_cas)
241	deallocate (t_prof_cas)
242	deallocate (th_prof_cas)
243	deallocate (thl_prof_cas)
244	deallocate (qt_prof_cas)
245	deallocate (qv_prof_cas)
246	deallocate (ql_prof_cas)
247	deallocate (qi_prof_cas)
248	deallocate (rh_prof_cas)
249	deallocate (rt_prof_cas)
250	deallocate (rv_prof_cas)
251	deallocate (rl_prof_cas)
252	deallocate (ri_prof_cas)
253	deallocate (u_prof_cas)
254	deallocate (v_prof_cas)
255	deallocate (w_prof_cas)
256	deallocate (omega_prof_cas)
257	deallocate (ug_prof_cas)
258	deallocate (vg_prof_cas)
259	deallocate (temp_nudg_prof_cas,th_nudg_prof_cas,thl_nudg_prof_cas)
260	deallocate (qt_nudg_prof_cas,qv_nudg_prof_cas,rt_nudg_prof_cas,rv_nudg_prof_cas)
261	deallocate (u_nudg_prof_cas,v_nudg_prof_cas)
262
263	END SUBROUTINE deallocate2_1D_cases
264
265
266	!=====================================================================
267	SUBROUTINE read_SCM(nid,nlevel,ntime, &
268	& zzs,pp,zzforc,ppforc,temp0,theta0,thetal0,qt0,qv0,ql0,qi0,rh0,rt0,rv0,rl0,ri0, &
269	& u0,v0,w,omega,ug,vg,uadv,vadv,tadv,thadv,thladv,qvadv,qtadv,rvadv,rtadv, &
270	& temp_nudg,th_nudg,thl_nudg,qv_nudg,qt_nudg,rv_nudg,rt_nudg,u_nudg,v_nudg, &
271	& trad,thrad,thlrad,tke,sens,flat,ts,ps,ps_forc,ustar, &
272	& wpthetap,wpqvp,wpqtp,wprvp,wprtp)
273
274	!program reading forcing of the case study
275	implicit none
276	#include "netcdf.inc"
277	#include "compar1d_std.h"
278
279	integer ntime,nlevel,k,t
280
281	real zzs,zzforc(nlevel,ntime)
282	real pp,ppforc(nlevel,ntime)
283	!profils initiaux
284	real temp0(nlevel),theta0(nlevel),thetal0(nlevel),tke(nlevel)
285	real temp(nlevel,ntime),qv(nlevel,ntime),ql(nlevel,ntime),qi(nlevel,ntime),u(nlevel,ntime),v(nlevel,ntime)
286	real qt0(nlevel),qv0(nlevel),ql0(nlevel),qi0(nlevel),u0(nlevel),v0(nlevel)
287	real rt0(nlevel),rv0(nlevel),rl0(nlevel),ri0(nlevel),rh0(nlevel)
288	real w(nlevel,ntime),omega(nlevel,ntime)
289	real ug(nlevel,ntime),vg(nlevel,ntime)
290	!advections et nudging
291	real uadv(nlevel,ntime),vadv(nlevel,ntime)
292	real tadv(nlevel,ntime),thadv(nlevel,ntime),thladv(nlevel,ntime)
293	real qtadv(nlevel,ntime),qvadv(nlevel,ntime)
294	real rtadv(nlevel,ntime),rvadv(nlevel,ntime)
295	real trad(nlevel,ntime),thrad(nlevel,ntime),thlrad(nlevel,ntime)
296	real temp_nudg(nlevel,ntime),th_nudg(nlevel,ntime),thl_nudg(nlevel,ntime)
297	real qv_nudg(nlevel,ntime),qt_nudg(nlevel,ntime)
298	real rv_nudg(nlevel,ntime),rt_nudg(nlevel,ntime)
299	real u_nudg(nlevel,ntime),v_nudg(nlevel,ntime)
300	! flux
301	real flat(ntime),sens(ntime),ustar(ntime)
302	real ts(ntime),ps(ntime),ps_forc(ntime)
303	real wpthetap(ntime),wpqvp(ntime),wpqtp(ntime),wprtp(ntime),wprvp(ntime)
304	real resul(nlevel,ntime),resul1(nlevel),resul2(ntime),resul3
305
306
307	integer nid, ierr,ierr1,ierr2,rid,i
308	integer nbvar3d
309	parameter(nbvar3d=55)
310	integer var3didin(nbvar3d),missing_var(nbvar3d)
311	character*14 name_var(1:nbvar3d)
312
313
314	data name_var/ &
315	! coordonnees pression (n niveaux) profils intiaux #1-#15
316	& 'qt','qv','ql','qi','rt','rv','rl','ri', &
317	& 'rh','temp','theta','thetal','u','v','tke', &
318	! coordonnees pression (n niveaux) + temps #16-#42
319	& 'height_forc','pressure_forc','w','omega','ug','vg','u_adv','v_adv', &
320	& 'temp_adv','theta_adv','thetal_adv','qt_adv','qv_adv','rt_adv','rv_adv', &
321	& 'temp_rad','theta_rad','thetal_rad','temp_nudging','theta_nudging','thetal_nudging', &
322	& 'qv_nudging','qt_nudging','rv_nudging','rt_nudging','u_nudging','v_nudging', &
323	! coordonnees temps #43-#53
324	& 'sfc_sens_flx','sfc_lat_flx','ts','ps','ps_forc','ustar', &
325	& 'wpthetap','wpqvp','wpqtp','wprtp','wprvp', &
326	! scalaires #54-55
327	& 'height','pressure'/
328
329	!-----------------------------------------------------------------------
330	! Checking availability of variable #i in the cas.nc file
331	! missing_var=1 if the variable is missing
332	!-----------------------------------------------------------------------
333
334	do i=1,nbvar3d
335	missing_var(i)=0.
336	ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
337	if(ierr/=NF_NOERR) then
338	print *,'Variable manquante dans cas.nc:',i,name_var(i)
339	ierr=NF_NOERR
340	missing_var(i)=1
341	else
342
343	!-----------------------------------------------------------------------
344	! Activating keys depending on the presence of specific variables in cas.nc
345	!-----------------------------------------------------------------------
346	if ( 1 == 1 ) THEN
347	if ( name_var(i) == 'temp_nudging' .and. nint(nudging_temp)==0) stop 'Nudging inconsistency temp'
348	if ( name_var(i) == 'theta_nudging' .and. nint(nudging_theta)==0) stop 'Nudging inconsistency theta'
349	if ( name_var(i) == 'thetal_nudging' .and. nint(nudging_thetal)==0) stop 'Nudging inconsistency thetal'
350	if ( name_var(i) == 'qv_nudging' .and. nint(nudging_qv)==0) stop 'Nudging inconsistency qv'
351	if ( name_var(i) == 'qt_nudging' .and. nint(nudging_qt)==0) stop 'Nudging inconsistency qt'
352	if ( name_var(i) == 'rv_nudging' .and. nint(nudging_rv)==0) stop 'Nudging inconsistency rv'
353	if ( name_var(i) == 'rt_nudging' .and. nint(nudging_rt)==0) stop 'Nudging inconsistency rt'
354	if ( name_var(i) == 'u_nudging' .and. nint(nudging_u)==0) stop 'Nudging inconsistency u'
355	if ( name_var(i) == 'v_nudging' .and. nint(nudging_u)==0) stop 'Nudging inconsistency v'
356	ELSE
357	print*,'GUIDAGE : CONSISTENCY CHECK DEACTIVATED FOR TESTS of SANDU/REF'
358	ENDIF
359
360	!-----------------------------------------------------------------------
361	! Reading 1D (N) vertical varialbes (nlevel,lat,lon)
362	!-----------------------------------------------------------------------
363	if(i.LE.15) then
364	#ifdef NC_DOUBLE
365	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul1)
366	#else
367	ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul1)
368	#endif
369	print *,'read2_cas(resul1), on a lu ',i,name_var(i)
370	if(ierr/=NF_NOERR) then
371	print *,'Pb a la lecture de cas.nc: ',name_var(i)
372	stop "getvarup"
373	endif
374	print*,'Lecture de la variable #i ',i,name_var(i),minval(resul1),maxval(resul1)
375
376	!-----------------------------------------------------------------------
377	! Reading 2D tim-vertical variables (time,nlevel,lat,lon)
378	! TBD : seems to be the same as above.
379	!-----------------------------------------------------------------------
380	else if(i.ge.16.and.i.LE.42) then
381	#ifdef NC_DOUBLE
382	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
383	#else
384	ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
385	#endif
386	print *,'read2_cas(resul), on a lu ',i,name_var(i)
387	if(ierr/=NF_NOERR) then
388	print *,'Pb a la lecture de cas.nc: ',name_var(i)
389	stop "getvarup"
390	endif
391	print*,'Lecture de la variable #i ',i,name_var(i),minval(resul),maxval(resul)
392
393	!-----------------------------------------------------------------------
394	! Reading 1D time variables (time,lat,lon)
395	!-----------------------------------------------------------------------
396	else if (i.gt.43.and.i.LE.53) then
397	#ifdef NC_DOUBLE
398	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul2)
399	#else
400	ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul2)
401	#endif
402	print *,'read2_cas(resul2), on a lu ',i,name_var(i)
403	if(ierr/=NF_NOERR) then
404	print *,'Pb a la lecture de cas.nc: ',name_var(i)
405	stop "getvarup"
406	endif
407	print*,'Lecture de la variable #i ',i,name_var(i),minval(resul2),maxval(resul2)
408	!-----------------------------------------------------------------------
409	! Reading scalar variables (t0,lat,lon)
410	!-----------------------------------------------------------------------
411	else
412	#ifdef NC_DOUBLE
413	ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul3)
414	#else
415	ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul3)
416	#endif
417	print *,'read2_cas(resul3), on a lu ',i,name_var(i)
418	if(ierr/=NF_NOERR) then
419	print *,'Pb a la lecture de cas.nc: ',name_var(i)
420	stop "getvarup"
421	endif
422	print*,'Lecture de la variable #i ',i,name_var(i),resul3
423	endif
424	endif
425
426	!-----------------------------------------------------------------------
427	! Attributing variables
428	!-----------------------------------------------------------------------
429	select case(i)
430	case(1) ; qt0 =resul1
431	case(2) ; qv0 =resul1
432	case(3) ; ql0 =resul1
433	case(4) ; qi0 =resul1
434	case(5) ; rt0 =resul1
435	case(6) ; rv0 =resul1
436	case(7) ; rl0 =resul1
437	case(8) ; ri0 =resul1
438	case(9) ; rh0 =resul1
439	case(10) ; temp0 =resul1
440	case(11) ; theta0 =resul1
441	case(12) ; thetal0 =resul1
442	case(13) ; u0 =resul1
443	case(14) ; v0 =resul1
444	case(15) ; tke =resul1
445	case(16) ; zzforc =resul ! donnees indexees en nlevel,time
446	case(17) ; ppforc =resul
447	case(18) ; w =resul
448	case(19) ; omega =resul
449	case(20) ; ug =resul
450	case(21) ; vg =resul
451	case(22) ; uadv =resul
452	case(23) ; vadv =resul
453	case(24) ; tadv =resul
454	case(25) ; thadv =resul
455	case(26) ; thladv =resul
456	case(27) ; qtadv =resul
457	case(28) ; qvadv =resul
458	case(29) ; rtadv =resul
459	case(30) ; rvadv =resul
460	case(31) ; trad =resul
461	case(32) ; thrad =resul
462	case(33) ; thlrad =resul
463	case(34) ; temp_nudg =resul
464	case(35) ; th_nudg =resul
465	case(36) ; thl_nudg =resul
466	case(37) ; qv_nudg =resul
467	case(38) ; qt_nudg =resul
468	case(39) ; rv_nudg =resul
469	case(40) ; rt_nudg =resul
470	case(41) ; u_nudg =resul
471	case(42) ; v_nudg =resul
472	case(43) ; sens =resul2 ! donnees indexees en time seulement
473	case(44) ; flat =resul2
474	case(45) ; ts =resul2
475	case(46) ; ps =resul2
476	case(47) ; ps_forc =resul2
477	case(48) ; ustar =resul2
478	case(49) ; wpthetap =resul2
479	case(50) ; wpqvp =resul2
480	case(51) ; wpqtp =resul2
481	case(52) ; wprvp =resul2
482	case(53) ; wprtp =resul2
483	case(54) ; zzs =resul3 ! scalaires
484	case(55) ; pp =resul3
485	end select
486	resul=0.
487	resul1=0.
488	resul2=0.
489	resul3=0.
490	enddo
491	! print*,'Lecture de la variable APRES ,sens ',minval(sens),maxval(sens)
492	! print*,'Lecture de la variable APRES ,flat ',minval(flat),maxval(flat)
493
494	!CR:ATTENTION EN ATTENTE DE REGLER LA QUESTION DU PAS DE TEMPS INITIAL
495	! do t=1,ntime
496	! do k=1,nlevel
497	! temp(k,t)=temp0(k)
498	! qv(k,t)=qv0(k)
499	! ql(k,t)=ql0(k)
500	! qi(k,t)=qi0(k)
501	! u(k,t)=u0(k)
502	! v(k,t)=v0(k)
503	! !tke(k,t)=tke0(k)
504	! enddo
505	! enddo
506	!!!! TRAVAIL : EN FONCTION DES DECISIONS SUR LA SPECIFICATION DE W
507	!!!omega=-vitwpresrg/(rd*temp)
508	!-----------------------------------------------------------------------
509
510	return
511	END SUBROUTINE read_SCM
512	!======================================================================
513
514	!======================================================================
515
516	!**********************************************************************************************
517	SUBROUTINE interp_case_time_std(day,day1,annee_ref &
518	& ,nt_cas,nlev_cas &
519	& ,ts_cas,ps_cas,ps_forc_cas,plev_cas,ppforc_cas,t_cas,th_cas,thl_cas &
520	& ,qt_cas,qv_cas,ql_cas,qi_cas &
521	& ,rt_cas,rv_cas,rl_cas,ri_cas,rh_cas &
522	& ,u_cas,v_cas,w_cas,omega_cas,ug_cas,vg_cas &
523	& ,temp_nudg_cas,th_nudg_cas,thl_nudg_cas,qt_nudg_cas,qv_nudg_cas &
524	& ,rt_nudg_cas,rv_nudg_cas,u_nudg_cas,v_nudg_cas &
525	& ,uadv_cas,vadv_cas,tadv_cas,thadv_cas,thladv_cas &
526	& ,qtadv_cas,qvadv_cas,rtadv_cas,rvadv_cas &
527	& ,trad_cas,thrad_cas,thlrad_cas &
528	& ,tke_cas,lat_cas,sens_cas,ustar_cas &
529	& ,wpthetap_cas,wpqtp_cas,wpqvp_cas,wprtp_cas,wprvp_cas &
530	!
531	& ,ts_prof_cas,ps_prof_cas,ps_forc_prof_cas,plev_prof_cas,pforc_prof_cas&
532	& ,t_prof_cas,th_prof_cas,thl_prof_cas &
533	& ,qt_prof_cas,qv_prof_cas,ql_prof_cas,qi_prof_cas &
534	& ,rt_prof_cas,rv_prof_cas,rl_prof_cas,ri_prof_cas,rh_prof_cas &
535	& ,u_prof_cas,v_prof_cas,w_prof_cas,omega_prof_cas &
536	& ,ug_prof_cas,vg_prof_cas &
537	& ,temp_nudg_prof_cas,th_nudg_prof_cas,thl_nudg_prof_cas &
538	& ,qt_nudg_prof_cas,qv_nudg_prof_cas &
539	& ,rt_nudg_prof_cas,rv_nudg_prof_cas,u_nudg_prof_cas,v_nudg_prof_cas &
540	& ,uadv_prof_cas,vadv_prof_cas,tadv_prof_cas,thadv_prof_cas,thladv_prof_cas&
541	& ,qtadv_prof_cas,qvadv_prof_cas,rtadv_prof_cas,rvadv_prof_cas &
542	& ,trad_prof_cas,thrad_prof_cas,thlrad_prof_cas &
543	& ,tke_prof_cas,lat_prof_cas,sens_prof_cas,ustar_prof_cas &
544	& ,wpthetap_prof_cas,wpqtp_prof_cas,wpqvp_prof_cas,wprtp_prof_cas,wprvp_prof_cas)
545
546
547	implicit none
548
549	!---------------------------------------------------------------------------------------
550	! Time interpolation of a 2D field to the timestep corresponding to day
551	!
552	! day: current julian day (e.g. 717538.2)
553	! day1: first day of the simulation
554	! nt_cas: total nb of data in the forcing
555	! pdt_cas: total time interval (in sec) between 2 forcing data
556	!---------------------------------------------------------------------------------------
557
558	#include "compar1d_std.h"
559	#include "date_cas.h"
560
561	! inputs:
562	integer annee_ref
563	integer nt_cas,nlev_cas
564	real day, day1,day_cas
565	real ts_cas(nt_cas),ps_cas(nt_cas),ps_forc_cas(nt_cas)
566	real plev_cas(nlev_cas,nt_cas),ppforc_cas(nt_cas)
567	real t_cas(nlev_cas,nt_cas),th_cas(nlev_cas,nt_cas),thl_cas(nlev_cas,nt_cas)
568	real qt_cas(nlev_cas,nt_cas),qv_cas(nlev_cas,nt_cas),ql_cas(nlev_cas,nt_cas),qi_cas(nlev_cas,nt_cas)
569	real rt_cas(nlev_cas,nt_cas),rv_cas(nlev_cas,nt_cas),rl_cas(nlev_cas,nt_cas),ri_cas(nlev_cas,nt_cas)
570	real rh_cas(nlev_cas,nt_cas),u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas)
571	real w_cas(nlev_cas,nt_cas),omega_cas(nlev_cas,nt_cas)
572	real ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas)
573	real temp_nudg_cas(nlev_cas,nt_cas),th_nudg_cas(nlev_cas,nt_cas),thl_nudg_cas(nlev_cas,nt_cas)
574	real qt_nudg_cas(nlev_cas,nt_cas),qv_nudg_cas(nlev_cas,nt_cas)
575	real rt_nudg_cas(nlev_cas,nt_cas),rv_nudg_cas(nlev_cas,nt_cas)
576	real u_nudg_cas(nlev_cas,nt_cas),v_nudg_cas(nlev_cas,nt_cas)
577	real uadv_cas(nlev_cas,nt_cas),vadv_cas(nlev_cas,nt_cas)
578	real tadv_cas(nlev_cas,nt_cas),thadv_cas(nlev_cas,nt_cas),thladv_cas(nlev_cas,nt_cas)
579	real qtadv_cas(nlev_cas,nt_cas),qvadv_cas(nlev_cas,nt_cas)
580	real rtadv_cas(nlev_cas,nt_cas),rvadv_cas(nlev_cas,nt_cas)
581	real trad_cas(nlev_cas,nt_cas),thrad_cas(nlev_cas,nt_cas),thlrad_cas(nlev_cas,nt_cas)
582	real lat_cas(nt_cas),sens_cas(nt_cas),tke_cas(nt_cas)
583	real wpthetap_cas(nt_cas),wpqtp_cas(nt_cas),wpqvp_cas(nt_cas)
584	real ustar_cas(nt_cas),wprtp_cas(nt_cas),wprvp_cas(nt_cas)
585
586	! output:
587	real plev_prof_cas(nlev_cas),pforc_prof_cas(nt_cas)
588	real t_prof_cas(nlev_cas),th_prof_cas(nlev_cas),thl_prof_cas(nlev_cas)
589	real qt_prof_cas(nlev_cas),qv_prof_cas(nlev_cas),ql_prof_cas(nlev_cas),qi_prof_cas(nlev_cas)
590	real rt_prof_cas(nlev_cas),rv_prof_cas(nlev_cas),rl_prof_cas(nlev_cas),ri_prof_cas(nlev_cas)
591	real rh_prof_cas(nlev_cas),u_prof_cas(nlev_cas),v_prof_cas(nlev_cas)
592	real w_prof_cas(nlev_cas),omega_prof_cas(nlev_cas)
593	real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas)
594	real temp_nudg_prof_cas(nlev_cas),th_nudg_prof_cas(nlev_cas),thl_nudg_prof_cas(nlev_cas)
595	real qt_nudg_prof_cas(nlev_cas),qv_nudg_prof_cas(nlev_cas)
596	real rt_nudg_prof_cas(nlev_cas),rv_nudg_prof_cas(nlev_cas)
597	real u_nudg_prof_cas(nlev_cas),v_nudg_prof_cas(nlev_cas)
598	real uadv_prof_cas(nlev_cas),vadv_prof_cas(nlev_cas)
599	real tadv_prof_cas(nlev_cas),thadv_prof_cas(nlev_cas),thladv_prof_cas(nlev_cas)
600	real qtadv_prof_cas(nlev_cas),qvadv_prof_cas(nlev_cas)
601	real rtadv_prof_cas(nlev_cas),rvadv_prof_cas(nlev_cas)
602	real trad_prof_cas(nlev_cas),thrad_prof_cas(nlev_cas),thlrad_prof_cas(nlev_cas)
603	real lat_prof_cas,sens_prof_cas,tke_prof_cas
604	real ts_prof_cas,ps_prof_cas,ps_forc_prof_cas
605	real wpthetap_prof_cas,wpqtp_prof_cas,wpqvp_prof_cas
606	real ustar_prof_cas,wprtp_prof_cas,wprvp_prof_cas
607
608	! local:
609	integer it_cas1, it_cas2,k
610	real timeit,time_cas1,time_cas2,frac
611
612	print*,'Check time',day1,day_ju_ini_cas,day_deb+1,pdt_cas
613	! do k=1,nlev_cas
614	! print*,'debut de interp_case_time, plev_cas=',k,plev_cas(k,1)
615	! enddo
616
617	! On teste si la date du cas AMMA est correcte.
618	! C est pour memoire car en fait les fichiers .def
619	! sont censes etre corrects.
620	! A supprimer a terme (MPL 20150623)
621	! if ((forcing_type.eq.10).and.(1.eq.0)) then
622	! Check that initial day of the simulation consistent with AMMA case:
623	! if (annee_ref.ne.2006) then
624	! print*,'Pour AMMA, annee_ref doit etre 2006'
625	! print*,'Changer annee_ref dans run.def'
626	! stop
627	! endif
628	! if (annee_ref.eq.2006 .and. day1.lt.day_cas) then
629	! print*,'AMMA a debute le 10 juillet 2006',day1,day_cas
630	! print*,'Changer dayref dans run.def'
631	! stop
632	! endif
633	! if (annee_ref.eq.2006 .and. day1.gt.day_cas+1) then
634	! print*,'AMMA a fini le 11 juillet'
635	! print*,'Changer dayref ou nday dans run.def'
636	! stop
637	! endif
638	! endif
639
640	! Determine timestep relative to the 1st day:
641	! timeit=(day-day1)*86400.
642	! if (annee_ref.eq.1992) then
643	! timeit=(day-day_cas)*86400.
644	! else
645	! timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
646	! endif
647	timeit=(day-day_ju_ini_cas)*86400
648	print *,'day=',day
649	print *,'day_ju_ini_cas=',day_ju_ini_cas
650	print *,'pdt_cas=',pdt_cas
651	print *,'timeit=',timeit
652	print *,'nt_cas=',nt_cas
653
654	! Determine the closest observation times:
655	! it_cas1=INT(timeit/pdt_cas)+1
656	! it_cas2=it_cas1 + 1
657	! time_cas1=(it_cas1-1)*pdt_cas
658	! time_cas2=(it_cas2-1)*pdt_cas
659
660	it_cas1=INT(timeit/pdt_cas)+1
661	IF (it_cas1 .EQ. nt_cas) THEN
662	it_cas2=it_cas1
663	ELSE
664	it_cas2=it_cas1 + 1
665	ENDIF
666	time_cas1=(it_cas1-1)*pdt_cas
667	time_cas2=(it_cas2-1)*pdt_cas
668	! print *,'timeit,pdt_cas,nt_cas=',timeit,pdt_cas,nt_cas
669	! print *,'it_cas1,it_cas2,time_cas1,time_cas2=',it_cas1,it_cas2,time_cas1,time_cas2
670
671	if (it_cas1 .gt. nt_cas) then
672	write(,) 'PB-stop: day, day_ju_ini_cas,it_cas1, it_cas2, timeit: ' &
673	& ,day,day_ju_ini_cas,it_cas1,it_cas2,timeit
674	stop
675	endif
676
677	! time interpolation:
678	IF (it_cas1 .EQ. it_cas2) THEN
679	frac=0.
680	ELSE
681	frac=(time_cas2-timeit)/(time_cas2-time_cas1)
682	frac=max(frac,0.0)
683	ENDIF
684
685	lat_prof_cas = lat_cas(it_cas2) &
686	& -frac*(lat_cas(it_cas2)-lat_cas(it_cas1))
687	sens_prof_cas = sens_cas(it_cas2) &
688	& -frac*(sens_cas(it_cas2)-sens_cas(it_cas1))
689	tke_prof_cas = tke_cas(it_cas2) &
690	& -frac*(tke_cas(it_cas2)-tke_cas(it_cas1))
691	ts_prof_cas = ts_cas(it_cas2) &
692	& -frac*(ts_cas(it_cas2)-ts_cas(it_cas1))
693	ps_prof_cas = ps_cas(it_cas2) &
694	& -frac*(ps_cas(it_cas2)-ps_cas(it_cas1))
695	ps_forc_prof_cas = ps_forc_cas(it_cas2) &
696	& -frac*(ps_forc_cas(it_cas2)-ps_forc_cas(it_cas1))
697	ustar_prof_cas = ustar_cas(it_cas2) &
698	& -frac*(ustar_cas(it_cas2)-ustar_cas(it_cas1))
699	wpthetap_prof_cas = wpthetap_cas(it_cas2) &
700	& -frac*(wpthetap_cas(it_cas2)-wpthetap_cas(it_cas1))
701	wpqtp_prof_cas = wpqtp_cas(it_cas2) &
702	& -frac*(wpqtp_cas(it_cas2)-wpqtp_cas(it_cas1))
703	wpqvp_prof_cas = wpqvp_cas(it_cas2) &
704	& -frac*(wpqvp_cas(it_cas2)-wpqvp_cas(it_cas1))
705	wprtp_prof_cas = wprtp_cas(it_cas2) &
706	& -frac*(wprtp_cas(it_cas2)-wprtp_cas(it_cas1))
707	wprvp_prof_cas = wprvp_cas(it_cas2) &
708	& -frac*(wprvp_cas(it_cas2)-wprvp_cas(it_cas1))
709
710	do k=1,nlev_cas
711	plev_prof_cas(k) = plev_cas(k,it_cas2) &
712	& -frac*(plev_cas(k,it_cas2)-plev_cas(k,it_cas1))
713	t_prof_cas(k) = t_cas(k,it_cas2) &
714	& -frac*(t_cas(k,it_cas2)-t_cas(k,it_cas1))
715	!print *,'k,frac,plev_cas1,plev_cas2=',k,frac,plev_cas(k,it_cas1),plev_cas(k,it_cas2)
716	th_prof_cas(k) = th_cas(k,it_cas2) &
717	& -frac*(th_cas(k,it_cas2)-th_cas(k,it_cas1))
718	thl_prof_cas(k) = thl_cas(k,it_cas2) &
719	& -frac*(thl_cas(k,it_cas2)-thl_cas(k,it_cas1))
720	qt_prof_cas(k) = qt_cas(k,it_cas2) &
721	& -frac*(qt_cas(k,it_cas2)-qt_cas(k,it_cas1))
722	qv_prof_cas(k) = qv_cas(k,it_cas2) &
723	& -frac*(qv_cas(k,it_cas2)-qv_cas(k,it_cas1))
724	ql_prof_cas(k) = ql_cas(k,it_cas2) &
725	& -frac*(ql_cas(k,it_cas2)-ql_cas(k,it_cas1))
726	qi_prof_cas(k) = qi_cas(k,it_cas2) &
727	& -frac*(qi_cas(k,it_cas2)-qi_cas(k,it_cas1))
728	rt_prof_cas(k) = rt_cas(k,it_cas2) &
729	& -frac*(rt_cas(k,it_cas2)-rt_cas(k,it_cas1))
730	rv_prof_cas(k) = rv_cas(k,it_cas2) &
731	& -frac*(rv_cas(k,it_cas2)-rv_cas(k,it_cas1))
732	rl_prof_cas(k) = rl_cas(k,it_cas2) &
733	& -frac*(rl_cas(k,it_cas2)-rl_cas(k,it_cas1))
734	ri_prof_cas(k) = ri_cas(k,it_cas2) &
735	& -frac*(ri_cas(k,it_cas2)-ri_cas(k,it_cas1))
736	rh_prof_cas(k) = rh_cas(k,it_cas2) &
737	& -frac*(rh_cas(k,it_cas2)-rh_cas(k,it_cas1))
738	u_prof_cas(k) = u_cas(k,it_cas2) &
739	& -frac*(u_cas(k,it_cas2)-u_cas(k,it_cas1))
740	v_prof_cas(k) = v_cas(k,it_cas2) &
741	& -frac*(v_cas(k,it_cas2)-v_cas(k,it_cas1))
742	w_prof_cas(k) = w_cas(k,it_cas2) &
743	& -frac*(w_cas(k,it_cas2)-w_cas(k,it_cas1))
744	omega_prof_cas(k) = omega_cas(k,it_cas2) &
745	& -frac*(omega_cas(k,it_cas2)-omega_cas(k,it_cas1))
746	ug_prof_cas(k) = ug_cas(k,it_cas2) &
747	& -frac*(ug_cas(k,it_cas2)-ug_cas(k,it_cas1))
748	vg_prof_cas(k) = vg_cas(k,it_cas2) &
749	& -frac*(vg_cas(k,it_cas2)-vg_cas(k,it_cas1))
750	temp_nudg_prof_cas(k) = temp_nudg_cas(k,it_cas2) &
751	& -frac*(temp_nudg_cas(k,it_cas2)-temp_nudg_cas(k,it_cas1))
752	th_nudg_prof_cas(k) = th_nudg_cas(k,it_cas2) &
753	& -frac*(th_nudg_cas(k,it_cas2)-th_nudg_cas(k,it_cas1))
754	thl_nudg_prof_cas(k) = thl_nudg_cas(k,it_cas2) &
755	& -frac*(thl_nudg_cas(k,it_cas2)-thl_nudg_cas(k,it_cas1))
756	qt_nudg_prof_cas(k) = qt_nudg_cas(k,it_cas2) &
757	& -frac*(qt_nudg_cas(k,it_cas2)-qt_nudg_cas(k,it_cas1))
758	qv_nudg_prof_cas(k) = qv_nudg_cas(k,it_cas2) &
759	& -frac*(qv_nudg_cas(k,it_cas2)-qv_nudg_cas(k,it_cas1))
760	rt_nudg_prof_cas(k) = rt_nudg_cas(k,it_cas2) &
761	& -frac*(rt_nudg_cas(k,it_cas2)-rt_nudg_cas(k,it_cas1))
762	rv_nudg_prof_cas(k) = rv_nudg_cas(k,it_cas2) &
763	& -frac*(rv_nudg_cas(k,it_cas2)-rv_nudg_cas(k,it_cas1))
764	u_nudg_prof_cas(k) = u_nudg_cas(k,it_cas2) &
765	& -frac*(u_nudg_cas(k,it_cas2)-u_nudg_cas(k,it_cas1))
766	v_nudg_prof_cas(k) = v_nudg_cas(k,it_cas2) &
767	& -frac*(v_nudg_cas(k,it_cas2)-v_nudg_cas(k,it_cas1))
768	uadv_prof_cas(k) = uadv_cas(k,it_cas2) &
769	& -frac*(uadv_cas(k,it_cas2)-uadv_cas(k,it_cas1))
770	vadv_prof_cas(k) = vadv_cas(k,it_cas2) &
771	& -frac*(vadv_cas(k,it_cas2)-vadv_cas(k,it_cas1))
772	tadv_prof_cas(k) = tadv_cas(k,it_cas2) &
773	& -frac*(tadv_cas(k,it_cas2)-tadv_cas(k,it_cas1))
774	thadv_prof_cas(k) = thadv_cas(k,it_cas2) &
775	& -frac*(thadv_cas(k,it_cas2)-thadv_cas(k,it_cas1))
776	thladv_prof_cas(k) = thladv_cas(k,it_cas2) &
777	& -frac*(thladv_cas(k,it_cas2)-thladv_cas(k,it_cas1))
778	qtadv_prof_cas(k) = qtadv_cas(k,it_cas2) &
779	& -frac*(qtadv_cas(k,it_cas2)-qtadv_cas(k,it_cas1))
780	qvadv_prof_cas(k) = qvadv_cas(k,it_cas2) &
781	& -frac*(qvadv_cas(k,it_cas2)-qvadv_cas(k,it_cas1))
782	rtadv_prof_cas(k) = rtadv_cas(k,it_cas2) &
783	& -frac*(rtadv_cas(k,it_cas2)-rtadv_cas(k,it_cas1))
784	rvadv_prof_cas(k) = rvadv_cas(k,it_cas2) &
785	& -frac*(rvadv_cas(k,it_cas2)-rvadv_cas(k,it_cas1))
786	trad_prof_cas(k) = trad_cas(k,it_cas2) &
787	& -frac*(trad_cas(k,it_cas2)-trad_cas(k,it_cas1))
788	thrad_prof_cas(k) = thrad_cas(k,it_cas2) &
789	& -frac*(thrad_cas(k,it_cas2)-thrad_cas(k,it_cas1))
790	thlrad_prof_cas(k) = thlrad_cas(k,it_cas2) &
791	& -frac*(thlrad_cas(k,it_cas2)-thlrad_cas(k,it_cas1))
792	enddo
793
794	return
795	END SUBROUTINE interp_case_time_std
796
797	!**********************************************************************************************
798	!=====================================================================
799	SUBROUTINE interp_case_vertical_std(nlev_cas &
800	& ,plev_prof_cas,t_prof_cas,th_prof_cas,thl_prof_cas &
801	& ,qt_prof_cas,qv_prof_cas,ql_prof_cas,qi_prof_cas &
802	& ,rt_prof_cas,rv_prof_cas,rl_prof_cas,ri_prof_cas,rh_prof_cas &
803	& ,u_prof_cas,v_prof_cas,w_prof_cas,omega_prof_cas &
804	& ,ug_prof_cas,vg_prof_cas &
805	& ,temp_nudg_prof_cas,th_nudg_prof_cas,thl_nudg_prof_cas &
806	& ,qt_nudg_prof_cas,qv_nudg_prof_cas &
807	& ,rt_nudg_prof_cas,rv_nudg_prof_cas,u_nudg_prof_cas,v_nudg_prof_cas &
808	& ,uadv_prof_cas,vadv_prof_cas,tadv_prof_cas,thadv_prof_cas,thladv_prof_cas &
809	& ,qtadv_prof_cas,qvadv_prof_cas,rtadv_prof_cas,rvadv_prof_cas &
810	& ,trad_prof_cas,thrad_prof_cas,thlrad_prof_cas &
811	!
812	& ,plev_mod_cas,t_mod_cas,th_mod_cas,thl_mod_cas &
813	& ,qt_mod_cas,qv_mod_cas,ql_mod_cas,qi_mod_cas &
814	& ,rt_mod_cas,rv_mod_cas,rl_mod_cas,ri_mod_cas,rh_mod_cas &
815	& ,u_mod_cas,v_mod_cas,w_mod_cas,omega_mod_cas &
816	& ,ug_mod_cas,vg_mod_cas &
817	& ,temp_nudg_mod_cas,th_nudg_mod_cas,thl_nudg_mod_cas &
818	& ,qt_nudg_mod_cas,qv_nudg_mod_cas &
819	& ,rt_nudg_mod_cas,rv_nudg_mod_cas,u_nudg_mod_cas,v_nudg_mod_cas &
820	& ,uadv_mod_cas,vadv_mod_cas,tadv_mod_cas,thadv_mod_cas,thladv_mod_cas &
821	& ,qtadv_mod_cas,qvadv_mod_cas,rtadv_mod_cas,rvadv_mod_cas &
822	& ,trad_mod_cas,thrad_mod_cas,thlrad_mod_cas)
823
824	implicit none
825
826	#include "YOMCST.h"
827	#include "dimensions.h"
828
829	!-------------------------------------------------------------------------
830	! Vertical interpolation of generic case forcing data onto mod_casel levels
831	!-------------------------------------------------------------------------
832
833	integer nlevmax
834	parameter (nlevmax=41)
835	integer nlev_cas,mxcalc
836	! real play(llm), plev_prof(nlevmax)
837	! real t_prof(nlevmax),q_prof(nlevmax)
838	! real u_prof(nlevmax),v_prof(nlevmax), w_prof(nlevmax)
839	! real ht_prof(nlevmax),vt_prof(nlevmax)
840	! real hq_prof(nlevmax),vq_prof(nlevmax)
841	real plev_prof_cas(nlev_cas)
842	real t_prof_cas(nlev_cas),th_prof_cas(nlev_cas),thl_prof_cas(nlev_cas)
843	real qt_prof_cas(nlev_cas),qv_prof_cas(nlev_cas),ql_prof_cas(nlev_cas),qi_prof_cas(nlev_cas)
844	real rt_prof_cas(nlev_cas),rv_prof_cas(nlev_cas),rl_prof_cas(nlev_cas),ri_prof_cas(nlev_cas)
845	real rh_prof_cas(nlev_cas)
846	real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas),w_prof_cas(nlev_cas),omega_prof_cas(nlev_cas)
847	real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas)
848	real temp_nudg_prof_cas(nlev_cas),th_nudg_prof_cas(nlev_cas),thl_nudg_prof_cas(nlev_cas)
849	real qt_nudg_prof_cas(nlev_cas),qv_nudg_prof_cas(nlev_cas)
850	real rt_nudg_prof_cas(nlev_cas),rv_nudg_prof_cas(nlev_cas)
851	real u_nudg_prof_cas(nlev_cas),v_nudg_prof_cas(nlev_cas)
852	real uadv_prof_cas(nlev_cas),vadv_prof_cas(nlev_cas)
853	real tadv_prof_cas(nlev_cas),thadv_prof_cas(nlev_cas),thladv_prof_cas(nlev_cas)
854	real qtadv_prof_cas(nlev_cas),qvadv_prof_cas(nlev_cas)
855	real rtadv_prof_cas(nlev_cas),rvadv_prof_cas(nlev_cas)
856	real trad_prof_cas(nlev_cas),thrad_prof_cas(nlev_cas),thlrad_prof_cas(nlev_cas)
857
858	real play(llm),plev_mod_cas(llm),t_mod_cas(llm),th_mod_cas(llm),thl_mod_cas(llm)
859	real qt_mod_cas(llm),qv_mod_cas(llm),ql_mod_cas(llm),qi_mod_cas(llm)
860	real rt_mod_cas(llm),rv_mod_cas(llm),rl_mod_cas(llm),ri_mod_cas(llm)
861	real rh_mod_cas(llm)
862	real u_mod_cas(llm),v_mod_cas(llm),w_mod_cas(llm),omega_mod_cas(llm)
863	real ug_mod_cas(llm),vg_mod_cas(llm)
864	real temp_nudg_mod_cas(llm),th_nudg_mod_cas(llm),thl_nudg_mod_cas(llm)
865	real qt_nudg_mod_cas(llm),qv_nudg_mod_cas(llm)
866	real rt_nudg_mod_cas(llm),rv_nudg_mod_cas(llm),u_nudg_mod_cas(llm),v_nudg_mod_cas(llm)
867	real uadv_mod_cas(llm),vadv_mod_cas(llm)
868	real tadv_mod_cas(llm),thadv_mod_cas(llm),thladv_mod_cas(llm)
869	real qtadv_mod_cas(llm),qvadv_mod_cas(llm)
870	real rtadv_mod_cas(llm),rvadv_mod_cas(llm)
871	real trad_mod_cas(llm),thrad_mod_cas(llm),thlrad_mod_cas(llm)
872
873	integer l,k,k1,k2
874	real frac,frac1,frac2,fact
875
876	! do l = 1, llm
877	! print *,'debut interp, play=',l,play(l)
878	! enddo
879	! do l = 1, nlev_cas
880	! print *,'debut interp, plev_prof_cas=',l,play(l),plev_prof_cas(l)
881	! enddo
882
883	do l = 1, llm
884
885	if (play(l).ge.plev_prof_cas(nlev_cas)) then
886
887	mxcalc=l
888	! print *,'debut interp, mxcalc=',mxcalc
889	k1=0
890	k2=0
891
892	if (play(l).le.plev_prof_cas(1)) then
893
894	do k = 1, nlev_cas-1
895	if (play(l).le.plev_prof_cas(k).and. play(l).gt.plev_prof_cas(k+1)) then
896	k1=k
897	k2=k+1
898	endif
899	enddo
900
901	if (k1.eq.0 .or. k2.eq.0) then
902	write(,) 'PB! k1, k2 = ',k1,k2
903	write(,) 'l,play(l) = ',l,play(l)/100
904	do k = 1, nlev_cas-1
905	write(,) 'k,plev_prof_cas(k) = ',k,plev_prof_cas(k)/100
906	enddo
907	endif
908
909	frac = (plev_prof_cas(k2)-play(l))/(plev_prof_cas(k2)-plev_prof_cas(k1))
910	t_mod_cas(l)= t_prof_cas(k2) - frac*(t_prof_cas(k2)-t_prof_cas(k1))
911	th_mod_cas(l)= th_prof_cas(k2) - frac*(th_prof_cas(k2)-th_prof_cas(k1))
912	if(th_mod_cas(l).NE.0) t_mod_cas(l)= th_mod_cas(l)(play(l)/100000.)*(RD/RCPD)
913	thl_mod_cas(l)= thl_prof_cas(k2) - frac*(thl_prof_cas(k2)-thl_prof_cas(k1))
914	qt_mod_cas(l)= qt_prof_cas(k2) - frac*(qt_prof_cas(k2)-qt_prof_cas(k1))
915	qv_mod_cas(l)= qv_prof_cas(k2) - frac*(qv_prof_cas(k2)-qv_prof_cas(k1))
916	ql_mod_cas(l)= ql_prof_cas(k2) - frac*(ql_prof_cas(k2)-ql_prof_cas(k1))
917	qi_mod_cas(l)= qi_prof_cas(k2) - frac*(qi_prof_cas(k2)-qi_prof_cas(k1))
918	rt_mod_cas(l)= rt_prof_cas(k2) - frac*(rt_prof_cas(k2)-rt_prof_cas(k1))
919	rv_mod_cas(l)= rv_prof_cas(k2) - frac*(rv_prof_cas(k2)-rv_prof_cas(k1))
920	rl_mod_cas(l)= rl_prof_cas(k2) - frac*(rl_prof_cas(k2)-rl_prof_cas(k1))
921	ri_mod_cas(l)= ri_prof_cas(k2) - frac*(ri_prof_cas(k2)-ri_prof_cas(k1))
922	rh_mod_cas(l)= rh_prof_cas(k2) - frac*(rh_prof_cas(k2)-rh_prof_cas(k1))
923	u_mod_cas(l)= u_prof_cas(k2) - frac*(u_prof_cas(k2)-u_prof_cas(k1))
924	v_mod_cas(l)= v_prof_cas(k2) - frac*(v_prof_cas(k2)-v_prof_cas(k1))
925	w_mod_cas(l)= w_prof_cas(k2) - frac*(w_prof_cas(k2)-w_prof_cas(k1))
926	omega_mod_cas(l)= omega_prof_cas(k2) - frac*(omega_prof_cas(k2)-omega_prof_cas(k1))
927	ug_mod_cas(l)= ug_prof_cas(k2) - frac*(ug_prof_cas(k2)-ug_prof_cas(k1))
928	vg_mod_cas(l)= vg_prof_cas(k2) - frac*(vg_prof_cas(k2)-vg_prof_cas(k1))
929	temp_nudg_mod_cas(l)= temp_nudg_prof_cas(k2) - frac*(temp_nudg_prof_cas(k2)-temp_nudg_prof_cas(k1))
930	th_nudg_mod_cas(l)= th_nudg_prof_cas(k2) - frac*(th_nudg_prof_cas(k2)-th_nudg_prof_cas(k1))
931	thl_nudg_mod_cas(l)= thl_nudg_prof_cas(k2) - frac*(thl_nudg_prof_cas(k2)-thl_nudg_prof_cas(k1))
932	qt_nudg_mod_cas(l)= qt_nudg_prof_cas(k2) - frac*(qt_nudg_prof_cas(k2)-qt_nudg_prof_cas(k1))
933	qv_nudg_mod_cas(l)= qv_nudg_prof_cas(k2) - frac*(qv_nudg_prof_cas(k2)-qv_nudg_prof_cas(k1))
934	rt_nudg_mod_cas(l)= rt_nudg_prof_cas(k2) - frac*(rt_nudg_prof_cas(k2)-rt_nudg_prof_cas(k1))
935	rv_nudg_mod_cas(l)= rv_nudg_prof_cas(k2) - frac*(rv_nudg_prof_cas(k2)-rv_nudg_prof_cas(k1))
936	u_nudg_mod_cas(l)= u_nudg_prof_cas(k2) - frac*(u_nudg_prof_cas(k2)-u_nudg_prof_cas(k1))
937	v_nudg_mod_cas(l)= v_nudg_prof_cas(k2) - frac*(v_nudg_prof_cas(k2)-v_nudg_prof_cas(k1))
938	uadv_mod_cas(l)= uadv_prof_cas(k2) - frac*(uadv_prof_cas(k2)-uadv_prof_cas(k1))
939	vadv_mod_cas(l)= vadv_prof_cas(k2) - frac*(vadv_prof_cas(k2)-vadv_prof_cas(k1))
940	tadv_mod_cas(l)= tadv_prof_cas(k2) - frac*(tadv_prof_cas(k2)-tadv_prof_cas(k1))
941	thadv_mod_cas(l)= thadv_prof_cas(k2) - frac*(thadv_prof_cas(k2)-thadv_prof_cas(k1))
942	thladv_mod_cas(l)= thladv_prof_cas(k2) - frac*(thladv_prof_cas(k2)-thladv_prof_cas(k1))
943	qtadv_mod_cas(l)= qtadv_prof_cas(k2) - frac*(qtadv_prof_cas(k2)-qtadv_prof_cas(k1))
944	qvadv_mod_cas(l)= qvadv_prof_cas(k2) - frac*(qvadv_prof_cas(k2)-qvadv_prof_cas(k1))
945	rtadv_mod_cas(l)= rtadv_prof_cas(k2) - frac*(rtadv_prof_cas(k2)-rtadv_prof_cas(k1))
946	rvadv_mod_cas(l)= rvadv_prof_cas(k2) - frac*(rvadv_prof_cas(k2)-rvadv_prof_cas(k1))
947	trad_mod_cas(l)= trad_prof_cas(k2) - frac*(trad_prof_cas(k2)-trad_prof_cas(k1))
948	thrad_mod_cas(l)= thrad_prof_cas(k2) - frac*(thrad_prof_cas(k2)-thrad_prof_cas(k1))
949	thlrad_mod_cas(l)= thlrad_prof_cas(k2) - frac*(thlrad_prof_cas(k2)-thlrad_prof_cas(k1))
950
951	else !play>plev_prof_cas(1)
952
953	k1=1
954	k2=2
955	print *,'interp_vert, k1,k2=',plev_prof_cas(k1),plev_prof_cas(k2)
956	frac1 = (play(l)-plev_prof_cas(k2))/(plev_prof_cas(k1)-plev_prof_cas(k2))
957	frac2 = (play(l)-plev_prof_cas(k1))/(plev_prof_cas(k1)-plev_prof_cas(k2))
958	t_mod_cas(l)= frac1t_prof_cas(k1) - frac2t_prof_cas(k2)
959	th_mod_cas(l)= frac1th_prof_cas(k1) - frac2th_prof_cas(k2)
960	if(th_mod_cas(l).NE.0) t_mod_cas(l)= th_mod_cas(l)(play(l)/100000.)*(RD/RCPD)
961	thl_mod_cas(l)= frac1thl_prof_cas(k1) - frac2thl_prof_cas(k2)
962	qt_mod_cas(l)= frac1qt_prof_cas(k1) - frac2qt_prof_cas(k2)
963	qv_mod_cas(l)= frac1qv_prof_cas(k1) - frac2qv_prof_cas(k2)
964	ql_mod_cas(l)= frac1ql_prof_cas(k1) - frac2ql_prof_cas(k2)
965	qi_mod_cas(l)= frac1qi_prof_cas(k1) - frac2qi_prof_cas(k2)
966	rt_mod_cas(l)= frac1rt_prof_cas(k1) - frac2rt_prof_cas(k2)
967	rv_mod_cas(l)= frac1rv_prof_cas(k1) - frac2rv_prof_cas(k2)
968	rl_mod_cas(l)= frac1rl_prof_cas(k1) - frac2rl_prof_cas(k2)
969	ri_mod_cas(l)= frac1ri_prof_cas(k1) - frac2ri_prof_cas(k2)
970	rh_mod_cas(l)= frac1rh_prof_cas(k1) - frac2rh_prof_cas(k2)
971	u_mod_cas(l)= frac1u_prof_cas(k1) - frac2u_prof_cas(k2)
972	v_mod_cas(l)= frac1v_prof_cas(k1) - frac2v_prof_cas(k2)
973	w_mod_cas(l)= frac1w_prof_cas(k1) - frac2w_prof_cas(k2)
974	omega_mod_cas(l)= frac1omega_prof_cas(k1) - frac2omega_prof_cas(k2)
975	ug_mod_cas(l)= frac1ug_prof_cas(k1) - frac2ug_prof_cas(k2)
976	vg_mod_cas(l)= frac1vg_prof_cas(k1) - frac2vg_prof_cas(k2)
977	temp_nudg_mod_cas(l)= frac1temp_nudg_prof_cas(k1) - frac2temp_nudg_prof_cas(k2)
978	th_nudg_mod_cas(l)= frac1th_nudg_prof_cas(k1) - frac2th_nudg_prof_cas(k2)
979	thl_nudg_mod_cas(l)= frac1thl_nudg_prof_cas(k1) - frac2thl_nudg_prof_cas(k2)
980	qt_nudg_mod_cas(l)= frac1qt_nudg_prof_cas(k1) - frac2qt_nudg_prof_cas(k2)
981	qv_nudg_mod_cas(l)= frac1qv_nudg_prof_cas(k1) - frac2qv_nudg_prof_cas(k2)
982	rt_nudg_mod_cas(l)= frac1rt_nudg_prof_cas(k1) - frac2rt_nudg_prof_cas(k2)
983	rv_nudg_mod_cas(l)= frac1rv_nudg_prof_cas(k1) - frac2rv_nudg_prof_cas(k2)
984	u_nudg_mod_cas(l)= frac1u_nudg_prof_cas(k1) - frac2u_nudg_prof_cas(k2)
985	v_nudg_mod_cas(l)= frac1v_nudg_prof_cas(k1) - frac2v_nudg_prof_cas(k2)
986	uadv_mod_cas(l)= frac1uadv_prof_cas(k1) - frac2uadv_prof_cas(k2)
987	vadv_mod_cas(l)= frac1vadv_prof_cas(k1) - frac2vadv_prof_cas(k2)
988	tadv_mod_cas(l)= frac1tadv_prof_cas(k1) - frac2tadv_prof_cas(k2)
989	thadv_mod_cas(l)= frac1thadv_prof_cas(k1) - frac2thadv_prof_cas(k2)
990	thladv_mod_cas(l)= frac1thladv_prof_cas(k1) - frac2thladv_prof_cas(k2)
991	qtadv_mod_cas(l)= frac1qtadv_prof_cas(k1) - frac2qtadv_prof_cas(k2)
992	qvadv_mod_cas(l)= frac1qvadv_prof_cas(k1) - frac2qvadv_prof_cas(k2)
993	rtadv_mod_cas(l)= frac1rtadv_prof_cas(k1) - frac2rtadv_prof_cas(k2)
994	rvadv_mod_cas(l)= frac1rvadv_prof_cas(k1) - frac2rvadv_prof_cas(k2)
995	trad_mod_cas(l)= frac1trad_prof_cas(k1) - frac2trad_prof_cas(k2)
996	thrad_mod_cas(l)= frac1thrad_prof_cas(k1) - frac2thrad_prof_cas(k2)
997	thlrad_mod_cas(l)= frac1thlrad_prof_cas(k1) - frac2thlrad_prof_cas(k2)
998
999	endif ! play.le.plev_prof_cas(1)
1000
1001	else ! above max altitude of forcing file
1002
1003	!jyg
1004	fact=20.*(plev_prof_cas(nlev_cas)-play(l))/plev_prof_cas(nlev_cas) !jyg
1005	fact = max(fact,0.) !jyg
1006	fact = exp(-fact) !jyg
1007	t_mod_cas(l)= t_prof_cas(nlev_cas) !jyg
1008	th_mod_cas(l)= th_prof_cas(nlev_cas) !jyg
1009	thl_mod_cas(l)= thl_prof_cas(nlev_cas) !jyg
1010	qt_mod_cas(l)= qt_prof_cas(nlev_cas)*fact !jyg
1011	qv_mod_cas(l)= qv_prof_cas(nlev_cas)*fact !jyg
1012	ql_mod_cas(l)= ql_prof_cas(nlev_cas)*fact !jyg
1013	qi_mod_cas(l)= qi_prof_cas(nlev_cas)*fact !jyg
1014	rt_mod_cas(l)= rt_prof_cas(nlev_cas)*fact !jyg
1015	rv_mod_cas(l)= rv_prof_cas(nlev_cas)*fact !jyg
1016	rl_mod_cas(l)= rl_prof_cas(nlev_cas)*fact !jyg
1017	ri_mod_cas(l)= ri_prof_cas(nlev_cas)*fact !jyg
1018	rh_mod_cas(l)= rh_prof_cas(nlev_cas)*fact !jyg
1019	u_mod_cas(l)= u_prof_cas(nlev_cas)*fact !jyg
1020	v_mod_cas(l)= v_prof_cas(nlev_cas)*fact !jyg
1021	w_mod_cas(l)= 0.0 !jyg
1022	omega_mod_cas(l)= 0.0 !jyg
1023	ug_mod_cas(l)= ug_prof_cas(nlev_cas) !jyg
1024	vg_mod_cas(l)= vg_prof_cas(nlev_cas) !jyg
1025	temp_nudg_mod_cas(l)= temp_nudg_prof_cas(nlev_cas) !jyg
1026	th_nudg_mod_cas(l)= th_nudg_prof_cas(nlev_cas) !jyg
1027	thl_nudg_mod_cas(l)= thl_nudg_prof_cas(nlev_cas) !jyg
1028	qt_nudg_mod_cas(l)= qt_nudg_prof_cas(nlev_cas) !jyg
1029	qv_nudg_mod_cas(l)= qv_nudg_prof_cas(nlev_cas) !jyg
1030	rt_nudg_mod_cas(l)= rt_nudg_prof_cas(nlev_cas) !jyg
1031	rv_nudg_mod_cas(l)= rv_nudg_prof_cas(nlev_cas) !jyg
1032	u_nudg_mod_cas(l)= u_nudg_prof_cas(nlev_cas) !jyg
1033	v_nudg_mod_cas(l)= v_nudg_prof_cas(nlev_cas) !jyg
1034	uadv_mod_cas(l)= uadv_prof_cas(nlev_cas) !jyg
1035	vadv_mod_cas(l)= vadv_prof_cas(nlev_cas) !jyg
1036	tadv_mod_cas(l)= tadv_prof_cas(nlev_cas) !jyg
1037	thadv_mod_cas(l)= thadv_prof_cas(nlev_cas) !jyg
1038	thladv_mod_cas(l)= thladv_prof_cas(nlev_cas) !jyg
1039	qtadv_mod_cas(l)= qtadv_prof_cas(nlev_cas) !jyg
1040	qvadv_mod_cas(l)= qvadv_prof_cas(nlev_cas) !jyg
1041	rtadv_mod_cas(l)= rtadv_prof_cas(nlev_cas) !jyg
1042	rvadv_mod_cas(l)= rvadv_prof_cas(nlev_cas) !jyg
1043	trad_mod_cas(l)= trad_prof_cas(nlev_cas)*fact !jyg
1044	thrad_mod_cas(l)= thrad_prof_cas(nlev_cas)*fact !jyg
1045	thlrad_mod_cas(l)= thlrad_prof_cas(nlev_cas)*fact !jyg
1046
1047	endif ! play
1048
1049	enddo ! l
1050
1051	return
1052	end SUBROUTINE interp_case_vertical_std
1053	!*****************************************************************************
1054
1055
1056
1057
1058
1059	END MODULE mod_1D_cases_read_std

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

Download in other formats:

Original Format