source: LMDZ6/trunk/libf/phylmd/dyn1d/mod_1D_cases_read_std.F90 @ 3575

Last change on this file since 3575 was 3541, checked in by fhourdin, 5 years ago

Gros nettoyage en cours sur le 1D.
Le nouveau lmdz1d.F90 est une coquille vide qui choisit entre
old_lmdz1d.F90 (l'ancien lmdz1d.F90) et scm.F90 (le nouveau au format standard).
Plusieur fichiers sont renommés old_truc, le truc étant au format standard,
nettoyé des anciens formats.
Le 1DUTILS.h est lui même séparé entre des routines génériques venant remplacer
notamment des routines de dyn3d (la vocation d'origine de 1DUTILS.h) et
les routiles de lecture spécifiques mises dans old_1DUTILS.h
On perdra un peu de l'utilister de svn au moment de cette grosse bascule.
Mais les old_ sont faits pour ne plus bouger, et les versions standard
sont en pleine évolution.
Fredho

File size: 45.4 KB
Line 
1!
2! $Id: mod_1D_cases_read.F90 2373 2015-10-13 17:28:01Z jyg $
3!
4MODULE 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
12
13!profils environnementaux
14        real, allocatable::  plev_cas(:,:),plevh_cas(:)
15        real, allocatable::  ap_cas(:),bp_cas(:)
16
17        real, allocatable::  z_cas(:,:),zh_cas(:)
18        real, allocatable::  t_cas(:,:),q_cas(:,:),qv_cas(:,:),ql_cas(:,:),qi_cas(:,:),rh_cas(:,:)
19        real, allocatable::  th_cas(:,:),thv_cas(:,:),thl_cas(:,:),rv_cas(:,:)
20        real, allocatable::  u_cas(:,:),v_cas(:,:),vitw_cas(:,:),omega_cas(:,:)
21
22!forcing
23        real, allocatable::  ht_cas(:,:),vt_cas(:,:),dt_cas(:,:),dtrad_cas(:,:)
24        real, allocatable::  hth_cas(:,:),vth_cas(:,:),dth_cas(:,:)
25        real, allocatable::  hq_cas(:,:),vq_cas(:,:),dq_cas(:,:)
26        real, allocatable::  hr_cas(:,:),vr_cas(:,:),dr_cas(:,:)
27        real, allocatable::  hu_cas(:,:),vu_cas(:,:),du_cas(:,:)
28        real, allocatable::  hv_cas(:,:),vv_cas(:,:),dv_cas(:,:)
29        real, allocatable::  ug_cas(:,:),vg_cas(:,:)
30        real, allocatable::  lat_cas(:),sens_cas(:),ts_cas(:),ps_cas(:),ustar_cas(:)
31        real, allocatable::  uw_cas(:,:),vw_cas(:,:),q1_cas(:,:),q2_cas(:,:),tke_cas(:)
32
33!champs interpoles
34        real, allocatable::  plev_prof_cas(:)
35        real, allocatable::  t_prof_cas(:)
36        real, allocatable::  theta_prof_cas(:)
37        real, allocatable::  thl_prof_cas(:)
38        real, allocatable::  thv_prof_cas(:)
39        real, allocatable::  q_prof_cas(:)
40        real, allocatable::  qv_prof_cas(:)
41        real, allocatable::  ql_prof_cas(:)
42        real, allocatable::  qi_prof_cas(:)
43        real, allocatable::  rh_prof_cas(:)
44        real, allocatable::  rv_prof_cas(:)
45        real, allocatable::  u_prof_cas(:)
46        real, allocatable::  v_prof_cas(:)       
47        real, allocatable::  vitw_prof_cas(:)
48        real, allocatable::  omega_prof_cas(:)
49        real, allocatable::  ug_prof_cas(:)
50        real, allocatable::  vg_prof_cas(:)
51        real, allocatable::  ht_prof_cas(:)
52        real, allocatable::  hth_prof_cas(:)
53        real, allocatable::  hq_prof_cas(:)
54        real, allocatable::  vt_prof_cas(:)
55        real, allocatable::  vth_prof_cas(:)
56        real, allocatable::  vq_prof_cas(:)
57        real, allocatable::  dt_prof_cas(:)
58        real, allocatable::  dth_prof_cas(:)
59        real, allocatable::  dtrad_prof_cas(:)
60        real, allocatable::  dq_prof_cas(:)
61        real, allocatable::  hu_prof_cas(:)
62        real, allocatable::  hv_prof_cas(:)
63        real, allocatable::  vu_prof_cas(:)
64        real, allocatable::  vv_prof_cas(:)
65        real, allocatable::  du_prof_cas(:)
66        real, allocatable::  dv_prof_cas(:)
67        real, allocatable::  uw_prof_cas(:)
68        real, allocatable::  vw_prof_cas(:)
69        real, allocatable::  q1_prof_cas(:)
70        real, allocatable::  q2_prof_cas(:)
71
72
73        real lat_prof_cas,sens_prof_cas,ts_prof_cas,ps_prof_cas,ustar_prof_cas,tke_prof_cas
74        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
75     
76
77
78CONTAINS
79
80
81!**********************************************************************************************
82SUBROUTINE read_SCM_cas
83      implicit none
84
85#include "netcdf.inc"
86#include "date_cas.h"
87
88      INTEGER nid,rid,ierr
89      INTEGER ii,jj,timeid
90      REAL, ALLOCATABLE :: time_val(:)
91
92      print*,'ON EST VRAIMENT LA'
93      fich_cas='cas.nc'
94      print*,'fich_cas ',fich_cas
95      ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid)
96      print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid
97      if (ierr.NE.NF_NOERR) then
98         write(*,*) 'ERROR: GROS Pb opening forcings nc file '
99         write(*,*) NF_STRERROR(ierr)
100         stop ""
101      endif
102!.......................................................................
103      ierr=NF_INQ_DIMID(nid,'lat',rid)
104      IF (ierr.NE.NF_NOERR) THEN
105         print*, 'Oh probleme lecture dimension lat'
106      ENDIF
107      ierr=NF_INQ_DIMLEN(nid,rid,ii)
108      print*,'OK1 read2: nid,rid,lat',nid,rid,ii
109!.......................................................................
110      ierr=NF_INQ_DIMID(nid,'lon',rid)
111      IF (ierr.NE.NF_NOERR) THEN
112         print*, 'Oh probleme lecture dimension lon'
113      ENDIF
114      ierr=NF_INQ_DIMLEN(nid,rid,jj)
115      print*,'OK2 read2: nid,rid,lat',nid,rid,jj
116!.......................................................................
117      ierr=NF_INQ_DIMID(nid,'lev',rid)
118      IF (ierr.NE.NF_NOERR) THEN
119         print*, 'Oh probleme lecture dimension nlev'
120      ENDIF
121      ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas)
122      print*,'OK3 read2: nid,rid,nlev_cas',nid,rid,nlev_cas
123      IF ( .NOT. ( nlev_cas > 10 .AND. nlev_cas < 1000 )) THEN
124              print*,'Valeur de nlev_cas peu probable'
125              STOP
126      ENDIF
127!.......................................................................
128      ierr=NF_INQ_DIMID(nid,'time',rid)
129      nt_cas=0
130      IF (ierr.NE.NF_NOERR) THEN
131        stop 'Oh probleme lecture dimension time'
132      ENDIF
133      ierr=NF_INQ_DIMLEN(nid,rid,nt_cas)
134      print*,'OK4 read2: nid,rid,nt_cas',nid,rid,nt_cas
135! Lecture de l'axe des temps
136      print*,'LECTURE DU TEMPS'
137      ierr=NF_INQ_VARID(nid,'time',timeid)
138         if(ierr/=NF_NOERR) then
139           print *,'Variable time manquante dans cas.nc:'
140           ierr=NF_NOERR
141         else
142                 allocate(time_val(nt_cas))
143#ifdef NC_DOUBLE
144         ierr = NF_GET_VAR_DOUBLE(nid,timeid,time_val)
145#else
146           ierr = NF_GET_VAR_REAL(nid,timeid,time_val)
147#endif
148           if(ierr/=NF_NOERR) then
149              print *,'Pb a la lecture de time cas.nc: '
150           endif
151   endif
152   IF (nt_cas>1) THEN
153           pdt_cas=time_val(2)-time_val(1)
154   ELSE
155           pdt_cas=0.
156   ENDIF
157
158
159!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
160!profils moyens:
161        allocate(plev_cas(nlev_cas,nt_cas),plevh_cas(nlev_cas+1))       
162        allocate(z_cas(nlev_cas,nt_cas),zh_cas(nlev_cas+1))
163        allocate(ap_cas(nlev_cas+1),bp_cas(nt_cas+1))
164        allocate(t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas),qv_cas(nlev_cas,nt_cas),ql_cas(nlev_cas,nt_cas), &
165             qi_cas(nlev_cas,nt_cas),rh_cas(nlev_cas,nt_cas))
166        allocate(th_cas(nlev_cas,nt_cas),thl_cas(nlev_cas,nt_cas),thv_cas(nlev_cas,nt_cas),rv_cas(nlev_cas,nt_cas))
167        allocate(u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas),vitw_cas(nlev_cas,nt_cas),omega_cas(nlev_cas,nt_cas))
168
169!forcing
170        allocate(ht_cas(nlev_cas,nt_cas),vt_cas(nlev_cas,nt_cas),dt_cas(nlev_cas,nt_cas),dtrad_cas(nlev_cas,nt_cas))
171        allocate(hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas),dq_cas(nlev_cas,nt_cas))
172        allocate(hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas),dth_cas(nlev_cas,nt_cas))
173        allocate(hr_cas(nlev_cas,nt_cas),vr_cas(nlev_cas,nt_cas),dr_cas(nlev_cas,nt_cas))
174        allocate(hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas),du_cas(nlev_cas,nt_cas))
175        allocate(hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas),dv_cas(nlev_cas,nt_cas))
176        allocate(ug_cas(nlev_cas,nt_cas))
177        allocate(vg_cas(nlev_cas,nt_cas))
178        allocate(lat_cas(nt_cas),sens_cas(nt_cas),ts_cas(nt_cas),ps_cas(nt_cas),ustar_cas(nt_cas),tke_cas(nt_cas))
179        allocate(uw_cas(nlev_cas,nt_cas),vw_cas(nlev_cas,nt_cas),q1_cas(nlev_cas,nt_cas),q2_cas(nlev_cas,nt_cas))
180
181
182
183!champs interpoles
184        allocate(plev_prof_cas(nlev_cas))
185        allocate(t_prof_cas(nlev_cas))
186        allocate(theta_prof_cas(nlev_cas))
187        allocate(thl_prof_cas(nlev_cas))
188        allocate(thv_prof_cas(nlev_cas))
189        allocate(q_prof_cas(nlev_cas))
190        allocate(qv_prof_cas(nlev_cas))
191        allocate(ql_prof_cas(nlev_cas))
192        allocate(qi_prof_cas(nlev_cas))
193        allocate(rh_prof_cas(nlev_cas))
194        allocate(rv_prof_cas(nlev_cas))
195        allocate(u_prof_cas(nlev_cas))
196        allocate(v_prof_cas(nlev_cas))
197        allocate(vitw_prof_cas(nlev_cas))
198        allocate(omega_prof_cas(nlev_cas))
199        allocate(ug_prof_cas(nlev_cas))
200        allocate(vg_prof_cas(nlev_cas))
201        allocate(ht_prof_cas(nlev_cas))
202        allocate(hth_prof_cas(nlev_cas))
203        allocate(hq_prof_cas(nlev_cas))
204        allocate(hu_prof_cas(nlev_cas))
205        allocate(hv_prof_cas(nlev_cas))
206        allocate(vt_prof_cas(nlev_cas))
207        allocate(vth_prof_cas(nlev_cas))
208        allocate(vq_prof_cas(nlev_cas))
209        allocate(vu_prof_cas(nlev_cas))
210        allocate(vv_prof_cas(nlev_cas))
211        allocate(dt_prof_cas(nlev_cas))
212        allocate(dth_prof_cas(nlev_cas))
213        allocate(dtrad_prof_cas(nlev_cas))
214        allocate(dq_prof_cas(nlev_cas))
215        allocate(du_prof_cas(nlev_cas))
216        allocate(dv_prof_cas(nlev_cas))
217        allocate(uw_prof_cas(nlev_cas))
218        allocate(vw_prof_cas(nlev_cas))
219        allocate(q1_prof_cas(nlev_cas))
220        allocate(q2_prof_cas(nlev_cas))
221
222        print*,'Allocations OK'
223        call read_SCM (nid,nlev_cas,nt_cas,                                                                     &
224     &     ap_cas,bp_cas,z_cas,plev_cas,zh_cas,plevh_cas,t_cas,th_cas,thv_cas,thl_cas,qv_cas,                    &
225     &     ql_cas,qi_cas,rh_cas,rv_cas,u_cas,v_cas,vitw_cas,omega_cas,ug_cas,vg_cas,du_cas,hu_cas,vu_cas,        &
226     &     dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas,dth_cas,hth_cas,vth_cas,               &
227     &     dr_cas,hr_cas,vr_cas,dtrad_cas,sens_cas,lat_cas,ts_cas,ps_cas,ustar_cas,tke_cas,                      &
228     &     uw_cas,vw_cas,q1_cas,q2_cas,orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,heat_rough, &
229     &     o3_cas,rugos_cas,clay_cas,sand_cas)
230        print*,'Read2 cas OK'
231        do ii=1,nlev_cas
232        print*,'apres read2_cas, plev_cas=',ii,plev_cas(ii,1)
233        enddo
234
235
236END SUBROUTINE read_SCM_cas
237
238
239!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
240SUBROUTINE deallocate2_1D_cases
241!profils environnementaux:
242        deallocate(plev_cas,plevh_cas)
243       
244        deallocate(z_cas,zh_cas)
245        deallocate(ap_cas,bp_cas)
246        deallocate(t_cas,q_cas,qv_cas,ql_cas,qi_cas,rh_cas)
247        deallocate(th_cas,thl_cas,thv_cas,rv_cas)
248        deallocate(u_cas,v_cas,vitw_cas,omega_cas)
249       
250!forcing
251        deallocate(ht_cas,vt_cas,dt_cas,dtrad_cas)
252        deallocate(hq_cas,vq_cas,dq_cas)
253        deallocate(hth_cas,vth_cas,dth_cas)
254        deallocate(hr_cas,vr_cas,dr_cas)
255        deallocate(hu_cas,vu_cas,du_cas)
256        deallocate(hv_cas,vv_cas,dv_cas)
257        deallocate(ug_cas)
258        deallocate(vg_cas)
259        deallocate(lat_cas,sens_cas,ts_cas,ps_cas,ustar_cas,tke_cas,uw_cas,vw_cas,q1_cas,q2_cas)
260
261!champs interpoles
262        deallocate(plev_prof_cas)
263        deallocate(t_prof_cas)
264        deallocate(theta_prof_cas)
265        deallocate(thl_prof_cas)
266        deallocate(thv_prof_cas)
267        deallocate(q_prof_cas)
268        deallocate(qv_prof_cas)
269        deallocate(ql_prof_cas)
270        deallocate(qi_prof_cas)
271        deallocate(rh_prof_cas)
272        deallocate(rv_prof_cas)
273        deallocate(u_prof_cas)
274        deallocate(v_prof_cas)
275        deallocate(vitw_prof_cas)
276        deallocate(omega_prof_cas)
277        deallocate(ug_prof_cas)
278        deallocate(vg_prof_cas)
279        deallocate(ht_prof_cas)
280        deallocate(hq_prof_cas)
281        deallocate(hu_prof_cas)
282        deallocate(hv_prof_cas)
283        deallocate(vt_prof_cas)
284        deallocate(vq_prof_cas)
285        deallocate(vu_prof_cas)
286        deallocate(vv_prof_cas)
287        deallocate(dt_prof_cas)
288        deallocate(dtrad_prof_cas)
289        deallocate(dq_prof_cas)
290        deallocate(du_prof_cas)
291        deallocate(dv_prof_cas)
292        deallocate(t_prof_cas)
293        deallocate(u_prof_cas)
294        deallocate(v_prof_cas)
295        deallocate(uw_prof_cas)
296        deallocate(vw_prof_cas)
297        deallocate(q1_prof_cas)
298        deallocate(q2_prof_cas)
299
300END SUBROUTINE deallocate2_1D_cases
301
302
303!=====================================================================
304      subroutine read_cas2(nid,nlevel,ntime                          &
305     &     ,zz,pp,temp,qv,rh,theta,rv,u,v,ug,vg,w,                   &
306     &     du,hu,vu,dv,hv,vv,dt,dtrad,ht,vt,dq,hq,vq,                     &
307     &     dth,hth,vth,dr,hr,vr,sens,flat,ts,ustar,uw,vw,q1,q2)
308
309!program reading forcing of the case study
310      implicit none
311#include "netcdf.inc"
312
313      integer ntime,nlevel
314
315      real zz(nlevel,ntime)
316      real pp(nlevel,ntime)
317      real temp(nlevel,ntime),qv(nlevel,ntime),rh(nlevel,ntime)
318      real theta(nlevel,ntime),rv(nlevel,ntime)
319      real u(nlevel,ntime)
320      real v(nlevel,ntime)
321      real ug(nlevel,ntime)
322      real vg(nlevel,ntime)
323      real w(nlevel,ntime)
324      real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
325      real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
326      real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
327      real dtrad(nlevel,ntime)
328      real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
329      real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime)
330      real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
331      real flat(ntime),sens(ntime),ts(ntime),ustar(ntime)
332      real uw(nlevel,ntime),vw(nlevel,ntime),q1(nlevel,ntime),q2(nlevel,ntime),resul(nlevel,ntime),resul1(ntime)
333
334
335      integer nid, ierr, ierr1,ierr2,rid,i
336      integer nbvar3d
337      parameter(nbvar3d=39)
338      integer var3didin(nbvar3d)
339      character*5 name_var(1:nbvar3d)
340      data name_var/'zz','pp','temp','qv','rh','theta','rv','u','v','ug','vg','w','advu','hu','vu',&
341     &'advv','hv','vv','advT','hT','vT','advq','hq','vq','advth','hth','vth','advr','hr','vr',&
342     &'radT','uw','vw','q1','q2','sens','flat','ts','ustar'/
343
344       do i=1,nbvar3d
345         print *,'Dans read_cas2, on va lire ',nid,i,name_var(i)
346       enddo
347       do i=1,nbvar3d
348         ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
349         print *,'ierr=',i,ierr,name_var(i),var3didin(i)
350         if(ierr/=NF_NOERR) then
351           print *,'Variable manquante dans cas.nc:',name_var(i)
352         endif
353       enddo
354       do i=1,nbvar3d
355         print *,'Dans read_cas2, on va lire ',var3didin(i),name_var(i)
356         if(i.LE.35) then
357#ifdef NC_DOUBLE
358         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
359#else
360         ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
361#endif
362         print *,'Dans read_cas2, on a lu ',ierr,var3didin(i),name_var(i)
363         if(ierr/=NF_NOERR) then
364            print *,'Pb a la lecture de cas.nc: ',name_var(i)
365            stop "getvarup"
366         endif
367         else
368#ifdef NC_DOUBLE
369         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul1)
370#else
371         ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul1)
372#endif
373         print *,'Dans read_cas2, on a lu ',ierr,var3didin(i),name_var(i)
374         if(ierr/=NF_NOERR) then
375            print *,'Pb a la lecture de cas.nc: ',name_var(i)
376            stop "getvarup"
377         endif
378         endif
379         select case(i)
380           case(1) ; zz=resul
381           case(2) ; pp=resul
382           case(3) ; temp=resul
383           case(4) ; qv=resul
384           case(5) ; rh=resul
385           case(6) ; theta=resul
386           case(7) ; rv=resul
387           case(8) ; u=resul
388           case(9) ; v=resul
389           case(10) ; ug=resul
390           case(11) ; vg=resul
391           case(12) ; w=resul
392           case(13) ; du=resul
393           case(14) ; hu=resul
394           case(15) ; vu=resul
395           case(16) ; dv=resul
396           case(17) ; hv=resul
397           case(18) ; vv=resul
398           case(19) ; dt=resul
399           case(20) ; ht=resul
400           case(21) ; vt=resul
401           case(22) ; dq=resul
402           case(23) ; hq=resul
403           case(24) ; vq=resul
404           case(25) ; dth=resul
405           case(26) ; hth=resul
406           case(27) ; vth=resul
407           case(28) ; dr=resul
408           case(29) ; hr=resul
409           case(30) ; vr=resul
410           case(31) ; dtrad=resul
411           case(32) ; uw=resul
412           case(33) ; vw=resul
413           case(34) ; q1=resul
414           case(35) ; q2=resul
415           case(36) ; sens=resul1
416           case(37) ; flat=resul1
417           case(38) ; ts=resul1
418           case(39) ; ustar=resul1
419         end select
420       enddo
421
422         return
423         end subroutine read_cas2
424!======================================================================
425      subroutine read2_cas(nid,nlevel,ntime,                                       &
426     &     ap,bp,zz,pp,zzh,pph,temp,theta,thv,thl,qv,ql,qi,rh,rv,u,v,vitw,omega,ug,vg,&
427     &     du,hu,vu,dv,hv,vv,dt,ht,vt,dq,hq,vq,                                    &
428     &     dth,hth,vth,dr,hr,vr,dtrad,sens,flat,ts,ps,ustar,tke,uw,vw,q1,q2,       &
429     &     orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,          &
430     &     heat_rough,o3_cas,rugos_cas,clay_cas,sand_cas)
431
432!program reading forcing of the case study
433      implicit none
434#include "netcdf.inc"
435
436      integer ntime,nlevel
437
438      real ap(nlevel+1),bp(nlevel+1)
439      real zz(nlevel,ntime),zzh(nlevel+1)
440      real pp(nlevel,ntime),pph(nlevel+1)
441      real temp(nlevel,ntime),qv(nlevel,ntime),ql(nlevel,ntime),qi(nlevel,ntime),rh(nlevel,ntime)
442      real theta(nlevel,ntime),thv(nlevel,ntime),thl(nlevel,ntime),rv(nlevel,ntime)
443      real u(nlevel,ntime),v(nlevel,ntime)
444      real ug(nlevel,ntime),vg(nlevel,ntime)
445      real vitw(nlevel,ntime),omega(nlevel,ntime)
446      real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
447      real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
448      real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
449      real dtrad(nlevel,ntime)
450      real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
451      real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime),hthl(nlevel,ntime)
452      real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
453      real flat(ntime),sens(ntime),ustar(ntime)
454      real uw(nlevel,ntime),vw(nlevel,ntime),q1(nlevel,ntime),q2(nlevel,ntime)
455      real ts(ntime),ps(ntime),tke(ntime)
456      real orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,heat_rough,o3_cas,rugos_cas,clay_cas,sand_cas
457      real apbp(nlevel+1),resul(nlevel,ntime),resul1(nlevel),resul2(ntime),resul3
458
459
460      integer nid, ierr,ierr1,ierr2,rid,i
461      integer nbvar3d
462      parameter(nbvar3d=62)
463      integer var3didin(nbvar3d),missing_var(nbvar3d)
464      character*12 name_var(1:nbvar3d)
465      data name_var/'coor_par_a','coor_par_b','height_h','pressure_h',&
466     &'w','omega','ug','vg','uadv','uadvh','uadvv','vadv','vadvh','vadvv','tadv','tadvh','tadvv',&
467     &'qadv','qadvh','qadvv','thadv','thadvh','thadvv','thladvh','radv','radvh','radvv','radcool','q1','q2','ustress','vstress', &
468     'rh',&
469     &'height_f','pressure_f','temp','theta','thv','thl','qv','ql','qi','rv','u','v',&
470     &'sfc_sens_flx','sfc_lat_flx','ts','ps','ustar','tke',&
471     &'orog','albedo','emiss','t_skin','q_skin','mom_rough','heat_rough','o3','rugos','clay','sand'/
472      do i=1,nbvar3d
473        missing_var(i)=0.
474      enddo
475
476!-----------------------------------------------------------------------
477       do i=1,nbvar3d
478         ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
479         if(ierr/=NF_NOERR) then
480           print *,'Variable manquante dans cas.nc:',i,name_var(i)
481           ierr=NF_NOERR
482           missing_var(i)=1
483         else
484!-----------------------------------------------------------------------
485           if(i.LE.4) then     ! Lecture des coord pression en (nlevelp1,lat,lon)
486#ifdef NC_DOUBLE
487           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),apbp)
488#else
489           ierr = NF_GET_VAR_REAL(nid,var3didin(i),apbp)
490#endif
491           print *,'read2_cas(apbp), on a lu ',i,name_var(i)
492           if(ierr/=NF_NOERR) then
493              print *,'Pb a la lecture de cas.nc: ',name_var(i)
494              stop "getvarup"
495           endif
496!-----------------------------------------------------------------------
497           else if(i.gt.4.and.i.LE.45) then   ! Lecture des variables en (time,nlevel,lat,lon)
498#ifdef NC_DOUBLE
499           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
500#else
501           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
502#endif
503           print *,'read2_cas(resul), on a lu ',i,name_var(i)
504           if(ierr/=NF_NOERR) then
505              print *,'Pb a la lecture de cas.nc: ',name_var(i)
506              stop "getvarup"
507           endif
508!-----------------------------------------------------------------------
509           else if (i.gt.45.and.i.LE.51) then   ! Lecture des variables en (time,lat,lon)
510#ifdef NC_DOUBLE
511           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul2)
512#else
513           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul2)
514#endif
515           print *,'read2_cas(resul2), on a lu ',i,name_var(i)
516           if(ierr/=NF_NOERR) then
517              print *,'Pb a la lecture de cas.nc: ',name_var(i)
518              stop "getvarup"
519           endif
520!-----------------------------------------------------------------------
521           else     ! Lecture des constantes (lat,lon)
522#ifdef NC_DOUBLE
523           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul3)
524#else
525           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul3)
526#endif
527           print *,'read2_cas(resul3), on a lu ',i,name_var(i)
528           if(ierr/=NF_NOERR) then
529              print *,'Pb a la lecture de cas.nc: ',name_var(i)
530              stop "getvarup"
531           endif
532           endif
533         endif
534!-----------------------------------------------------------------------
535         select case(i)
536           case(1) ; ap=apbp       ! donnees indexees en nlevel+1
537           case(2) ; bp=apbp
538           case(3) ; zzh=apbp
539           case(4) ; pph=apbp
540           case(5) ; vitw=resul    ! donnees indexees en nlevel,time
541           case(6) ; omega=resul
542           case(7) ; ug=resul
543           case(8) ; vg=resul
544           case(9) ; du=resul
545           case(10) ; hu=resul
546           case(11) ; vu=resul
547           case(12) ; dv=resul
548           case(13) ; hv=resul
549           case(14) ; vv=resul
550           case(15) ; dt=resul
551           case(16) ; ht=resul
552           case(17) ; vt=resul
553           case(18) ; dq=resul
554           case(19) ; hq=resul
555           case(20) ; vq=resul
556           case(21) ; dth=resul
557           case(22) ; hth=resul
558           case(23) ; vth=resul
559           case(24) ; hthl=resul
560           case(25) ; dr=resul
561           case(26) ; hr=resul
562           case(27) ; vr=resul
563           case(28) ; dtrad=resul
564           case(29) ; q1=resul
565           case(30) ; q2=resul
566           case(31) ; uw=resul
567           case(32) ; vw=resul
568           case(33) ; rh=resul
569           case(34) ; zz=resul      ! donnees en time,nlevel pour profil initial
570           case(35) ; pp=resul
571           case(36) ; temp=resul
572           case(37) ; theta=resul
573           case(38) ; thv=resul
574           case(39) ; thl=resul
575           case(40) ; qv=resul
576           case(41) ; ql=resul
577           case(42) ; qi=resul
578           case(43) ; rv=resul
579           case(44) ; u=resul
580           case(45) ; v=resul
581           case(46) ; sens=resul2   ! donnees indexees en time
582           case(47) ; flat=resul2
583           case(48) ; ts=resul2
584           case(49) ; ps=resul2
585           case(50) ; ustar=resul2
586           case(51) ; tke=resul2
587           case(52) ; orog_cas=resul3      ! constantes
588           case(53) ; albedo_cas=resul3
589           case(54) ; emiss_cas=resul3
590           case(55) ; t_skin_cas=resul3
591           case(56) ; q_skin_cas=resul3
592           case(57) ; mom_rough=resul3
593           case(58) ; heat_rough=resul3
594           case(59) ; o3_cas=resul3       
595           case(60) ; rugos_cas=resul3
596           case(61) ; clay_cas=resul3
597           case(62) ; sand_cas=resul3
598         end select
599         resul=0.
600         resul1=0.
601         resul2=0.
602         resul3=0.
603       enddo
604!-----------------------------------------------------------------------
605
606         return
607         end subroutine read2_cas
608
609!======================================================================
610      subroutine read_SCM(nid,nlevel,ntime,                                       &
611     &     ap,bp,zz,pp,zzh,pph,temp,theta,thv,thl,qv,ql,qi,rh,rv,u,v,vitw,omega,ug,vg,&
612     &     du,hu,vu,dv,hv,vv,dt,ht,vt,dq,hq,vq,                                    &
613     &     dth,hth,vth,dr,hr,vr,dtrad,sens,flat,ts,ps,ustar,tke,uw,vw,q1,q2,       &
614     &     orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,          &
615     &     heat_rough,o3_cas,rugos_cas,clay_cas,sand_cas)
616
617!program reading forcing of the case study
618      implicit none
619#include "netcdf.inc"
620
621      integer ntime,nlevel,k,t
622
623      real ap(nlevel+1),bp(nlevel+1)
624      real zz(nlevel,ntime),zzh(nlevel+1)
625      real pp(nlevel,ntime),pph(nlevel+1)
626!profils initiaux
627      real temp0(nlevel),qv0(nlevel),ql0(nlevel),qi0(nlevel),u0(nlevel),v0(nlevel),tke0(nlevel)
628      real pp0(nlevel)   
629      real temp(nlevel,ntime),qv(nlevel,ntime),ql(nlevel,ntime),qi(nlevel,ntime),rh(nlevel,ntime)
630      real theta(nlevel,ntime),thv(nlevel,ntime),thl(nlevel,ntime),rv(nlevel,ntime)
631      real u(nlevel,ntime),v(nlevel,ntime),tke(nlevel,ntime)
632      real ug(nlevel,ntime),vg(nlevel,ntime)
633      real vitw(nlevel,ntime),omega(nlevel,ntime)
634      real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
635      real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
636      real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
637      real dtrad(nlevel,ntime)
638      real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
639      real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime),hthl(nlevel,ntime)
640      real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
641      real flat(ntime),sens(ntime),ustar(ntime)
642      real uw(nlevel,ntime),vw(nlevel,ntime),q1(nlevel,ntime),q2(nlevel,ntime)
643      real ts(ntime),ps(ntime)
644      real orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,heat_rough,o3_cas,rugos_cas,clay_cas,sand_cas
645      real apbp(nlevel+1),resul(nlevel,ntime),resul1(nlevel),resul2(ntime),resul3
646
647
648      integer nid, ierr,ierr1,ierr2,rid,i
649      integer nbvar3d
650      parameter(nbvar3d=70)
651      integer var3didin(nbvar3d),missing_var(nbvar3d)
652      character*13 name_var(1:nbvar3d)
653      data name_var/'coor_par_a','coor_par_b','height_h','pressure_h',&
654     &'temp','qv','ql','qi','u','v','tke','pressure',&
655     &'w','omega','ug','vg','uadv','uadvh','uadvv','vadv','vadvh','vadvv','tadv','tadvh','tadvv',&
656     &'qvadv','qvadvh','qvadvv','thadv','thadvh','thadvv','thladvh','radv','radvh','radvv','radcool','q1','q2','ustress','vstress', &
657     'rh',&
658     &'height_f','pressure_forc','tempt','theta','thv','thl','qvt','qlt','qit','rv','ut','vt','tket',&
659     &'sfc_sens_flx','sfc_lat_flx','ts','ps','ustar',&
660     &'orog','albedo','emiss','t_skin','q_skin','mom_rough','heat_rough','o3','rugos','clay','sand'/
661      do i=1,nbvar3d
662        missing_var(i)=0.
663      enddo
664
665!-----------------------------------------------------------------------
666
667     print*,'ON EST LA'
668       do i=1,nbvar3d
669         ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
670         if(ierr/=NF_NOERR) then
671           print *,'Variable manquante dans cas.nc:',i,name_var(i)
672           ierr=NF_NOERR
673           missing_var(i)=1
674         else
675!-----------------------------------------------------------------------
676           if(i.LE.4) then     ! Lecture des coord pression en (nlevelp1,lat,lon)
677#ifdef NC_DOUBLE
678           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),apbp)
679#else
680           ierr = NF_GET_VAR_REAL(nid,var3didin(i),apbp)
681#endif
682           print *,'read2_cas(apbp), on a lu ',i,name_var(i)
683           if(ierr/=NF_NOERR) then
684              print *,'Pb a la lecture de cas.nc: ',name_var(i)
685              stop "getvarup"
686           endif
687!-----------------------------------------------------------------------
688           else if(i.gt.4.and.i.LE.12) then   ! Lecture des variables en (time,nlevel,lat,lon)
689#ifdef NC_DOUBLE
690           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul1)
691#else
692           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul1)
693#endif
694           print *,'read2_cas(resul1), on a lu ',i,name_var(i)
695           if(ierr/=NF_NOERR) then
696              print *,'Pb a la lecture de cas.nc: ',name_var(i)
697              stop "getvarup"
698           endif
699         print*,'Lecture de la variable #i ',i,name_var(i),minval(resul1),maxval(resul1)
700!-----------------------------------------------------------------------
701           else if(i.gt.12.and.i.LE.54) then   ! Lecture des variables en (time,nlevel,lat,lon)
702#ifdef NC_DOUBLE
703           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
704#else
705           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
706#endif
707           print *,'read2_cas(resul), on a lu ',i,name_var(i)
708           if(ierr/=NF_NOERR) then
709              print *,'Pb a la lecture de cas.nc: ',name_var(i)
710              stop "getvarup"
711           endif
712         print*,'Lecture de la variable #i ',i,name_var(i),minval(resul),maxval(resul)
713!-----------------------------------------------------------------------
714           else if (i.gt.54.and.i.LE.65) then   ! Lecture des variables en (time,lat,lon)
715#ifdef NC_DOUBLE
716           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul2)
717#else
718           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul2)
719#endif
720           print *,'read2_cas(resul2), on a lu ',i,name_var(i)
721           if(ierr/=NF_NOERR) then
722              print *,'Pb a la lecture de cas.nc: ',name_var(i)
723              stop "getvarup"
724           endif
725         print*,'Lecture de la variable #i  ',i,name_var(i),minval(resul2),maxval(resul2)
726!-----------------------------------------------------------------------
727           else     ! Lecture des constantes (lat,lon)
728#ifdef NC_DOUBLE
729           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul3)
730#else
731           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul3)
732#endif
733           print *,'read2_cas(resul3), on a lu ',i,name_var(i)
734           if(ierr/=NF_NOERR) then
735              print *,'Pb a la lecture de cas.nc: ',name_var(i)
736              stop "getvarup"
737           endif
738         print*,'Lecture de la variable #i ',i,name_var(i),resul3
739           endif
740         endif
741!-----------------------------------------------------------------------
742         select case(i)
743         !case(1) ; ap=apbp       ! donnees indexees en nlevel+1
744         ! case(2) ; bp=apbp
745           case(3) ; zzh=apbp
746           case(4) ; pph=apbp
747           case(5) ; temp0=resul1    ! donnees initiales
748           case(6) ; qv0=resul1
749           case(7) ; ql0=resul1
750           case(8) ; qi0=resul1
751           case(9) ; u0=resul1
752           case(10) ; v0=resul1
753           case(11) ; tke0=resul1
754           case(12) ; pp0=resul1
755           case(13) ; vitw=resul    ! donnees indexees en nlevel,time
756           case(14) ; omega=resul
757           case(15) ; ug=resul
758           case(16) ; vg=resul
759           case(17) ; du=resul
760           case(18) ; hu=resul
761           case(19) ; vu=resul
762           case(20) ; dv=resul
763           case(21) ; hv=resul
764           case(22) ; vv=resul
765           case(23) ; dt=resul
766           case(24) ; ht=resul
767           case(25) ; vt=resul
768           case(26) ; dq=resul
769           case(27) ; hq=resul
770           case(28) ; vq=resul
771           case(29) ; dth=resul
772           case(30) ; hth=resul
773           case(31) ; vth=resul
774           case(32) ; hthl=resul
775           case(33) ; dr=resul
776           case(34) ; hr=resul
777           case(35) ; vr=resul
778           case(36) ; dtrad=resul
779           case(37) ; q1=resul
780           case(38) ; q2=resul
781           case(39) ; uw=resul
782           case(40) ; vw=resul
783           case(41) ; rh=resul
784           case(42) ; zz=resul      ! donnees en time,nlevel pour profil initial
785           case(43) ; pp=resul
786           case(44) ; temp=resul
787           case(45) ; theta=resul
788           case(46) ; thv=resul
789           case(47) ; thl=resul
790           case(48) ; qv=resul
791           case(49) ; ql=resul
792           case(50) ; qi=resul
793           case(51) ; rv=resul
794           case(52) ; u=resul
795           case(53) ; v=resul
796           case(54) ; tke=resul
797           case(55) ; sens=resul2   ! donnees indexees en time
798           case(56) ; flat=resul2
799           case(57) ; ts=resul2
800           case(58) ; ps=resul2
801           case(59) ; ustar=resul2
802           case(60) ; orog_cas=resul3      ! constantes
803           case(61) ; albedo_cas=resul3
804           case(62) ; emiss_cas=resul3
805           case(63) ; t_skin_cas=resul3
806           case(64) ; q_skin_cas=resul3
807           case(65) ; mom_rough=resul3
808           case(66) ; heat_rough=resul3
809           case(67) ; o3_cas=resul3       
810           case(68) ; rugos_cas=resul3
811           case(69) ; clay_cas=resul3
812           case(70) ; sand_cas=resul3
813         end select
814         resul=0.
815         resul1=0.
816         resul2=0.
817         resul3=0.
818       enddo
819         print*,'Lecture de la variable APRES ,sens ',minval(sens),maxval(sens)
820         print*,'Lecture de la variable APRES ,flat ',minval(flat),maxval(flat)
821
822!CR:ATTENTION EN ATTENTE DE REGLER LA QUESTION DU PAS DE TEMPS INITIAL
823       do t=1,ntime
824          do k=1,nlevel
825             temp(k,t)=temp0(k)
826             qv(k,t)=qv0(k)
827             ql(k,t)=ql0(k)
828             qi(k,t)=qi0(k)
829             u(k,t)=u0(k)
830             v(k,t)=v0(k)
831             tke(k,t)=tke0(k)
832          enddo
833       enddo
834       !!!! TRAVAIL : EN FONCTION DES DECISIONS SUR LA SPECIFICATION DE W
835       !!!omega=-vitw*pres*rg/(rd*temp)
836!-----------------------------------------------------------------------
837
838         return
839         end subroutine read_SCM
840!======================================================================
841
842!======================================================================
843
844!**********************************************************************************************
845        SUBROUTINE interp_case_time_std(day,day1,annee_ref                           &
846!    &         ,year_cas,day_cas,nt_cas,pdt_forc,nlev_cas                         &
847     &         ,nt_cas,nlev_cas                                                   &
848     &         ,ts_cas,ps_cas,plev_cas,t_cas,theta_cas,thv_cas,thl_cas            &
849     &         ,qv_cas,ql_cas,qi_cas,u_cas,v_cas                                  &
850     &         ,ug_cas,vg_cas,vitw_cas,omega_cas,du_cas,hu_cas,vu_cas             &
851     &         ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas               &
852     &         ,dq_cas,hq_cas,vq_cas,dth_cas,hth_cas,vth_cas                      &
853     &         ,lat_cas,sens_cas,ustar_cas                                        &
854     &         ,uw_cas,vw_cas,q1_cas,q2_cas,tke_cas                               &
855!
856     &         ,ts_prof_cas,plev_prof_cas,t_prof_cas,theta_prof_cas               &
857     &         ,thv_prof_cas,thl_prof_cas,qv_prof_cas,ql_prof_cas,qi_prof_cas     &
858     &         ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas                     &
859     &         ,vitw_prof_cas,omega_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas  &
860     &         ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas                   &
861     &         ,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas                &
862     &         ,hq_prof_cas,vq_prof_cas,dth_prof_cas,hth_prof_cas,vth_prof_cas    &
863     &         ,lat_prof_cas,sens_prof_cas                                        &
864     &         ,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas,tke_prof_cas)
865         
866
867        implicit none
868
869!---------------------------------------------------------------------------------------
870! Time interpolation of a 2D field to the timestep corresponding to day
871!
872! day: current julian day (e.g. 717538.2)
873! day1: first day of the simulation
874! nt_cas: total nb of data in the forcing
875! pdt_cas: total time interval (in sec) between 2 forcing data
876!---------------------------------------------------------------------------------------
877
878#include "compar1d.h"
879#include "date_cas.h"
880
881! inputs:
882        integer annee_ref
883        integer nt_cas,nlev_cas
884        real day, day1,day_cas
885        real ts_cas(nt_cas),ps_cas(nt_cas)
886        real plev_cas(nlev_cas,nt_cas)
887        real t_cas(nlev_cas,nt_cas),theta_cas(nlev_cas,nt_cas),thv_cas(nlev_cas,nt_cas),thl_cas(nlev_cas,nt_cas)
888        real qv_cas(nlev_cas,nt_cas),ql_cas(nlev_cas,nt_cas),qi_cas(nlev_cas,nt_cas)
889        real u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas)
890        real ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas)
891        real vitw_cas(nlev_cas,nt_cas),omega_cas(nlev_cas,nt_cas)
892        real du_cas(nlev_cas,nt_cas),hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas)
893        real dv_cas(nlev_cas,nt_cas),hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas)
894        real dt_cas(nlev_cas,nt_cas),ht_cas(nlev_cas,nt_cas),vt_cas(nlev_cas,nt_cas)
895        real dth_cas(nlev_cas,nt_cas),hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas)
896        real dtrad_cas(nlev_cas,nt_cas)
897        real dq_cas(nlev_cas,nt_cas),hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas)
898        real lat_cas(nt_cas),sens_cas(nt_cas),tke_cas(nt_cas)
899        real ustar_cas(nt_cas),uw_cas(nlev_cas,nt_cas),vw_cas(nlev_cas,nt_cas)
900        real q1_cas(nlev_cas,nt_cas),q2_cas(nlev_cas,nt_cas)
901
902! outputs:
903        real plev_prof_cas(nlev_cas)
904        real t_prof_cas(nlev_cas),theta_prof_cas(nlev_cas),thl_prof_cas(nlev_cas),thv_prof_cas(nlev_cas)
905        real qv_prof_cas(nlev_cas),ql_prof_cas(nlev_cas),qi_prof_cas(nlev_cas)
906        real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas)
907        real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas)
908        real vitw_prof_cas(nlev_cas),omega_prof_cas(nlev_cas)
909        real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas)
910        real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas)
911        real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas)
912        real dth_prof_cas(nlev_cas),hth_prof_cas(nlev_cas),vth_prof_cas(nlev_cas)
913        real dtrad_prof_cas(nlev_cas)
914        real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas)
915        real lat_prof_cas,sens_prof_cas,tke_prof_cas,ts_prof_cas,ustar_prof_cas
916        real uw_prof_cas(nlev_cas),vw_prof_cas(nlev_cas),q1_prof_cas(nlev_cas),q2_prof_cas(nlev_cas)
917! local:
918        integer it_cas1, it_cas2,k
919        real timeit,time_cas1,time_cas2,frac
920
921
922        print*,'Check time',day1,day_ju_ini_cas,day_deb+1,pdt_cas
923!       do k=1,nlev_cas
924!       print*,'debut de interp2_case_time, plev_cas=',k,plev_cas(k,1)
925!       enddo
926
927! On teste si la date du cas AMMA est correcte.
928! C est pour memoire car en fait les fichiers .def
929! sont censes etre corrects.
930! A supprimer a terme (MPL 20150623)
931!     if ((forcing_type.eq.10).and.(1.eq.0)) then
932! Check that initial day of the simulation consistent with AMMA case:
933!      if (annee_ref.ne.2006) then
934!       print*,'Pour AMMA, annee_ref doit etre 2006'
935!       print*,'Changer annee_ref dans run.def'
936!       stop
937!      endif
938!      if (annee_ref.eq.2006 .and. day1.lt.day_cas) then
939!       print*,'AMMA a debute le 10 juillet 2006',day1,day_cas
940!       print*,'Changer dayref dans run.def'
941!       stop
942!      endif
943!      if (annee_ref.eq.2006 .and. day1.gt.day_cas+1) then
944!       print*,'AMMA a fini le 11 juillet'
945!       print*,'Changer dayref ou nday dans run.def'
946!       stop
947!      endif
948!      endif
949
950! Determine timestep relative to the 1st day:
951!       timeit=(day-day1)*86400.
952!       if (annee_ref.eq.1992) then
953!        timeit=(day-day_cas)*86400.
954!       else
955!        timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
956!       endif
957      timeit=(day-day_ju_ini_cas)*86400
958      print *,'day=',day
959      print *,'day_ju_ini_cas=',day_ju_ini_cas
960      print *,'pdt_cas=',pdt_cas
961      print *,'timeit=',timeit
962      print *,'nt_cas=',nt_cas
963
964! Determine the closest observation times:
965!       it_cas1=INT(timeit/pdt_cas)+1
966!       it_cas2=it_cas1 + 1
967!       time_cas1=(it_cas1-1)*pdt_cas
968!       time_cas2=(it_cas2-1)*pdt_cas
969
970       it_cas1=INT(timeit/pdt_cas)+1
971       IF (it_cas1 .EQ. nt_cas) THEN
972       it_cas2=it_cas1
973       ELSE
974       it_cas2=it_cas1 + 1
975       ENDIF
976       time_cas1=(it_cas1-1)*pdt_cas
977       time_cas2=(it_cas2-1)*pdt_cas
978      print *,'timeit,pdt_cas,nt_cas=',timeit,pdt_cas,nt_cas
979      print *,'it_cas1,it_cas2,time_cas1,time_cas2=',it_cas1,it_cas2,time_cas1,time_cas2
980
981       if (it_cas1 .gt. nt_cas) then
982        write(*,*) 'PB-stop: day, day_ju_ini_cas,it_cas1, it_cas2, timeit: '            &
983     &        ,day,day_ju_ini_cas,it_cas1,it_cas2,timeit
984        stop
985       endif
986
987! time interpolation:
988       IF (it_cas1 .EQ. it_cas2) THEN
989          frac=0.
990       ELSE
991          frac=(time_cas2-timeit)/(time_cas2-time_cas1)
992          frac=max(frac,0.0)
993       ENDIF
994
995       lat_prof_cas = lat_cas(it_cas2)                                   &
996     &          -frac*(lat_cas(it_cas2)-lat_cas(it_cas1))
997       sens_prof_cas = sens_cas(it_cas2)                                 &
998     &          -frac*(sens_cas(it_cas2)-sens_cas(it_cas1))
999       tke_prof_cas = tke_cas(it_cas2)                                   &
1000     &          -frac*(tke_cas(it_cas2)-tke_cas(it_cas1))
1001       ts_prof_cas = ts_cas(it_cas2)                                     &
1002     &          -frac*(ts_cas(it_cas2)-ts_cas(it_cas1))
1003       ustar_prof_cas = ustar_cas(it_cas2)                               &
1004     &          -frac*(ustar_cas(it_cas2)-ustar_cas(it_cas1))
1005
1006       do k=1,nlev_cas
1007        plev_prof_cas(k) = plev_cas(k,it_cas2)                           &     
1008     &          -frac*(plev_cas(k,it_cas2)-plev_cas(k,it_cas1))
1009        t_prof_cas(k) = t_cas(k,it_cas2)                                 &       
1010     &          -frac*(t_cas(k,it_cas2)-t_cas(k,it_cas1))
1011        print *,'k,frac,plev_cas1,plev_cas2=',k,frac,plev_cas(k,it_cas1),plev_cas(k,it_cas2)
1012        theta_prof_cas(k) = theta_cas(k,it_cas2)                         &                     
1013     &          -frac*(theta_cas(k,it_cas2)-theta_cas(k,it_cas1))
1014        thv_prof_cas(k) = thv_cas(k,it_cas2)                             &         
1015     &          -frac*(thv_cas(k,it_cas2)-thv_cas(k,it_cas1))
1016        thl_prof_cas(k) = thl_cas(k,it_cas2)                             &             
1017     &          -frac*(thl_cas(k,it_cas2)-thl_cas(k,it_cas1))
1018        qv_prof_cas(k) = qv_cas(k,it_cas2)                               &
1019     &          -frac*(qv_cas(k,it_cas2)-qv_cas(k,it_cas1))
1020        ql_prof_cas(k) = ql_cas(k,it_cas2)                               &
1021     &          -frac*(ql_cas(k,it_cas2)-ql_cas(k,it_cas1))
1022        qi_prof_cas(k) = qi_cas(k,it_cas2)                               &
1023     &          -frac*(qi_cas(k,it_cas2)-qi_cas(k,it_cas1))
1024        u_prof_cas(k) = u_cas(k,it_cas2)                                 &
1025     &          -frac*(u_cas(k,it_cas2)-u_cas(k,it_cas1))
1026        v_prof_cas(k) = v_cas(k,it_cas2)                                 &
1027     &          -frac*(v_cas(k,it_cas2)-v_cas(k,it_cas1))
1028        ug_prof_cas(k) = ug_cas(k,it_cas2)                               &
1029     &          -frac*(ug_cas(k,it_cas2)-ug_cas(k,it_cas1))
1030        vg_prof_cas(k) = vg_cas(k,it_cas2)                               &
1031     &          -frac*(vg_cas(k,it_cas2)-vg_cas(k,it_cas1))
1032        vitw_prof_cas(k) = vitw_cas(k,it_cas2)                           &
1033     &          -frac*(vitw_cas(k,it_cas2)-vitw_cas(k,it_cas1))
1034        omega_prof_cas(k) = omega_cas(k,it_cas2)                         &
1035     &          -frac*(omega_cas(k,it_cas2)-omega_cas(k,it_cas1))
1036        du_prof_cas(k) = du_cas(k,it_cas2)                               &
1037     &          -frac*(du_cas(k,it_cas2)-du_cas(k,it_cas1))
1038        hu_prof_cas(k) = hu_cas(k,it_cas2)                               &
1039     &          -frac*(hu_cas(k,it_cas2)-hu_cas(k,it_cas1))
1040        vu_prof_cas(k) = vu_cas(k,it_cas2)                               &
1041     &          -frac*(vu_cas(k,it_cas2)-vu_cas(k,it_cas1))
1042        dv_prof_cas(k) = dv_cas(k,it_cas2)                               &
1043     &          -frac*(dv_cas(k,it_cas2)-dv_cas(k,it_cas1))
1044        hv_prof_cas(k) = hv_cas(k,it_cas2)                               &
1045     &          -frac*(hv_cas(k,it_cas2)-hv_cas(k,it_cas1))
1046        vv_prof_cas(k) = vv_cas(k,it_cas2)                               &
1047     &          -frac*(vv_cas(k,it_cas2)-vv_cas(k,it_cas1))
1048        dt_prof_cas(k) = dt_cas(k,it_cas2)                               &
1049     &          -frac*(dt_cas(k,it_cas2)-dt_cas(k,it_cas1))
1050        ht_prof_cas(k) = ht_cas(k,it_cas2)                               &
1051     &          -frac*(ht_cas(k,it_cas2)-ht_cas(k,it_cas1))
1052        vt_prof_cas(k) = vt_cas(k,it_cas2)                               &
1053     &          -frac*(vt_cas(k,it_cas2)-vt_cas(k,it_cas1))
1054        dth_prof_cas(k) = dth_cas(k,it_cas2)                             &
1055     &          -frac*(dth_cas(k,it_cas2)-dth_cas(k,it_cas1))
1056        hth_prof_cas(k) = hth_cas(k,it_cas2)                             &
1057     &          -frac*(hth_cas(k,it_cas2)-hth_cas(k,it_cas1))
1058        vth_prof_cas(k) = vth_cas(k,it_cas2)                             &
1059     &          -frac*(vth_cas(k,it_cas2)-vth_cas(k,it_cas1))
1060        dtrad_prof_cas(k) = dtrad_cas(k,it_cas2)                         &
1061     &          -frac*(dtrad_cas(k,it_cas2)-dtrad_cas(k,it_cas1))
1062        dq_prof_cas(k) = dq_cas(k,it_cas2)                               &
1063     &          -frac*(dq_cas(k,it_cas2)-dq_cas(k,it_cas1))
1064        hq_prof_cas(k) = hq_cas(k,it_cas2)                               &
1065     &          -frac*(hq_cas(k,it_cas2)-hq_cas(k,it_cas1))
1066        vq_prof_cas(k) = vq_cas(k,it_cas2)                               &
1067     &          -frac*(vq_cas(k,it_cas2)-vq_cas(k,it_cas1))
1068       uw_prof_cas(k) = uw_cas(k,it_cas2)                                &
1069     &          -frac*(uw_cas(k,it_cas2)-uw_cas(k,it_cas1))
1070       vw_prof_cas(k) = vw_cas(k,it_cas2)                                &
1071     &          -frac*(vw_cas(k,it_cas2)-vw_cas(k,it_cas1))
1072       q1_prof_cas(k) = q1_cas(k,it_cas2)                                &
1073     &          -frac*(q1_cas(k,it_cas2)-q1_cas(k,it_cas1))
1074       q2_prof_cas(k) = q2_cas(k,it_cas2)                                &
1075     &          -frac*(q2_cas(k,it_cas2)-q2_cas(k,it_cas1))
1076        enddo
1077
1078        return
1079        END SUBROUTINE interp_case_time_std
1080
1081!**********************************************************************************************
1082
1083END MODULE mod_1D_cases_read_std
Note: See TracBrowser for help on using the repository browser.