source: LMDZ6/branches/Amaury_dev/libf/phylmd/dyn1d/mod_1D_cases_read2.F90

Last change on this file was 5295, checked in by abarral, 2 months ago

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