source: LMDZ6/trunk/libf/phylmd/dyn1d/mod_1D_cases_read2.F90 @ 4273

Last change on this file since 4273 was 4273, checked in by lguez, 2 years ago

Bug fix: split too long lines

The Fortran 2003 standard says a line may contain no more than 132
characters.
File size: 60.5 KB
Line 
1!
2! $Id: mod_1D_cases_read.F90 2373 2015-10-13 17:28:01Z jyg $
3!
4MODULE mod_1D_cases_read2
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  SUBROUTINE read_1D_cas
81    implicit none
82
83#include "netcdf.inc"
84
85    INTEGER nid,rid,ierr
86    INTEGER ii,jj
87
88    fich_cas='setup/cas.nc'
89    print*,'fich_cas ',fich_cas
90    ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid)
91    print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid
92    if (ierr.NE.NF_NOERR) then
93       write(*,*) 'ERROR: GROS Pb opening forcings nc file '
94       write(*,*) NF_STRERROR(ierr)
95       stop ""
96    endif
97    !.......................................................................
98    ierr=NF_INQ_DIMID(nid,'lat',rid)
99    IF (ierr.NE.NF_NOERR) THEN
100       print*, 'Oh probleme lecture dimension lat'
101    ENDIF
102    ierr=NF_INQ_DIMLEN(nid,rid,ii)
103    print*,'OK1 nid,rid,lat',nid,rid,ii
104    !.......................................................................
105    ierr=NF_INQ_DIMID(nid,'lon',rid)
106    IF (ierr.NE.NF_NOERR) THEN
107       print*, 'Oh probleme lecture dimension lon'
108    ENDIF
109    ierr=NF_INQ_DIMLEN(nid,rid,jj)
110    print*,'OK2 nid,rid,lat',nid,rid,jj
111    !.......................................................................
112    ierr=NF_INQ_DIMID(nid,'lev',rid)
113    IF (ierr.NE.NF_NOERR) THEN
114       print*, 'Oh probleme lecture dimension zz'
115    ENDIF
116    ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas)
117    print*,'OK3 nid,rid,nlev_cas',nid,rid,nlev_cas
118    !.......................................................................
119    ierr=NF_INQ_DIMID(nid,'time',rid)
120    print*,'nid,rid',nid,rid
121    nt_cas=0
122    IF (ierr.NE.NF_NOERR) THEN
123       stop 'probleme lecture dimension sens'
124    ENDIF
125    ierr=NF_INQ_DIMLEN(nid,rid,nt_cas)
126    print*,'OK4 nid,rid,nt_cas',nid,rid,nt_cas
127
128!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
129    !profils moyens:
130    allocate(plev_cas(nlev_cas,nt_cas))       
131    allocate(z_cas(nlev_cas,nt_cas))
132    allocate(t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas),rh_cas(nlev_cas,nt_cas))
133    allocate(th_cas(nlev_cas,nt_cas),rv_cas(nlev_cas,nt_cas))
134    allocate(u_cas(nlev_cas,nt_cas))
135    allocate(v_cas(nlev_cas,nt_cas))
136
137    !forcing
138    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))
139    allocate(hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas),dq_cas(nlev_cas,nt_cas))
140    allocate(hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas),dth_cas(nlev_cas,nt_cas))
141    allocate(hr_cas(nlev_cas,nt_cas),vr_cas(nlev_cas,nt_cas),dr_cas(nlev_cas,nt_cas))
142    allocate(hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas),du_cas(nlev_cas,nt_cas))
143    allocate(hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas),dv_cas(nlev_cas,nt_cas))
144    allocate(vitw_cas(nlev_cas,nt_cas))
145    allocate(ug_cas(nlev_cas,nt_cas))
146    allocate(vg_cas(nlev_cas,nt_cas))
147    allocate(lat_cas(nt_cas),sens_cas(nt_cas),ts_cas(nt_cas),ps_cas(nt_cas),ustar_cas(nt_cas))
148    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))
149
150
151    !champs interpoles
152    allocate(plev_prof_cas(nlev_cas))
153    allocate(t_prof_cas(nlev_cas))
154    allocate(q_prof_cas(nlev_cas))
155    allocate(u_prof_cas(nlev_cas))
156    allocate(v_prof_cas(nlev_cas))
157
158    allocate(vitw_prof_cas(nlev_cas))
159    allocate(ug_prof_cas(nlev_cas))
160    allocate(vg_prof_cas(nlev_cas))
161    allocate(ht_prof_cas(nlev_cas))
162    allocate(hq_prof_cas(nlev_cas))
163    allocate(hu_prof_cas(nlev_cas))
164    allocate(hv_prof_cas(nlev_cas))
165    allocate(vt_prof_cas(nlev_cas))
166    allocate(vq_prof_cas(nlev_cas))
167    allocate(vu_prof_cas(nlev_cas))
168    allocate(vv_prof_cas(nlev_cas))
169    allocate(dt_prof_cas(nlev_cas))
170    allocate(dtrad_prof_cas(nlev_cas))
171    allocate(dq_prof_cas(nlev_cas))
172    allocate(du_prof_cas(nlev_cas))
173    allocate(dv_prof_cas(nlev_cas))
174    allocate(uw_prof_cas(nlev_cas))
175    allocate(vw_prof_cas(nlev_cas))
176    allocate(q1_prof_cas(nlev_cas))
177    allocate(q2_prof_cas(nlev_cas))
178
179    print*,'Allocations OK'
180    call read_cas2(nid,nlev_cas,nt_cas                                       &
181         ,z_cas,plev_cas,t_cas,q_cas,rh_cas,th_cas,rv_cas,u_cas,v_cas         &
182         ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas,dv_cas,hv_cas,vv_cas    &
183         ,dt_cas,dtrad_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas                 &
184         ,dth_cas,hth_cas,vth_cas,dr_cas,hr_cas,vr_cas,sens_cas,lat_cas,ts_cas&
185         ,ustar_cas,uw_cas,vw_cas,q1_cas,q2_cas)
186    print*,'Read cas OK'
187
188
189  END SUBROUTINE read_1D_cas
190  !**********************************************************************************************
191  SUBROUTINE read2_1D_cas
192    implicit none
193
194#include "netcdf.inc"
195
196    INTEGER nid,rid,ierr
197    INTEGER ii,jj
198
199    fich_cas='setup/cas.nc'
200    print*,'fich_cas ',fich_cas
201    ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid)
202    print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid
203    if (ierr.NE.NF_NOERR) then
204       write(*,*) 'ERROR: GROS Pb opening forcings nc file '
205       write(*,*) NF_STRERROR(ierr)
206       stop ""
207    endif
208    !.......................................................................
209    ierr=NF_INQ_DIMID(nid,'lat',rid)
210    IF (ierr.NE.NF_NOERR) THEN
211       print*, 'Oh probleme lecture dimension lat'
212    ENDIF
213    ierr=NF_INQ_DIMLEN(nid,rid,ii)
214    print*,'OK1 read2: nid,rid,lat',nid,rid,ii
215    !.......................................................................
216    ierr=NF_INQ_DIMID(nid,'lon',rid)
217    IF (ierr.NE.NF_NOERR) THEN
218       print*, 'Oh probleme lecture dimension lon'
219    ENDIF
220    ierr=NF_INQ_DIMLEN(nid,rid,jj)
221    print*,'OK2 read2: nid,rid,lat',nid,rid,jj
222    !.......................................................................
223    ierr=NF_INQ_DIMID(nid,'nlev',rid)
224    IF (ierr.NE.NF_NOERR) THEN
225       print*, 'Oh probleme lecture dimension nlev'
226    ENDIF
227    ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas)
228    print*,'OK3 read2: nid,rid,nlev_cas',nid,rid,nlev_cas
229    !.......................................................................
230    ierr=NF_INQ_DIMID(nid,'time',rid)
231    nt_cas=0
232    IF (ierr.NE.NF_NOERR) THEN
233       stop 'Oh probleme lecture dimension time'
234    ENDIF
235    ierr=NF_INQ_DIMLEN(nid,rid,nt_cas)
236    print*,'OK4 read2: nid,rid,nt_cas',nid,rid,nt_cas
237
238!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
239    !profils moyens:
240    allocate(plev_cas(nlev_cas,nt_cas),plevh_cas(nlev_cas+1))       
241    allocate(z_cas(nlev_cas,nt_cas),zh_cas(nlev_cas+1))
242    allocate(ap_cas(nlev_cas+1),bp_cas(nt_cas+1))
243    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), &
244         qi_cas(nlev_cas,nt_cas),rh_cas(nlev_cas,nt_cas))
245    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))
246    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))
247
248    !forcing
249    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))
250    allocate(hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas),dq_cas(nlev_cas,nt_cas))
251    allocate(hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas),dth_cas(nlev_cas,nt_cas))
252    allocate(hr_cas(nlev_cas,nt_cas),vr_cas(nlev_cas,nt_cas),dr_cas(nlev_cas,nt_cas))
253    allocate(hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas),du_cas(nlev_cas,nt_cas))
254    allocate(hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas),dv_cas(nlev_cas,nt_cas))
255    allocate(ug_cas(nlev_cas,nt_cas))
256    allocate(vg_cas(nlev_cas,nt_cas))
257    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))
258    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))
259
260
261
262    !champs interpoles
263    allocate(plev_prof_cas(nlev_cas))
264    allocate(t_prof_cas(nlev_cas))
265    allocate(theta_prof_cas(nlev_cas))
266    allocate(thl_prof_cas(nlev_cas))
267    allocate(thv_prof_cas(nlev_cas))
268    allocate(q_prof_cas(nlev_cas))
269    allocate(qv_prof_cas(nlev_cas))
270    allocate(ql_prof_cas(nlev_cas))
271    allocate(qi_prof_cas(nlev_cas))
272    allocate(rh_prof_cas(nlev_cas))
273    allocate(rv_prof_cas(nlev_cas))
274    allocate(u_prof_cas(nlev_cas))
275    allocate(v_prof_cas(nlev_cas))
276    allocate(vitw_prof_cas(nlev_cas))
277    allocate(omega_prof_cas(nlev_cas))
278    allocate(ug_prof_cas(nlev_cas))
279    allocate(vg_prof_cas(nlev_cas))
280    allocate(ht_prof_cas(nlev_cas))
281    allocate(hth_prof_cas(nlev_cas))
282    allocate(hq_prof_cas(nlev_cas))
283    allocate(hu_prof_cas(nlev_cas))
284    allocate(hv_prof_cas(nlev_cas))
285    allocate(vt_prof_cas(nlev_cas))
286    allocate(vth_prof_cas(nlev_cas))
287    allocate(vq_prof_cas(nlev_cas))
288    allocate(vu_prof_cas(nlev_cas))
289    allocate(vv_prof_cas(nlev_cas))
290    allocate(dt_prof_cas(nlev_cas))
291    allocate(dth_prof_cas(nlev_cas))
292    allocate(dtrad_prof_cas(nlev_cas))
293    allocate(dq_prof_cas(nlev_cas))
294    allocate(du_prof_cas(nlev_cas))
295    allocate(dv_prof_cas(nlev_cas))
296    allocate(uw_prof_cas(nlev_cas))
297    allocate(vw_prof_cas(nlev_cas))
298    allocate(q1_prof_cas(nlev_cas))
299    allocate(q2_prof_cas(nlev_cas))
300
301    print*,'Allocations OK'
302    call read2_cas (nid,nlev_cas,nt_cas,                                                                     &
303         ap_cas,bp_cas,z_cas,plev_cas,zh_cas,plevh_cas,t_cas,th_cas,thv_cas,thl_cas,qv_cas,                    &
304         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,        &
305         dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas,dth_cas,hth_cas,vth_cas,               &
306         dr_cas,hr_cas,vr_cas,dtrad_cas,sens_cas,lat_cas,ts_cas,ps_cas,ustar_cas,tke_cas,                      &
307         uw_cas,vw_cas,q1_cas,q2_cas,orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,heat_rough, &
308         o3_cas,rugos_cas,clay_cas,sand_cas)
309    print*,'Read2 cas OK'
310    do ii=1,nlev_cas
311       print*,'apres read2_cas, plev_cas=',ii,plev_cas(ii,1)
312    enddo
313
314
315  END SUBROUTINE read2_1D_cas
316
317  !**********************************************************************************************
318  SUBROUTINE old_read_SCM_cas
319    implicit none
320
321#include "netcdf.inc"
322#include "date_cas.h"
323
324    INTEGER nid,rid,ierr
325    INTEGER ii,jj,timeid
326    REAL, ALLOCATABLE :: time_val(:)
327
328    fich_cas='cas.nc'
329    print*,'fich_cas ',fich_cas
330    ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid)
331    print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid
332    if (ierr.NE.NF_NOERR) then
333       write(*,*) 'ERROR: GROS Pb opening forcings nc file '
334       write(*,*) NF_STRERROR(ierr)
335       stop ""
336    endif
337    !.......................................................................
338    ierr=NF_INQ_DIMID(nid,'lat',rid)
339    IF (ierr.NE.NF_NOERR) THEN
340       print*, 'Oh probleme lecture dimension lat'
341    ENDIF
342    ierr=NF_INQ_DIMLEN(nid,rid,ii)
343    print*,'OK1 read2: nid,rid,lat',nid,rid,ii
344    !.......................................................................
345    ierr=NF_INQ_DIMID(nid,'lon',rid)
346    IF (ierr.NE.NF_NOERR) THEN
347       print*, 'Oh probleme lecture dimension lon'
348    ENDIF
349    ierr=NF_INQ_DIMLEN(nid,rid,jj)
350    print*,'OK2 read2: nid,rid,lat',nid,rid,jj
351    !.......................................................................
352    ierr=NF_INQ_DIMID(nid,'lev',rid)
353    IF (ierr.NE.NF_NOERR) THEN
354       print*, 'Oh probleme lecture dimension nlev'
355    ENDIF
356    ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas)
357    print*,'OK3 read2: nid,rid,nlev_cas',nid,rid,nlev_cas
358    IF ( .NOT. ( nlev_cas > 10 .AND. nlev_cas < 1000 )) THEN
359       print*,'Valeur de nlev_cas peu probable'
360       STOP
361    ENDIF
362    !.......................................................................
363    ierr=NF_INQ_DIMID(nid,'time',rid)
364    nt_cas=0
365    IF (ierr.NE.NF_NOERR) THEN
366       stop 'Oh probleme lecture dimension time'
367    ENDIF
368    ierr=NF_INQ_DIMLEN(nid,rid,nt_cas)
369    print*,'OK4 read2: nid,rid,nt_cas',nid,rid,nt_cas
370    ! Lecture de l'axe des temps
371    print*,'LECTURE DU TEMPS'
372    ierr=NF_INQ_VARID(nid,'time',timeid)
373    if(ierr/=NF_NOERR) then
374       print *,'Variable time manquante dans cas.nc:'
375       ierr=NF_NOERR
376    else
377       allocate(time_val(nt_cas))
378#ifdef NC_DOUBLE
379       ierr = NF_GET_VAR_DOUBLE(nid,timeid,time_val)
380#else
381       ierr = NF_GET_VAR_REAL(nid,timeid,time_val)
382#endif
383       if(ierr/=NF_NOERR) then
384          print *,'Pb a la lecture de time cas.nc: '
385       endif
386    endif
387    IF (nt_cas>1) THEN
388       pdt_cas=time_val(2)-time_val(1)
389    ELSE
390       pdt_cas=0.
391    ENDIF
392
393
394!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
395    !profils moyens:
396    allocate(plev_cas(nlev_cas,nt_cas),plevh_cas(nlev_cas+1))       
397    allocate(z_cas(nlev_cas,nt_cas),zh_cas(nlev_cas+1))
398    allocate(ap_cas(nlev_cas+1),bp_cas(nt_cas+1))
399    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), &
400         qi_cas(nlev_cas,nt_cas),rh_cas(nlev_cas,nt_cas))
401    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))
402    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))
403
404    !forcing
405    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))
406    allocate(hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas),dq_cas(nlev_cas,nt_cas))
407    allocate(hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas),dth_cas(nlev_cas,nt_cas))
408    allocate(hr_cas(nlev_cas,nt_cas),vr_cas(nlev_cas,nt_cas),dr_cas(nlev_cas,nt_cas))
409    allocate(hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas),du_cas(nlev_cas,nt_cas))
410    allocate(hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas),dv_cas(nlev_cas,nt_cas))
411    allocate(ug_cas(nlev_cas,nt_cas))
412    allocate(vg_cas(nlev_cas,nt_cas))
413    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))
414    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))
415
416
417
418    !champs interpoles
419    allocate(plev_prof_cas(nlev_cas))
420    allocate(t_prof_cas(nlev_cas))
421    allocate(theta_prof_cas(nlev_cas))
422    allocate(thl_prof_cas(nlev_cas))
423    allocate(thv_prof_cas(nlev_cas))
424    allocate(q_prof_cas(nlev_cas))
425    allocate(qv_prof_cas(nlev_cas))
426    allocate(ql_prof_cas(nlev_cas))
427    allocate(qi_prof_cas(nlev_cas))
428    allocate(rh_prof_cas(nlev_cas))
429    allocate(rv_prof_cas(nlev_cas))
430    allocate(u_prof_cas(nlev_cas))
431    allocate(v_prof_cas(nlev_cas))
432    allocate(vitw_prof_cas(nlev_cas))
433    allocate(omega_prof_cas(nlev_cas))
434    allocate(ug_prof_cas(nlev_cas))
435    allocate(vg_prof_cas(nlev_cas))
436    allocate(ht_prof_cas(nlev_cas))
437    allocate(hth_prof_cas(nlev_cas))
438    allocate(hq_prof_cas(nlev_cas))
439    allocate(hu_prof_cas(nlev_cas))
440    allocate(hv_prof_cas(nlev_cas))
441    allocate(vt_prof_cas(nlev_cas))
442    allocate(vth_prof_cas(nlev_cas))
443    allocate(vq_prof_cas(nlev_cas))
444    allocate(vu_prof_cas(nlev_cas))
445    allocate(vv_prof_cas(nlev_cas))
446    allocate(dt_prof_cas(nlev_cas))
447    allocate(dth_prof_cas(nlev_cas))
448    allocate(dtrad_prof_cas(nlev_cas))
449    allocate(dq_prof_cas(nlev_cas))
450    allocate(du_prof_cas(nlev_cas))
451    allocate(dv_prof_cas(nlev_cas))
452    allocate(uw_prof_cas(nlev_cas))
453    allocate(vw_prof_cas(nlev_cas))
454    allocate(q1_prof_cas(nlev_cas))
455    allocate(q2_prof_cas(nlev_cas))
456
457    print*,'Allocations OK'
458    call old_read_SCM (nid,nlev_cas,nt_cas,                                                                     &
459         ap_cas,bp_cas,z_cas,plev_cas,zh_cas,plevh_cas,t_cas,th_cas,thv_cas,thl_cas,qv_cas,                    &
460         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,        &
461         dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas,dth_cas,hth_cas,vth_cas,               &
462         dr_cas,hr_cas,vr_cas,dtrad_cas,sens_cas,lat_cas,ts_cas,ps_cas,ustar_cas,tke_cas,                      &
463         uw_cas,vw_cas,q1_cas,q2_cas,orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,heat_rough, &
464         o3_cas,rugos_cas,clay_cas,sand_cas)
465    print*,'Read2 cas OK'
466    do ii=1,nlev_cas
467       print*,'apres read2_cas, plev_cas=',ii,plev_cas(ii,1)
468    enddo
469
470
471  END SUBROUTINE old_read_SCM_cas
472
473
474!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
475  SUBROUTINE deallocate2_1D_cases
476    !profils environnementaux:
477    deallocate(plev_cas,plevh_cas)
478
479    deallocate(z_cas,zh_cas)
480    deallocate(ap_cas,bp_cas)
481    deallocate(t_cas,q_cas,qv_cas,ql_cas,qi_cas,rh_cas)
482    deallocate(th_cas,thl_cas,thv_cas,rv_cas)
483    deallocate(u_cas,v_cas,vitw_cas,omega_cas)
484
485    !forcing
486    deallocate(ht_cas,vt_cas,dt_cas,dtrad_cas)
487    deallocate(hq_cas,vq_cas,dq_cas)
488    deallocate(hth_cas,vth_cas,dth_cas)
489    deallocate(hr_cas,vr_cas,dr_cas)
490    deallocate(hu_cas,vu_cas,du_cas)
491    deallocate(hv_cas,vv_cas,dv_cas)
492    deallocate(ug_cas)
493    deallocate(vg_cas)
494    deallocate(lat_cas,sens_cas,ts_cas,ps_cas,ustar_cas,tke_cas,uw_cas,vw_cas,q1_cas,q2_cas)
495
496    !champs interpoles
497    deallocate(plev_prof_cas)
498    deallocate(t_prof_cas)
499    deallocate(theta_prof_cas)
500    deallocate(thl_prof_cas)
501    deallocate(thv_prof_cas)
502    deallocate(q_prof_cas)
503    deallocate(qv_prof_cas)
504    deallocate(ql_prof_cas)
505    deallocate(qi_prof_cas)
506    deallocate(rh_prof_cas)
507    deallocate(rv_prof_cas)
508    deallocate(u_prof_cas)
509    deallocate(v_prof_cas)
510    deallocate(vitw_prof_cas)
511    deallocate(omega_prof_cas)
512    deallocate(ug_prof_cas)
513    deallocate(vg_prof_cas)
514    deallocate(ht_prof_cas)
515    deallocate(hq_prof_cas)
516    deallocate(hu_prof_cas)
517    deallocate(hv_prof_cas)
518    deallocate(vt_prof_cas)
519    deallocate(vq_prof_cas)
520    deallocate(vu_prof_cas)
521    deallocate(vv_prof_cas)
522    deallocate(dt_prof_cas)
523    deallocate(dtrad_prof_cas)
524    deallocate(dq_prof_cas)
525    deallocate(du_prof_cas)
526    deallocate(dv_prof_cas)
527    deallocate(t_prof_cas)
528    deallocate(u_prof_cas)
529    deallocate(v_prof_cas)
530    deallocate(uw_prof_cas)
531    deallocate(vw_prof_cas)
532    deallocate(q1_prof_cas)
533    deallocate(q2_prof_cas)
534
535  END SUBROUTINE deallocate2_1D_cases
536
537
538END MODULE mod_1D_cases_read2
539!=====================================================================
540subroutine read_cas2(nid,nlevel,ntime                          &
541     ,zz,pp,temp,qv,rh,theta,rv,u,v,ug,vg,w,                   &
542     du,hu,vu,dv,hv,vv,dt,dtrad,ht,vt,dq,hq,vq,                     &
543     dth,hth,vth,dr,hr,vr,sens,flat,ts,ustar,uw,vw,q1,q2)
544
545  !program reading forcing of the case study
546  implicit none
547#include "netcdf.inc"
548
549  integer ntime,nlevel
550
551  real zz(nlevel,ntime)
552  real pp(nlevel,ntime)
553  real temp(nlevel,ntime),qv(nlevel,ntime),rh(nlevel,ntime)
554  real theta(nlevel,ntime),rv(nlevel,ntime)
555  real u(nlevel,ntime)
556  real v(nlevel,ntime)
557  real ug(nlevel,ntime)
558  real vg(nlevel,ntime)
559  real w(nlevel,ntime)
560  real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
561  real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
562  real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
563  real dtrad(nlevel,ntime)
564  real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
565  real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime)
566  real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
567  real flat(ntime),sens(ntime),ts(ntime),ustar(ntime)
568  real uw(nlevel,ntime),vw(nlevel,ntime),q1(nlevel,ntime),q2(nlevel,ntime),resul(nlevel,ntime),resul1(ntime)
569
570
571  integer nid, ierr, ierr1,ierr2,rid,i
572  integer nbvar3d
573  parameter(nbvar3d=39)
574  integer var3didin(nbvar3d)
575  character*5 name_var(1:nbvar3d)
576  data name_var/'zz','pp','temp','qv','rh','theta','rv','u','v','ug','vg','w','advu','hu','vu',&
577       'advv','hv','vv','advT','hT','vT','advq','hq','vq','advth','hth','vth','advr','hr','vr',&
578       'radT','uw','vw','q1','q2','sens','flat','ts','ustar'/
579
580  do i=1,nbvar3d
581     print *,'Dans read_cas2, on va lire ',nid,i,name_var(i)
582  enddo
583  do i=1,nbvar3d
584     ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
585     print *,'ierr=',i,ierr,name_var(i),var3didin(i)
586     if(ierr/=NF_NOERR) then
587        print *,'Variable manquante dans cas.nc:',name_var(i)
588     endif
589  enddo
590  do i=1,nbvar3d
591     print *,'Dans read_cas2, on va lire ',var3didin(i),name_var(i)
592     if(i.LE.35) then
593#ifdef NC_DOUBLE
594        ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
595#else
596        ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
597#endif
598        print *,'Dans read_cas2, on a lu ',ierr,var3didin(i),name_var(i)
599        if(ierr/=NF_NOERR) then
600           print *,'Pb a la lecture de cas.nc: ',name_var(i)
601           stop "getvarup"
602        endif
603     else
604#ifdef NC_DOUBLE
605        ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul1)
606#else
607        ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul1)
608#endif
609        print *,'Dans read_cas2, on a lu ',ierr,var3didin(i),name_var(i)
610        if(ierr/=NF_NOERR) then
611           print *,'Pb a la lecture de cas.nc: ',name_var(i)
612           stop "getvarup"
613        endif
614     endif
615     select case(i)
616     case(1) ; zz=resul
617     case(2) ; pp=resul
618     case(3) ; temp=resul
619     case(4) ; qv=resul
620     case(5) ; rh=resul
621     case(6) ; theta=resul
622     case(7) ; rv=resul
623     case(8) ; u=resul
624     case(9) ; v=resul
625     case(10) ; ug=resul
626     case(11) ; vg=resul
627     case(12) ; w=resul
628     case(13) ; du=resul
629     case(14) ; hu=resul
630     case(15) ; vu=resul
631     case(16) ; dv=resul
632     case(17) ; hv=resul
633     case(18) ; vv=resul
634     case(19) ; dt=resul
635     case(20) ; ht=resul
636     case(21) ; vt=resul
637     case(22) ; dq=resul
638     case(23) ; hq=resul
639     case(24) ; vq=resul
640     case(25) ; dth=resul
641     case(26) ; hth=resul
642     case(27) ; vth=resul
643     case(28) ; dr=resul
644     case(29) ; hr=resul
645     case(30) ; vr=resul
646     case(31) ; dtrad=resul
647     case(32) ; uw=resul
648     case(33) ; vw=resul
649     case(34) ; q1=resul
650     case(35) ; q2=resul
651     case(36) ; sens=resul1
652     case(37) ; flat=resul1
653     case(38) ; ts=resul1
654     case(39) ; ustar=resul1
655     end select
656  enddo
657
658  return
659end subroutine read_cas2
660!======================================================================
661subroutine read2_cas(nid,nlevel,ntime,                                       &
662     ap,bp,zz,pp,zzh,pph,temp,theta,thv,thl,qv,ql,qi,rh,rv,u,v,vitw,omega,ug,vg,&
663     du,hu,vu,dv,hv,vv,dt,ht,vt,dq,hq,vq,                                    &
664     dth,hth,vth,dr,hr,vr,dtrad,sens,flat,ts,ps,ustar,tke,uw,vw,q1,q2,       &
665     orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,          &
666     heat_rough,o3_cas,rugos_cas,clay_cas,sand_cas)
667
668  !program reading forcing of the case study
669  implicit none
670#include "netcdf.inc"
671
672  integer ntime,nlevel
673
674  real ap(nlevel+1),bp(nlevel+1)
675  real zz(nlevel,ntime),zzh(nlevel+1)
676  real pp(nlevel,ntime),pph(nlevel+1)
677  real temp(nlevel,ntime),qv(nlevel,ntime),ql(nlevel,ntime),qi(nlevel,ntime),rh(nlevel,ntime)
678  real theta(nlevel,ntime),thv(nlevel,ntime),thl(nlevel,ntime),rv(nlevel,ntime)
679  real u(nlevel,ntime),v(nlevel,ntime)
680  real ug(nlevel,ntime),vg(nlevel,ntime)
681  real vitw(nlevel,ntime),omega(nlevel,ntime)
682  real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
683  real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
684  real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
685  real dtrad(nlevel,ntime)
686  real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
687  real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime),hthl(nlevel,ntime)
688  real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
689  real flat(ntime),sens(ntime),ustar(ntime)
690  real uw(nlevel,ntime),vw(nlevel,ntime),q1(nlevel,ntime),q2(nlevel,ntime)
691  real ts(ntime),ps(ntime),tke(ntime)
692  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
693  real apbp(nlevel+1),resul(nlevel,ntime),resul1(nlevel),resul2(ntime),resul3
694
695
696  integer nid, ierr,ierr1,ierr2,rid,i
697  integer nbvar3d
698  parameter(nbvar3d=62)
699  integer var3didin(nbvar3d),missing_var(nbvar3d)
700  character*12 name_var(1:nbvar3d)
701  data name_var/'coor_par_a','coor_par_b','height_h','pressure_h',&
702       'w','omega','ug','vg','uadv','uadvh','uadvv','vadv','vadvh','vadvv','tadv','tadvh','tadvv',&
703       'qadv','qadvh','qadvv','thadv','thadvh','thadvv','thladvh','radv','radvh','radvv','radcool','q1','q2','ustress','vstress', &
704       'rh',&
705       'height_f','pressure_f','temp','theta','thv','thl','qv','ql','qi','rv','u','v',&
706       'sfc_sens_flx','sfc_lat_flx','ts','ps','ustar','tke',&
707       'orog','albedo','emiss','t_skin','q_skin','mom_rough','heat_rough','o3','rugos','clay','sand'/
708  do i=1,nbvar3d
709     missing_var(i)=0.
710  enddo
711
712  !-----------------------------------------------------------------------
713  do i=1,nbvar3d
714     ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
715     if(ierr/=NF_NOERR) then
716        print *,'Variable manquante dans cas.nc:',i,name_var(i)
717        ierr=NF_NOERR
718        missing_var(i)=1
719     else
720        !-----------------------------------------------------------------------
721        if(i.LE.4) then     ! Lecture des coord pression en (nlevelp1,lat,lon)
722#ifdef NC_DOUBLE
723           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),apbp)
724#else
725           ierr = NF_GET_VAR_REAL(nid,var3didin(i),apbp)
726#endif
727           print *,'read2_cas(apbp), on a lu ',i,name_var(i)
728           if(ierr/=NF_NOERR) then
729              print *,'Pb a la lecture de cas.nc: ',name_var(i)
730              stop "getvarup"
731           endif
732           !-----------------------------------------------------------------------
733        else if(i.gt.4.and.i.LE.45) then   ! Lecture des variables en (time,nlevel,lat,lon)
734#ifdef NC_DOUBLE
735           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
736#else
737           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
738#endif
739           print *,'read2_cas(resul), on a lu ',i,name_var(i)
740           if(ierr/=NF_NOERR) then
741              print *,'Pb a la lecture de cas.nc: ',name_var(i)
742              stop "getvarup"
743           endif
744           !-----------------------------------------------------------------------
745        else if (i.gt.45.and.i.LE.51) then   ! Lecture des variables en (time,lat,lon)
746#ifdef NC_DOUBLE
747           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul2)
748#else
749           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul2)
750#endif
751           print *,'read2_cas(resul2), on a lu ',i,name_var(i)
752           if(ierr/=NF_NOERR) then
753              print *,'Pb a la lecture de cas.nc: ',name_var(i)
754              stop "getvarup"
755           endif
756           !-----------------------------------------------------------------------
757        else     ! Lecture des constantes (lat,lon)
758#ifdef NC_DOUBLE
759           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul3)
760#else
761           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul3)
762#endif
763           print *,'read2_cas(resul3), on a lu ',i,name_var(i)
764           if(ierr/=NF_NOERR) then
765              print *,'Pb a la lecture de cas.nc: ',name_var(i)
766              stop "getvarup"
767           endif
768        endif
769     endif
770     !-----------------------------------------------------------------------
771     select case(i)
772     case(1) ; ap=apbp       ! donnees indexees en nlevel+1
773     case(2) ; bp=apbp
774     case(3) ; zzh=apbp
775     case(4) ; pph=apbp
776     case(5) ; vitw=resul    ! donnees indexees en nlevel,time
777     case(6) ; omega=resul
778     case(7) ; ug=resul
779     case(8) ; vg=resul
780     case(9) ; du=resul
781     case(10) ; hu=resul
782     case(11) ; vu=resul
783     case(12) ; dv=resul
784     case(13) ; hv=resul
785     case(14) ; vv=resul
786     case(15) ; dt=resul
787     case(16) ; ht=resul
788     case(17) ; vt=resul
789     case(18) ; dq=resul
790     case(19) ; hq=resul
791     case(20) ; vq=resul
792     case(21) ; dth=resul
793     case(22) ; hth=resul
794     case(23) ; vth=resul
795     case(24) ; hthl=resul
796     case(25) ; dr=resul
797     case(26) ; hr=resul
798     case(27) ; vr=resul
799     case(28) ; dtrad=resul
800     case(29) ; q1=resul
801     case(30) ; q2=resul
802     case(31) ; uw=resul
803     case(32) ; vw=resul
804     case(33) ; rh=resul
805     case(34) ; zz=resul      ! donnees en time,nlevel pour profil initial
806     case(35) ; pp=resul
807     case(36) ; temp=resul
808     case(37) ; theta=resul
809     case(38) ; thv=resul
810     case(39) ; thl=resul
811     case(40) ; qv=resul
812     case(41) ; ql=resul
813     case(42) ; qi=resul
814     case(43) ; rv=resul
815     case(44) ; u=resul
816     case(45) ; v=resul
817     case(46) ; sens=resul2   ! donnees indexees en time
818     case(47) ; flat=resul2
819     case(48) ; ts=resul2
820     case(49) ; ps=resul2
821     case(50) ; ustar=resul2
822     case(51) ; tke=resul2
823     case(52) ; orog_cas=resul3      ! constantes
824     case(53) ; albedo_cas=resul3
825     case(54) ; emiss_cas=resul3
826     case(55) ; t_skin_cas=resul3
827     case(56) ; q_skin_cas=resul3
828     case(57) ; mom_rough=resul3
829     case(58) ; heat_rough=resul3
830     case(59) ; o3_cas=resul3       
831     case(60) ; rugos_cas=resul3
832     case(61) ; clay_cas=resul3
833     case(62) ; sand_cas=resul3
834     end select
835     resul=0.
836     resul1=0.
837     resul2=0.
838     resul3=0.
839  enddo
840  !-----------------------------------------------------------------------
841
842
843  return
844end subroutine read2_cas
845
846!======================================================================
847subroutine old_read_SCM(nid,nlevel,ntime,                                       &
848     ap,bp,zz,pp,zzh,pph,temp,theta,thv,thl,qv,ql,qi,rh,rv,u,v,vitw,omega,ug,vg,&
849     du,hu,vu,dv,hv,vv,dt,ht,vt,dq,hq,vq,                                    &
850     dth,hth,vth,dr,hr,vr,dtrad,sens,flat,ts,ps,ustar,tke,uw,vw,q1,q2,       &
851     orog_cas,albedo_cas,emiss_cas,t_skin_cas,q_skin_cas,mom_rough,          &
852     heat_rough,o3_cas,rugos_cas,clay_cas,sand_cas)
853
854  !program reading forcing of the case study
855  implicit none
856#include "netcdf.inc"
857
858  integer ntime,nlevel,k,t
859
860  real ap(nlevel+1),bp(nlevel+1)
861  real zz(nlevel,ntime),zzh(nlevel+1)
862  real pp(nlevel,ntime),pph(nlevel+1)
863  !profils initiaux
864  real temp0(nlevel),qv0(nlevel),ql0(nlevel),qi0(nlevel),u0(nlevel),v0(nlevel),tke0(nlevel)
865  real pp0(nlevel)   
866  real temp(nlevel,ntime),qv(nlevel,ntime),ql(nlevel,ntime),qi(nlevel,ntime),rh(nlevel,ntime)
867  real theta(nlevel,ntime),thv(nlevel,ntime),thl(nlevel,ntime),rv(nlevel,ntime)
868  real u(nlevel,ntime),v(nlevel,ntime),tke(nlevel,ntime)
869  real ug(nlevel,ntime),vg(nlevel,ntime)
870  real vitw(nlevel,ntime),omega(nlevel,ntime)
871  real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
872  real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
873  real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
874  real dtrad(nlevel,ntime)
875  real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
876  real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime),hthl(nlevel,ntime)
877  real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
878  real flat(ntime),sens(ntime),ustar(ntime)
879  real uw(nlevel,ntime),vw(nlevel,ntime),q1(nlevel,ntime),q2(nlevel,ntime)
880  real ts(ntime),ps(ntime)
881  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
882  real apbp(nlevel+1),resul(nlevel,ntime),resul1(nlevel),resul2(ntime),resul3
883
884
885  integer nid, ierr,ierr1,ierr2,rid,i
886  integer nbvar3d
887  parameter(nbvar3d=70)
888  integer var3didin(nbvar3d),missing_var(nbvar3d)
889  character*13 name_var(1:nbvar3d)
890  data name_var/'coor_par_a','coor_par_b','height_h','pressure_h',&
891       'temp','qv','ql','qi','u','v','tke','pressure',&
892       'w','omega','ug','vg','uadv','uadvh','uadvv','vadv','vadvh','vadvv','tadv','tadvh','tadvv',&
893       'qvadv','qvadvh','qvadvv','thadv','thadvh','thadvv','thladvh','radv','radvh','radvv','radcool','q1','q2','ustress', &
894       'vstress','rh',&
895       'height_f','pressure_forc','tempt','theta','thv','thl','qvt','qlt','qit','rv','ut','vt','tket',&
896       'sfc_sens_flx','sfc_lat_flx','ts','ps','ustar',&
897       'orog','albedo','emiss','t_skin','q_skin','mom_rough','heat_rough','o3','rugos','clay','sand'/
898  do i=1,nbvar3d
899     missing_var(i)=0.
900  enddo
901
902  !-----------------------------------------------------------------------
903
904  print*,'ON EST LA'
905  do i=1,nbvar3d
906     ierr=NF_INQ_VARID(nid,name_var(i),var3didin(i))
907     if(ierr/=NF_NOERR) then
908        print *,'Variable manquante dans cas.nc:',i,name_var(i)
909        ierr=NF_NOERR
910        missing_var(i)=1
911     else
912        !-----------------------------------------------------------------------
913        if(i.LE.4) then     ! Lecture des coord pression en (nlevelp1,lat,lon)
914#ifdef NC_DOUBLE
915           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),apbp)
916#else
917           ierr = NF_GET_VAR_REAL(nid,var3didin(i),apbp)
918#endif
919           print *,'read2_cas(apbp), on a lu ',i,name_var(i)
920           if(ierr/=NF_NOERR) then
921              print *,'Pb a la lecture de cas.nc: ',name_var(i)
922              stop "getvarup"
923           endif
924           !-----------------------------------------------------------------------
925        else if(i.gt.4.and.i.LE.12) then   ! Lecture des variables en (time,nlevel,lat,lon)
926#ifdef NC_DOUBLE
927           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul1)
928#else
929           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul1)
930#endif
931           print *,'read2_cas(resul1), on a lu ',i,name_var(i)
932           if(ierr/=NF_NOERR) then
933              print *,'Pb a la lecture de cas.nc: ',name_var(i)
934              stop "getvarup"
935           endif
936           print*,'Lecture de la variable #i ',i,name_var(i),minval(resul1),maxval(resul1)
937           !-----------------------------------------------------------------------
938        else if(i.gt.12.and.i.LE.54) then   ! Lecture des variables en (time,nlevel,lat,lon)
939#ifdef NC_DOUBLE
940           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul)
941#else
942           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul)
943#endif
944           print *,'read2_cas(resul), on a lu ',i,name_var(i)
945           if(ierr/=NF_NOERR) then
946              print *,'Pb a la lecture de cas.nc: ',name_var(i)
947              stop "getvarup"
948           endif
949           print*,'Lecture de la variable #i ',i,name_var(i),minval(resul),maxval(resul)
950           !-----------------------------------------------------------------------
951        else if (i.gt.54.and.i.LE.65) then   ! Lecture des variables en (time,lat,lon)
952#ifdef NC_DOUBLE
953           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul2)
954#else
955           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul2)
956#endif
957           print *,'read2_cas(resul2), on a lu ',i,name_var(i)
958           if(ierr/=NF_NOERR) then
959              print *,'Pb a la lecture de cas.nc: ',name_var(i)
960              stop "getvarup"
961           endif
962           print*,'Lecture de la variable #i  ',i,name_var(i),minval(resul2),maxval(resul2)
963           !-----------------------------------------------------------------------
964        else     ! Lecture des constantes (lat,lon)
965#ifdef NC_DOUBLE
966           ierr = NF_GET_VAR_DOUBLE(nid,var3didin(i),resul3)
967#else
968           ierr = NF_GET_VAR_REAL(nid,var3didin(i),resul3)
969#endif
970           print *,'read2_cas(resul3), on a lu ',i,name_var(i)
971           if(ierr/=NF_NOERR) then
972              print *,'Pb a la lecture de cas.nc: ',name_var(i)
973              stop "getvarup"
974           endif
975           print*,'Lecture de la variable #i ',i,name_var(i),resul3
976        endif
977     endif
978     !-----------------------------------------------------------------------
979     select case(i)
980        !case(1) ; ap=apbp       ! donnees indexees en nlevel+1
981        ! case(2) ; bp=apbp
982     case(3) ; zzh=apbp
983     case(4) ; pph=apbp
984     case(5) ; temp0=resul1    ! donnees initiales
985     case(6) ; qv0=resul1
986     case(7) ; ql0=resul1
987     case(8) ; qi0=resul1
988     case(9) ; u0=resul1
989     case(10) ; v0=resul1
990     case(11) ; tke0=resul1
991     case(12) ; pp0=resul1
992     case(13) ; vitw=resul    ! donnees indexees en nlevel,time
993     case(14) ; omega=resul
994     case(15) ; ug=resul
995     case(16) ; vg=resul
996     case(17) ; du=resul
997     case(18) ; hu=resul
998     case(19) ; vu=resul
999     case(20) ; dv=resul
1000     case(21) ; hv=resul
1001     case(22) ; vv=resul
1002     case(23) ; dt=resul
1003     case(24) ; ht=resul
1004     case(25) ; vt=resul
1005     case(26) ; dq=resul
1006     case(27) ; hq=resul
1007     case(28) ; vq=resul
1008     case(29) ; dth=resul
1009     case(30) ; hth=resul
1010     case(31) ; vth=resul
1011     case(32) ; hthl=resul
1012     case(33) ; dr=resul
1013     case(34) ; hr=resul
1014     case(35) ; vr=resul
1015     case(36) ; dtrad=resul
1016     case(37) ; q1=resul
1017     case(38) ; q2=resul
1018     case(39) ; uw=resul
1019     case(40) ; vw=resul
1020     case(41) ; rh=resul
1021     case(42) ; zz=resul      ! donnees en time,nlevel pour profil initial
1022     case(43) ; pp=resul
1023     case(44) ; temp=resul
1024     case(45) ; theta=resul
1025     case(46) ; thv=resul
1026     case(47) ; thl=resul
1027     case(48) ; qv=resul
1028     case(49) ; ql=resul
1029     case(50) ; qi=resul
1030     case(51) ; rv=resul
1031     case(52) ; u=resul
1032     case(53) ; v=resul
1033     case(54) ; tke=resul
1034     case(55) ; sens=resul2   ! donnees indexees en time
1035     case(56) ; flat=resul2
1036     case(57) ; ts=resul2
1037     case(58) ; ps=resul2
1038     case(59) ; ustar=resul2
1039     case(60) ; orog_cas=resul3      ! constantes
1040     case(61) ; albedo_cas=resul3
1041     case(62) ; emiss_cas=resul3
1042     case(63) ; t_skin_cas=resul3
1043     case(64) ; q_skin_cas=resul3
1044     case(65) ; mom_rough=resul3
1045     case(66) ; heat_rough=resul3
1046     case(67) ; o3_cas=resul3       
1047     case(68) ; rugos_cas=resul3
1048     case(69) ; clay_cas=resul3
1049     case(70) ; sand_cas=resul3
1050     end select
1051     resul=0.
1052     resul1=0.
1053     resul2=0.
1054     resul3=0.
1055  enddo
1056  print*,'Lecture de la variable APRES ,sens ',minval(sens),maxval(sens)
1057  print*,'Lecture de la variable APRES ,flat ',minval(flat),maxval(flat)
1058
1059  !CR:ATTENTION EN ATTENTE DE REGLER LA QUESTION DU PAS DE TEMPS INITIAL
1060  do t=1,ntime
1061     do k=1,nlevel
1062        temp(k,t)=temp0(k)
1063        qv(k,t)=qv0(k)
1064        ql(k,t)=ql0(k)
1065        qi(k,t)=qi0(k)
1066        u(k,t)=u0(k)
1067        v(k,t)=v0(k)
1068        tke(k,t)=tke0(k)
1069     enddo
1070  enddo
1071  !-----------------------------------------------------------------------
1072
1073  return
1074end subroutine old_read_SCM
1075!======================================================================
1076
1077!======================================================================
1078SUBROUTINE interp_case_time2(day,day1,annee_ref                &
1079     !    &         ,year_cas,day_cas,nt_cas,pdt_forc,nlev_cas      &
1080     ,nt_cas,nlev_cas                                       &
1081     ,ts_cas,ps_cas,plev_cas,t_cas,q_cas,u_cas,v_cas               &
1082     ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas           &
1083     ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas   &
1084     ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas,ustar_cas       &
1085     ,uw_cas,vw_cas,q1_cas,q2_cas                           &
1086     ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas       &
1087     ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas         &
1088     ,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas     &
1089     ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas       &
1090     ,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas    &
1091     ,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas    &
1092     ,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas)
1093
1094
1095  implicit none
1096
1097  !---------------------------------------------------------------------------------------
1098  ! Time interpolation of a 2D field to the timestep corresponding to day
1099  !
1100  ! day: current julian day (e.g. 717538.2)
1101  ! day1: first day of the simulation
1102  ! nt_cas: total nb of data in the forcing
1103  ! pdt_cas: total time interval (in sec) between 2 forcing data
1104  !---------------------------------------------------------------------------------------
1105
1106#include "compar1d.h"
1107#include "date_cas.h"
1108
1109  ! inputs:
1110  integer annee_ref
1111  integer nt_cas,nlev_cas
1112  real day, day1,day_cas
1113  real ts_cas(nt_cas),ps_cas(nt_cas)
1114  real plev_cas(nlev_cas,nt_cas)
1115  real t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas)
1116  real u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas)
1117  real ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas)
1118  real vitw_cas(nlev_cas,nt_cas)
1119  real du_cas(nlev_cas,nt_cas),hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas)
1120  real dv_cas(nlev_cas,nt_cas),hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas)
1121  real dt_cas(nlev_cas,nt_cas),ht_cas(nlev_cas,nt_cas),vt_cas(nlev_cas,nt_cas)
1122  real dtrad_cas(nlev_cas,nt_cas)
1123  real dq_cas(nlev_cas,nt_cas),hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas)
1124  real lat_cas(nt_cas)
1125  real sens_cas(nt_cas)
1126  real ustar_cas(nt_cas),uw_cas(nlev_cas,nt_cas),vw_cas(nlev_cas,nt_cas)
1127  real q1_cas(nlev_cas,nt_cas),q2_cas(nlev_cas,nt_cas)
1128
1129  ! outputs:
1130  real plev_prof_cas(nlev_cas)
1131  real t_prof_cas(nlev_cas),q_prof_cas(nlev_cas)
1132  real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas)
1133  real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas)
1134  real vitw_prof_cas(nlev_cas)
1135  real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas)
1136  real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas)
1137  real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas)
1138  real dtrad_prof_cas(nlev_cas)
1139  real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas)
1140  real lat_prof_cas,sens_prof_cas,ts_prof_cas,ustar_prof_cas
1141  real uw_prof_cas(nlev_cas),vw_prof_cas(nlev_cas),q1_prof_cas(nlev_cas),q2_prof_cas(nlev_cas)
1142  ! local:
1143  integer it_cas1, it_cas2,k
1144  real timeit,time_cas1,time_cas2,frac
1145
1146
1147  print*,'Check time',day1,day_ju_ini_cas,day_deb+1,pdt_cas
1148
1149  ! On teste si la date du cas AMMA est correcte.
1150  ! C est pour memoire car en fait les fichiers .def
1151  ! sont censes etre corrects.
1152  ! A supprimer a terme (MPL 20150623)
1153  !     if ((forcing_type.eq.10).and.(1.eq.0)) then
1154  ! Check that initial day of the simulation consistent with AMMA case:
1155  !      if (annee_ref.ne.2006) then
1156  !       print*,'Pour AMMA, annee_ref doit etre 2006'
1157  !       print*,'Changer annee_ref dans run.def'
1158  !       stop
1159  !      endif
1160  !      if (annee_ref.eq.2006 .and. day1.lt.day_cas) then
1161  !       print*,'AMMA a debute le 10 juillet 2006',day1,day_cas
1162  !       print*,'Changer dayref dans run.def'
1163  !       stop
1164  !      endif
1165  !      if (annee_ref.eq.2006 .and. day1.gt.day_cas+1) then
1166  !       print*,'AMMA a fini le 11 juillet'
1167  !       print*,'Changer dayref ou nday dans run.def'
1168  !       stop
1169  !      endif
1170  !      endif
1171
1172  ! Determine timestep relative to the 1st day:
1173  !       timeit=(day-day1)*86400.
1174  !       if (annee_ref.eq.1992) then
1175  !        timeit=(day-day_cas)*86400.
1176  !       else
1177  !        timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
1178  !       endif
1179  timeit=(day-day_ju_ini_cas)*86400
1180  print *,'day=',day
1181  print *,'day_ju_ini_cas=',day_ju_ini_cas
1182  print *,'pdt_cas=',pdt_cas
1183  print *,'timeit=',timeit
1184  print *,'nt_cas=',nt_cas
1185
1186  ! Determine the closest observation times:
1187  !       it_cas1=INT(timeit/pdt_cas)+1
1188  !       it_cas2=it_cas1 + 1
1189  !       time_cas1=(it_cas1-1)*pdt_cas
1190  !       time_cas2=(it_cas2-1)*pdt_cas
1191
1192  it_cas1=INT(timeit/pdt_cas)+1
1193  IF (it_cas1 .EQ. nt_cas) THEN
1194     it_cas2=it_cas1
1195  ELSE
1196     it_cas2=it_cas1 + 1
1197  ENDIF
1198  time_cas1=(it_cas1-1)*pdt_cas
1199  time_cas2=(it_cas2-1)*pdt_cas
1200  print *,'it_cas1,it_cas2,time_cas1,time_cas2=',it_cas1,it_cas2,time_cas1,time_cas2
1201
1202  if (it_cas1 .gt. nt_cas) then
1203     write(*,*) 'PB-stop: day, day_ju_ini_cas,it_cas1, it_cas2, timeit: '            &
1204          ,day,day_ju_ini_cas,it_cas1,it_cas2,timeit
1205     stop
1206  endif
1207
1208  ! time interpolation:
1209  IF (it_cas1 .EQ. it_cas2) THEN
1210     frac=0.
1211  ELSE
1212     frac=(time_cas2-timeit)/(time_cas2-time_cas1)
1213     frac=max(frac,0.0)
1214  ENDIF
1215
1216  lat_prof_cas = lat_cas(it_cas2)                                       &
1217       -frac*(lat_cas(it_cas2)-lat_cas(it_cas1))
1218  sens_prof_cas = sens_cas(it_cas2)                                     &
1219       -frac*(sens_cas(it_cas2)-sens_cas(it_cas1))
1220  ts_prof_cas = ts_cas(it_cas2)                                         &
1221       -frac*(ts_cas(it_cas2)-ts_cas(it_cas1))
1222  ustar_prof_cas = ustar_cas(it_cas2)                                   &
1223       -frac*(ustar_cas(it_cas2)-ustar_cas(it_cas1))
1224
1225  do k=1,nlev_cas
1226     plev_prof_cas(k) = plev_cas(k,it_cas2)                               &
1227          -frac*(plev_cas(k,it_cas2)-plev_cas(k,it_cas1))
1228     t_prof_cas(k) = t_cas(k,it_cas2)                               &
1229          -frac*(t_cas(k,it_cas2)-t_cas(k,it_cas1))
1230     q_prof_cas(k) = q_cas(k,it_cas2)                               &
1231          -frac*(q_cas(k,it_cas2)-q_cas(k,it_cas1))
1232     u_prof_cas(k) = u_cas(k,it_cas2)                               &
1233          -frac*(u_cas(k,it_cas2)-u_cas(k,it_cas1))
1234     v_prof_cas(k) = v_cas(k,it_cas2)                               &
1235          -frac*(v_cas(k,it_cas2)-v_cas(k,it_cas1))
1236     ug_prof_cas(k) = ug_cas(k,it_cas2)                               &
1237          -frac*(ug_cas(k,it_cas2)-ug_cas(k,it_cas1))
1238     vg_prof_cas(k) = vg_cas(k,it_cas2)                               &
1239          -frac*(vg_cas(k,it_cas2)-vg_cas(k,it_cas1))
1240     vitw_prof_cas(k) = vitw_cas(k,it_cas2)                               &
1241          -frac*(vitw_cas(k,it_cas2)-vitw_cas(k,it_cas1))
1242     du_prof_cas(k) = du_cas(k,it_cas2)                                   &
1243          -frac*(du_cas(k,it_cas2)-du_cas(k,it_cas1))
1244     hu_prof_cas(k) = hu_cas(k,it_cas2)                                   &
1245          -frac*(hu_cas(k,it_cas2)-hu_cas(k,it_cas1))
1246     vu_prof_cas(k) = vu_cas(k,it_cas2)                                   &
1247          -frac*(vu_cas(k,it_cas2)-vu_cas(k,it_cas1))
1248     dv_prof_cas(k) = dv_cas(k,it_cas2)                                   &
1249          -frac*(dv_cas(k,it_cas2)-dv_cas(k,it_cas1))
1250     hv_prof_cas(k) = hv_cas(k,it_cas2)                                   &
1251          -frac*(hv_cas(k,it_cas2)-hv_cas(k,it_cas1))
1252     vv_prof_cas(k) = vv_cas(k,it_cas2)                                   &
1253          -frac*(vv_cas(k,it_cas2)-vv_cas(k,it_cas1))
1254     dt_prof_cas(k) = dt_cas(k,it_cas2)                                   &
1255          -frac*(dt_cas(k,it_cas2)-dt_cas(k,it_cas1))
1256     ht_prof_cas(k) = ht_cas(k,it_cas2)                                   &
1257          -frac*(ht_cas(k,it_cas2)-ht_cas(k,it_cas1))
1258     vt_prof_cas(k) = vt_cas(k,it_cas2)                                   &
1259          -frac*(vt_cas(k,it_cas2)-vt_cas(k,it_cas1))
1260     dtrad_prof_cas(k) = dtrad_cas(k,it_cas2)                                   &
1261          -frac*(dtrad_cas(k,it_cas2)-dtrad_cas(k,it_cas1))
1262     dq_prof_cas(k) = dq_cas(k,it_cas2)                                   &
1263          -frac*(dq_cas(k,it_cas2)-dq_cas(k,it_cas1))
1264     hq_prof_cas(k) = hq_cas(k,it_cas2)                                   &
1265          -frac*(hq_cas(k,it_cas2)-hq_cas(k,it_cas1))
1266     vq_prof_cas(k) = vq_cas(k,it_cas2)                                   &
1267          -frac*(vq_cas(k,it_cas2)-vq_cas(k,it_cas1))
1268     uw_prof_cas(k) = uw_cas(k,it_cas2)                                   &
1269          -frac*(uw_cas(k,it_cas2)-uw_cas(k,it_cas1))
1270     vw_prof_cas(k) = vw_cas(k,it_cas2)                                   &
1271          -frac*(vw_cas(k,it_cas2)-vw_cas(k,it_cas1))
1272     q1_prof_cas(k) = q1_cas(k,it_cas2)                                   &
1273          -frac*(q1_cas(k,it_cas2)-q1_cas(k,it_cas1))
1274     q2_prof_cas(k) = q2_cas(k,it_cas2)                                   &
1275          -frac*(q2_cas(k,it_cas2)-q2_cas(k,it_cas1))
1276  enddo
1277
1278  return
1279END SUBROUTINE interp_case_time2
1280
1281!**********************************************************************************************
1282SUBROUTINE interp2_case_time(day,day1,annee_ref                           &
1283     !    &         ,year_cas,day_cas,nt_cas,pdt_forc,nlev_cas                         &
1284     ,nt_cas,nlev_cas                                                   &
1285     ,ts_cas,ps_cas,plev_cas,t_cas,theta_cas,thv_cas,thl_cas            &
1286     ,qv_cas,ql_cas,qi_cas,u_cas,v_cas                                  &
1287     ,ug_cas,vg_cas,vitw_cas,omega_cas,du_cas,hu_cas,vu_cas             &
1288     ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas               &
1289     ,dq_cas,hq_cas,vq_cas,dth_cas,hth_cas,vth_cas                      &
1290     ,lat_cas,sens_cas,ustar_cas                                        &
1291     ,uw_cas,vw_cas,q1_cas,q2_cas,tke_cas                               &
1292     !
1293     ,ts_prof_cas,plev_prof_cas,t_prof_cas,theta_prof_cas               &
1294     ,thv_prof_cas,thl_prof_cas,qv_prof_cas,ql_prof_cas,qi_prof_cas     &
1295     ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas                     &
1296     ,vitw_prof_cas,omega_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas  &
1297     ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas                   &
1298     ,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas                &
1299     ,hq_prof_cas,vq_prof_cas,dth_prof_cas,hth_prof_cas,vth_prof_cas    &
1300     ,lat_prof_cas,sens_prof_cas                                        &
1301     ,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas,tke_prof_cas)
1302
1303
1304  implicit none
1305
1306  !---------------------------------------------------------------------------------------
1307  ! Time interpolation of a 2D field to the timestep corresponding to day
1308  !
1309  ! day: current julian day (e.g. 717538.2)
1310  ! day1: first day of the simulation
1311  ! nt_cas: total nb of data in the forcing
1312  ! pdt_cas: total time interval (in sec) between 2 forcing data
1313  !---------------------------------------------------------------------------------------
1314
1315#include "compar1d.h"
1316#include "date_cas.h"
1317
1318  ! inputs:
1319  integer annee_ref
1320  integer nt_cas,nlev_cas
1321  real day, day1,day_cas
1322  real ts_cas(nt_cas),ps_cas(nt_cas)
1323  real plev_cas(nlev_cas,nt_cas)
1324  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)
1325  real qv_cas(nlev_cas,nt_cas),ql_cas(nlev_cas,nt_cas),qi_cas(nlev_cas,nt_cas)
1326  real u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas)
1327  real ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas)
1328  real vitw_cas(nlev_cas,nt_cas),omega_cas(nlev_cas,nt_cas)
1329  real du_cas(nlev_cas,nt_cas),hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas)
1330  real dv_cas(nlev_cas,nt_cas),hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas)
1331  real dt_cas(nlev_cas,nt_cas),ht_cas(nlev_cas,nt_cas),vt_cas(nlev_cas,nt_cas)
1332  real dth_cas(nlev_cas,nt_cas),hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas)
1333  real dtrad_cas(nlev_cas,nt_cas)
1334  real dq_cas(nlev_cas,nt_cas),hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas)
1335  real lat_cas(nt_cas),sens_cas(nt_cas),tke_cas(nt_cas)
1336  real ustar_cas(nt_cas),uw_cas(nlev_cas,nt_cas),vw_cas(nlev_cas,nt_cas)
1337  real q1_cas(nlev_cas,nt_cas),q2_cas(nlev_cas,nt_cas)
1338
1339  ! outputs:
1340  real plev_prof_cas(nlev_cas)
1341  real t_prof_cas(nlev_cas),theta_prof_cas(nlev_cas),thl_prof_cas(nlev_cas),thv_prof_cas(nlev_cas)
1342  real qv_prof_cas(nlev_cas),ql_prof_cas(nlev_cas),qi_prof_cas(nlev_cas)
1343  real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas)
1344  real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas)
1345  real vitw_prof_cas(nlev_cas),omega_prof_cas(nlev_cas)
1346  real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas)
1347  real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas)
1348  real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas)
1349  real dth_prof_cas(nlev_cas),hth_prof_cas(nlev_cas),vth_prof_cas(nlev_cas)
1350  real dtrad_prof_cas(nlev_cas)
1351  real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas)
1352  real lat_prof_cas,sens_prof_cas,tke_prof_cas,ts_prof_cas,ustar_prof_cas
1353  real uw_prof_cas(nlev_cas),vw_prof_cas(nlev_cas),q1_prof_cas(nlev_cas),q2_prof_cas(nlev_cas)
1354  ! local:
1355  integer it_cas1, it_cas2,k
1356  real timeit,time_cas1,time_cas2,frac
1357
1358
1359  print*,'Check time',day1,day_ju_ini_cas,day_deb+1,pdt_cas
1360  !       do k=1,nlev_cas
1361  !       print*,'debut de interp2_case_time, plev_cas=',k,plev_cas(k,1)
1362  !       enddo
1363
1364  ! On teste si la date du cas AMMA est correcte.
1365  ! C est pour memoire car en fait les fichiers .def
1366  ! sont censes etre corrects.
1367  ! A supprimer a terme (MPL 20150623)
1368  !     if ((forcing_type.eq.10).and.(1.eq.0)) then
1369  ! Check that initial day of the simulation consistent with AMMA case:
1370  !      if (annee_ref.ne.2006) then
1371  !       print*,'Pour AMMA, annee_ref doit etre 2006'
1372  !       print*,'Changer annee_ref dans run.def'
1373  !       stop
1374  !      endif
1375  !      if (annee_ref.eq.2006 .and. day1.lt.day_cas) then
1376  !       print*,'AMMA a debute le 10 juillet 2006',day1,day_cas
1377  !       print*,'Changer dayref dans run.def'
1378  !       stop
1379  !      endif
1380  !      if (annee_ref.eq.2006 .and. day1.gt.day_cas+1) then
1381  !       print*,'AMMA a fini le 11 juillet'
1382  !       print*,'Changer dayref ou nday dans run.def'
1383  !       stop
1384  !      endif
1385  !      endif
1386
1387  ! Determine timestep relative to the 1st day:
1388  !       timeit=(day-day1)*86400.
1389  !       if (annee_ref.eq.1992) then
1390  !        timeit=(day-day_cas)*86400.
1391  !       else
1392  !        timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
1393  !       endif
1394  timeit=(day-day_ju_ini_cas)*86400
1395  print *,'day=',day
1396  print *,'day_ju_ini_cas=',day_ju_ini_cas
1397  print *,'pdt_cas=',pdt_cas
1398  print *,'timeit=',timeit
1399  print *,'nt_cas=',nt_cas
1400
1401  ! Determine the closest observation times:
1402  !       it_cas1=INT(timeit/pdt_cas)+1
1403  !       it_cas2=it_cas1 + 1
1404  !       time_cas1=(it_cas1-1)*pdt_cas
1405  !       time_cas2=(it_cas2-1)*pdt_cas
1406
1407  it_cas1=INT(timeit/pdt_cas)+1
1408  IF (it_cas1 .EQ. nt_cas) THEN
1409     it_cas2=it_cas1
1410  ELSE
1411     it_cas2=it_cas1 + 1
1412  ENDIF
1413  time_cas1=(it_cas1-1)*pdt_cas
1414  time_cas2=(it_cas2-1)*pdt_cas
1415  print *,'timeit,pdt_cas,nt_cas=',timeit,pdt_cas,nt_cas
1416  print *,'it_cas1,it_cas2,time_cas1,time_cas2=',it_cas1,it_cas2,time_cas1,time_cas2
1417
1418  if (it_cas1 .gt. nt_cas) then
1419     write(*,*) 'PB-stop: day, day_ju_ini_cas,it_cas1, it_cas2, timeit: '            &
1420          ,day,day_ju_ini_cas,it_cas1,it_cas2,timeit
1421     stop
1422  endif
1423
1424  ! time interpolation:
1425  IF (it_cas1 .EQ. it_cas2) THEN
1426     frac=0.
1427  ELSE
1428     frac=(time_cas2-timeit)/(time_cas2-time_cas1)
1429     frac=max(frac,0.0)
1430  ENDIF
1431
1432  lat_prof_cas = lat_cas(it_cas2)                                   &
1433       -frac*(lat_cas(it_cas2)-lat_cas(it_cas1))
1434  sens_prof_cas = sens_cas(it_cas2)                                 &
1435       -frac*(sens_cas(it_cas2)-sens_cas(it_cas1))
1436  tke_prof_cas = tke_cas(it_cas2)                                   &
1437       -frac*(tke_cas(it_cas2)-tke_cas(it_cas1))
1438  ts_prof_cas = ts_cas(it_cas2)                                     &
1439       -frac*(ts_cas(it_cas2)-ts_cas(it_cas1))
1440  ustar_prof_cas = ustar_cas(it_cas2)                               &
1441       -frac*(ustar_cas(it_cas2)-ustar_cas(it_cas1))
1442
1443  do k=1,nlev_cas
1444     plev_prof_cas(k) = plev_cas(k,it_cas2)                           &     
1445          -frac*(plev_cas(k,it_cas2)-plev_cas(k,it_cas1))
1446     t_prof_cas(k) = t_cas(k,it_cas2)                                 &       
1447          -frac*(t_cas(k,it_cas2)-t_cas(k,it_cas1))
1448     print *,'k,frac,plev_cas1,plev_cas2=',k,frac,plev_cas(k,it_cas1),plev_cas(k,it_cas2)
1449     theta_prof_cas(k) = theta_cas(k,it_cas2)                         &                     
1450          -frac*(theta_cas(k,it_cas2)-theta_cas(k,it_cas1))
1451     thv_prof_cas(k) = thv_cas(k,it_cas2)                             &         
1452          -frac*(thv_cas(k,it_cas2)-thv_cas(k,it_cas1))
1453     thl_prof_cas(k) = thl_cas(k,it_cas2)                             &             
1454          -frac*(thl_cas(k,it_cas2)-thl_cas(k,it_cas1))
1455     qv_prof_cas(k) = qv_cas(k,it_cas2)                               &
1456          -frac*(qv_cas(k,it_cas2)-qv_cas(k,it_cas1))
1457     ql_prof_cas(k) = ql_cas(k,it_cas2)                               &
1458          -frac*(ql_cas(k,it_cas2)-ql_cas(k,it_cas1))
1459     qi_prof_cas(k) = qi_cas(k,it_cas2)                               &
1460          -frac*(qi_cas(k,it_cas2)-qi_cas(k,it_cas1))
1461     u_prof_cas(k) = u_cas(k,it_cas2)                                 &
1462          -frac*(u_cas(k,it_cas2)-u_cas(k,it_cas1))
1463     v_prof_cas(k) = v_cas(k,it_cas2)                                 &
1464          -frac*(v_cas(k,it_cas2)-v_cas(k,it_cas1))
1465     ug_prof_cas(k) = ug_cas(k,it_cas2)                               &
1466          -frac*(ug_cas(k,it_cas2)-ug_cas(k,it_cas1))
1467     vg_prof_cas(k) = vg_cas(k,it_cas2)                               &
1468          -frac*(vg_cas(k,it_cas2)-vg_cas(k,it_cas1))
1469     vitw_prof_cas(k) = vitw_cas(k,it_cas2)                           &
1470          -frac*(vitw_cas(k,it_cas2)-vitw_cas(k,it_cas1))
1471     omega_prof_cas(k) = omega_cas(k,it_cas2)                         &
1472          -frac*(omega_cas(k,it_cas2)-omega_cas(k,it_cas1))
1473     du_prof_cas(k) = du_cas(k,it_cas2)                               &
1474          -frac*(du_cas(k,it_cas2)-du_cas(k,it_cas1))
1475     hu_prof_cas(k) = hu_cas(k,it_cas2)                               &
1476          -frac*(hu_cas(k,it_cas2)-hu_cas(k,it_cas1))
1477     vu_prof_cas(k) = vu_cas(k,it_cas2)                               &
1478          -frac*(vu_cas(k,it_cas2)-vu_cas(k,it_cas1))
1479     dv_prof_cas(k) = dv_cas(k,it_cas2)                               &
1480          -frac*(dv_cas(k,it_cas2)-dv_cas(k,it_cas1))
1481     hv_prof_cas(k) = hv_cas(k,it_cas2)                               &
1482          -frac*(hv_cas(k,it_cas2)-hv_cas(k,it_cas1))
1483     vv_prof_cas(k) = vv_cas(k,it_cas2)                               &
1484          -frac*(vv_cas(k,it_cas2)-vv_cas(k,it_cas1))
1485     dt_prof_cas(k) = dt_cas(k,it_cas2)                               &
1486          -frac*(dt_cas(k,it_cas2)-dt_cas(k,it_cas1))
1487     ht_prof_cas(k) = ht_cas(k,it_cas2)                               &
1488          -frac*(ht_cas(k,it_cas2)-ht_cas(k,it_cas1))
1489     vt_prof_cas(k) = vt_cas(k,it_cas2)                               &
1490          -frac*(vt_cas(k,it_cas2)-vt_cas(k,it_cas1))
1491     dth_prof_cas(k) = dth_cas(k,it_cas2)                             &
1492          -frac*(dth_cas(k,it_cas2)-dth_cas(k,it_cas1))
1493     hth_prof_cas(k) = hth_cas(k,it_cas2)                             &
1494          -frac*(hth_cas(k,it_cas2)-hth_cas(k,it_cas1))
1495     vth_prof_cas(k) = vth_cas(k,it_cas2)                             &
1496          -frac*(vth_cas(k,it_cas2)-vth_cas(k,it_cas1))
1497     dtrad_prof_cas(k) = dtrad_cas(k,it_cas2)                         &
1498          -frac*(dtrad_cas(k,it_cas2)-dtrad_cas(k,it_cas1))
1499     dq_prof_cas(k) = dq_cas(k,it_cas2)                               &
1500          -frac*(dq_cas(k,it_cas2)-dq_cas(k,it_cas1))
1501     hq_prof_cas(k) = hq_cas(k,it_cas2)                               &
1502          -frac*(hq_cas(k,it_cas2)-hq_cas(k,it_cas1))
1503     vq_prof_cas(k) = vq_cas(k,it_cas2)                               &
1504          -frac*(vq_cas(k,it_cas2)-vq_cas(k,it_cas1))
1505     uw_prof_cas(k) = uw_cas(k,it_cas2)                                &
1506          -frac*(uw_cas(k,it_cas2)-uw_cas(k,it_cas1))
1507     vw_prof_cas(k) = vw_cas(k,it_cas2)                                &
1508          -frac*(vw_cas(k,it_cas2)-vw_cas(k,it_cas1))
1509     q1_prof_cas(k) = q1_cas(k,it_cas2)                                &
1510          -frac*(q1_cas(k,it_cas2)-q1_cas(k,it_cas1))
1511     q2_prof_cas(k) = q2_cas(k,it_cas2)                                &
1512          -frac*(q2_cas(k,it_cas2)-q2_cas(k,it_cas1))
1513  enddo
1514
1515  return
1516END SUBROUTINE interp2_case_time
1517
1518!**********************************************************************************************
1519
Note: See TracBrowser for help on using the repository browser.