source: LMDZ6/branches/Amaury_dev/libf/phylmd/dyn1d/mod_1D_cases_read_std.F90 @ 5158

Last change on this file since 5158 was 5158, checked in by abarral, 7 weeks ago

Add missing klon on strataer_emiss_mod.F90
Correct various missing explicit declarations
Replace tabs by spaces (tabs are not part of the fortran charset)
Continue cleaning modules
Removed unused arguments and variables
File size: 59.4 KB
Line 
1! $Id: mod_1D_cases_read.F90 2373 2015-10-13 17:28:01Z jyg $
2
3MODULE mod_1D_cases_read_std
4  USE netcdf, ONLY: nf90_noerr, nf90_inq_varid, nf90_inq_dimid, nf90_inquire_dimension, nf90_open, nf90_nowrite, &
5          nf90_strerror, nf90_get_var
6
7  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
8  !Declarations specifiques au cas standard
9  CHARACTER*80 :: fich_cas
10  ! Discr?tisation
11  INTEGER nlev_cas, nt_cas
12
13
14  !profils environnementaux
15  REAL, ALLOCATABLE :: plev_cas(:, :), plevh_cas(:)
16  REAL, ALLOCATABLE :: ap_cas(:), bp_cas(:)
17
18  REAL, ALLOCATABLE :: z_cas(:, :), zh_cas(:)
19  REAL, ALLOCATABLE :: t_cas(:, :), q_cas(:, :), qv_cas(:, :), ql_cas(:, :), qi_cas(:, :), rh_cas(:, :)
20  REAL, ALLOCATABLE :: th_cas(:, :), thv_cas(:, :), thl_cas(:, :), rv_cas(:, :)
21  REAL, ALLOCATABLE :: u_cas(:, :), v_cas(:, :), vitw_cas(:, :), omega_cas(:, :), tke_cas(:, :)
22
23  !forcing
24  REAL, ALLOCATABLE :: ht_cas(:, :), vt_cas(:, :), dt_cas(:, :), dtrad_cas(:, :)
25  REAL, ALLOCATABLE :: hth_cas(:, :), vth_cas(:, :), dth_cas(:, :)
26  REAL, ALLOCATABLE :: hq_cas(:, :), vq_cas(:, :), dq_cas(:, :)
27  REAL, ALLOCATABLE :: hr_cas(:, :), vr_cas(:, :), dr_cas(:, :)
28  REAL, ALLOCATABLE :: hu_cas(:, :), vu_cas(:, :), du_cas(:, :)
29  REAL, ALLOCATABLE :: hv_cas(:, :), vv_cas(:, :), dv_cas(:, :)
30  REAL, ALLOCATABLE :: ug_cas(:, :), vg_cas(:, :)
31  REAL, ALLOCATABLE :: temp_nudg_cas(:, :), qv_nudg_cas(:, :), u_nudg_cas(:, :), v_nudg_cas(:, :)
32  REAL, ALLOCATABLE :: invtau_temp_nudg_cas(:, :), invtau_qv_nudg_cas(:, :), invtau_u_nudg_cas(:, :), invtau_v_nudg_cas(:, :)
33  REAL, ALLOCATABLE :: lat_cas(:), sens_cas(:), tskin_cas(:), ts_cas(:), ps_cas(:), ustar_cas(:)
34  REAL, ALLOCATABLE :: uw_cas(:, :), vw_cas(:, :), q1_cas(:, :), q2_cas(:, :), tkes_cas(:)
35
36  !champs interpoles
37  REAL, ALLOCATABLE :: plev_prof_cas(:)
38  REAL, ALLOCATABLE :: t_prof_cas(:)
39  REAL, ALLOCATABLE :: theta_prof_cas(:)
40  REAL, ALLOCATABLE :: thl_prof_cas(:)
41  REAL, ALLOCATABLE :: thv_prof_cas(:)
42  REAL, ALLOCATABLE :: q_prof_cas(:)
43  REAL, ALLOCATABLE :: qv_prof_cas(:)
44  REAL, ALLOCATABLE :: ql_prof_cas(:)
45  REAL, ALLOCATABLE :: qi_prof_cas(:)
46  REAL, ALLOCATABLE :: rh_prof_cas(:)
47  REAL, ALLOCATABLE :: rv_prof_cas(:)
48  REAL, ALLOCATABLE :: u_prof_cas(:)
49  REAL, ALLOCATABLE :: v_prof_cas(:)
50  REAL, ALLOCATABLE :: vitw_prof_cas(:)
51  REAL, ALLOCATABLE :: omega_prof_cas(:)
52  REAL, ALLOCATABLE :: tke_prof_cas(:)
53  REAL, ALLOCATABLE :: ug_prof_cas(:)
54  REAL, ALLOCATABLE :: vg_prof_cas(:)
55  REAL, ALLOCATABLE :: temp_nudg_prof_cas(:), qv_nudg_prof_cas(:), u_nudg_prof_cas(:), v_nudg_prof_cas(:)
56  REAL, ALLOCATABLE :: invtau_temp_nudg_prof_cas(:), invtau_qv_nudg_prof_cas(:), invtau_u_nudg_prof_cas(:), invtau_v_nudg_prof_cas(:)
57
58  REAL, ALLOCATABLE :: ht_prof_cas(:)
59  REAL, ALLOCATABLE :: hth_prof_cas(:)
60  REAL, ALLOCATABLE :: hq_prof_cas(:)
61  REAL, ALLOCATABLE :: vt_prof_cas(:)
62  REAL, ALLOCATABLE :: vth_prof_cas(:)
63  REAL, ALLOCATABLE :: vq_prof_cas(:)
64  REAL, ALLOCATABLE :: dt_prof_cas(:)
65  REAL, ALLOCATABLE :: dth_prof_cas(:)
66  REAL, ALLOCATABLE :: dtrad_prof_cas(:)
67  REAL, ALLOCATABLE :: dq_prof_cas(:)
68  REAL, ALLOCATABLE :: hu_prof_cas(:)
69  REAL, ALLOCATABLE :: hv_prof_cas(:)
70  REAL, ALLOCATABLE :: vu_prof_cas(:)
71  REAL, ALLOCATABLE :: vv_prof_cas(:)
72  REAL, ALLOCATABLE :: du_prof_cas(:)
73  REAL, ALLOCATABLE :: dv_prof_cas(:)
74  REAL, ALLOCATABLE :: uw_prof_cas(:)
75  REAL, ALLOCATABLE :: vw_prof_cas(:)
76  REAL, ALLOCATABLE :: q1_prof_cas(:)
77  REAL, ALLOCATABLE :: q2_prof_cas(:)
78
79  REAL o3_cas, lat_prof_cas, sens_prof_cas, ts_prof_cas, tskin_prof_cas, ps_prof_cas, ustar_prof_cas, tkes_prof_cas
80  REAL orog_cas, albedo_cas, emiss_cas, q_skin_cas, mom_rough, heat_rough, rugos_cas, sand_cas, clay_cas
81
82
83CONTAINS
84
85
86  !**********************************************************************************************
87  SUBROUTINE read_SCM_cas
88    USE lmdz_date_cas, ONLY: year_ini_cas, mth_ini_cas, day_deb, heure_ini_cas, pdt_cas, day_ju_ini_cas
89
90    IMPLICIT NONE
91
92    INTEGER nid, rid, ierr
93    INTEGER ii, jj, timeid
94    REAL, ALLOCATABLE :: time_val(:)
95
96    fich_cas = 'cas.nc'
97    PRINT*, 'fich_cas ', fich_cas
98    ierr = nf90_open(fich_cas, nf90_nowrite, nid)
99    PRINT*, 'fich_cas,nf90_nowrite,nid ', fich_cas, nf90_nowrite, nid
100    IF (ierr/=nf90_noerr) THEN
101      WRITE(*, *) 'ERROR: GROS Pb opening forcings nc file '
102      WRITE(*, *) nf90_strerror(ierr)
103      stop ""
104    endif
105    !.......................................................................
106    ierr = nf90_inq_dimid(nid, 'lat', rid)
107    IF (ierr/=nf90_noerr) THEN
108      PRINT*, 'Oh probleme lecture dimension lat'
109    ENDIF
110    ierr = nf90_inquire_dimension(nid, rid, len = ii)
111    PRINT*, 'OK1 read_SCM_cas: nid,rid,lat', nid, rid, ii
112    !.......................................................................
113    ierr = nf90_inq_dimid(nid, 'lon', rid)
114    IF (ierr/=nf90_noerr) THEN
115      PRINT*, 'Oh probleme lecture dimension lon'
116    ENDIF
117    ierr = nf90_inquire_dimension(nid, rid, len = jj)
118    PRINT*, 'OK2 read_SCM_cas: nid,rid,lat', nid, rid, jj
119    !.......................................................................
120    ierr = nf90_inq_dimid(nid, 'lev', rid)
121    IF (ierr/=nf90_noerr) THEN
122      PRINT*, 'Oh probleme lecture dimension nlev'
123    ENDIF
124    ierr = nf90_inquire_dimension(nid, rid, len = nlev_cas)
125    PRINT*, 'OK3 read_SCM_cas: nid,rid,nlev_cas', nid, rid, nlev_cas
126    IF (.NOT. (nlev_cas > 10 .AND. nlev_cas < 200000)) THEN
127      PRINT*, 'Valeur de nlev_cas peu probable'
128      STOP
129    ENDIF
130    !.......................................................................
131    ierr = nf90_inq_dimid(nid, 'time', rid)
132    nt_cas = 0
133    IF (ierr/=nf90_noerr) THEN
134      stop 'Oh probleme lecture dimension time'
135    ENDIF
136    ierr = nf90_inquire_dimension(nid, rid, len = nt_cas)
137    PRINT*, 'OK4 read_SCM_cas: nid,rid,nt_cas', nid, rid, nt_cas
138    ! Lecture de l'axe des temps
139    PRINT*, 'LECTURE DU TEMPS'
140    ierr = nf90_inq_varid(nid, 'time', timeid)
141    IF(ierr/=nf90_noerr) THEN
142      print *, 'Variable time manquante dans cas.nc:'
143      ierr = nf90_noerr
144    else
145      allocate(time_val(nt_cas))
146      ierr = nf90_get_var(nid, timeid, time_val)
147      IF(ierr/=nf90_noerr) THEN
148        print *, 'A Pb a la lecture de time cas.nc: '
149      endif
150    endif
151    IF (nt_cas>1) THEN
152      pdt_cas = time_val(2) - time_val(1)
153    ELSE
154      pdt_cas = 0.
155    ENDIF
156
157
158    !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
159    !profils moyens:
160    allocate(plev_cas(nlev_cas, nt_cas), plevh_cas(nlev_cas + 1))
161    allocate(z_cas(nlev_cas, nt_cas), zh_cas(nlev_cas + 1))
162    allocate(ap_cas(nlev_cas + 1), bp_cas(nlev_cas + 1))
163    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), &
164            qi_cas(nlev_cas, nt_cas), rh_cas(nlev_cas, nt_cas))
165    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))
166    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))
167    allocate(tke_cas(nlev_cas, nt_cas))
168    !forcing
169    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))
170    allocate(hq_cas(nlev_cas, nt_cas), vq_cas(nlev_cas, nt_cas), dq_cas(nlev_cas, nt_cas))
171    allocate(hth_cas(nlev_cas, nt_cas), vth_cas(nlev_cas, nt_cas), dth_cas(nlev_cas, nt_cas))
172    allocate(hr_cas(nlev_cas, nt_cas), vr_cas(nlev_cas, nt_cas), dr_cas(nlev_cas, nt_cas))
173    allocate(hu_cas(nlev_cas, nt_cas), vu_cas(nlev_cas, nt_cas), du_cas(nlev_cas, nt_cas))
174    allocate(hv_cas(nlev_cas, nt_cas), vv_cas(nlev_cas, nt_cas), dv_cas(nlev_cas, nt_cas))
175    allocate(ug_cas(nlev_cas, nt_cas), vg_cas(nlev_cas, nt_cas))
176    allocate(temp_nudg_cas(nlev_cas, nt_cas), qv_nudg_cas(nlev_cas, nt_cas))
177    allocate(u_nudg_cas(nlev_cas, nt_cas), v_nudg_cas(nlev_cas, nt_cas))
178    allocate(invtau_temp_nudg_cas(nlev_cas, nt_cas), invtau_qv_nudg_cas(nlev_cas, nt_cas))
179    allocate(invtau_u_nudg_cas(nlev_cas, nt_cas), invtau_v_nudg_cas(nlev_cas, nt_cas))
180    allocate(lat_cas(nt_cas), sens_cas(nt_cas), ts_cas(nt_cas), tskin_cas(nt_cas), ps_cas(nt_cas), ustar_cas(nt_cas), tkes_cas(nt_cas))
181    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))
182
183
184
185    !champs interpoles
186    allocate(plev_prof_cas(nlev_cas))
187    allocate(t_prof_cas(nlev_cas))
188    allocate(theta_prof_cas(nlev_cas))
189    allocate(thl_prof_cas(nlev_cas))
190    allocate(thv_prof_cas(nlev_cas))
191    allocate(q_prof_cas(nlev_cas))
192    allocate(qv_prof_cas(nlev_cas))
193    allocate(ql_prof_cas(nlev_cas))
194    allocate(qi_prof_cas(nlev_cas))
195    allocate(rh_prof_cas(nlev_cas))
196    allocate(rv_prof_cas(nlev_cas))
197    allocate(u_prof_cas(nlev_cas))
198    allocate(v_prof_cas(nlev_cas))
199    allocate(vitw_prof_cas(nlev_cas))
200    allocate(omega_prof_cas(nlev_cas))
201    allocate(tke_prof_cas(nlev_cas))
202    allocate(ug_prof_cas(nlev_cas))
203    allocate(vg_prof_cas(nlev_cas))
204    allocate(temp_nudg_prof_cas(nlev_cas), qv_nudg_prof_cas(nlev_cas))
205    allocate(u_nudg_prof_cas(nlev_cas), v_nudg_prof_cas(nlev_cas))
206    allocate(invtau_temp_nudg_prof_cas(nlev_cas), invtau_qv_nudg_prof_cas(nlev_cas))
207    allocate(invtau_u_nudg_prof_cas(nlev_cas), invtau_v_nudg_prof_cas(nlev_cas))
208    allocate(ht_prof_cas(nlev_cas))
209    allocate(hth_prof_cas(nlev_cas))
210    allocate(hq_prof_cas(nlev_cas))
211    allocate(hu_prof_cas(nlev_cas))
212    allocate(hv_prof_cas(nlev_cas))
213    allocate(vt_prof_cas(nlev_cas))
214    allocate(vth_prof_cas(nlev_cas))
215    allocate(vq_prof_cas(nlev_cas))
216    allocate(vu_prof_cas(nlev_cas))
217    allocate(vv_prof_cas(nlev_cas))
218    allocate(dt_prof_cas(nlev_cas))
219    allocate(dth_prof_cas(nlev_cas))
220    allocate(dtrad_prof_cas(nlev_cas))
221    allocate(dq_prof_cas(nlev_cas))
222    allocate(du_prof_cas(nlev_cas))
223    allocate(dv_prof_cas(nlev_cas))
224    allocate(uw_prof_cas(nlev_cas))
225    allocate(vw_prof_cas(nlev_cas))
226    allocate(q1_prof_cas(nlev_cas))
227    allocate(q2_prof_cas(nlev_cas))
228
229    PRINT*, 'Allocations OK'
230    CALL read_SCM (nid, nlev_cas, nt_cas, &
231            ap_cas, bp_cas, z_cas, plev_cas, zh_cas, plevh_cas, t_cas, th_cas, thv_cas, thl_cas, qv_cas, &
232            ql_cas, qi_cas, rh_cas, rv_cas, u_cas, v_cas, vitw_cas, omega_cas, tke_cas, ug_cas, vg_cas, &
233            temp_nudg_cas, qv_nudg_cas, u_nudg_cas, v_nudg_cas, &
234            invtau_temp_nudg_cas, invtau_qv_nudg_cas, invtau_u_nudg_cas, invtau_v_nudg_cas, &
235            du_cas, hu_cas, vu_cas, &
236            dv_cas, hv_cas, vv_cas, dt_cas, ht_cas, vt_cas, dq_cas, hq_cas, vq_cas, dth_cas, hth_cas, vth_cas, &
237            dr_cas, hr_cas, vr_cas, dtrad_cas, sens_cas, lat_cas, ts_cas, tskin_cas, ps_cas, ustar_cas, tkes_cas, &
238            uw_cas, vw_cas, q1_cas, q2_cas, orog_cas, albedo_cas, emiss_cas, q_skin_cas, mom_rough, heat_rough, &
239            o3_cas, rugos_cas, clay_cas, sand_cas)
240    PRINT*, 'read_SCM cas OK'
241    DO ii = 1, nlev_cas
242      PRINT*, 'apres read_SCM_cas, plev_cas=', ii, plev_cas(ii, 1)
243      !PRINT*,'apres read_SCM, plev_cas=',ii,omega_cas(ii,nt_cas/2+1)
244    enddo
245
246  END SUBROUTINE read_SCM_cas
247
248
249  !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
250  SUBROUTINE deallocate2_1D_cases
251    !profils environnementaux:
252    deallocate(plev_cas, plevh_cas)
253
254    deallocate(z_cas, zh_cas)
255    deallocate(ap_cas, bp_cas)
256    deallocate(t_cas, q_cas, qv_cas, ql_cas, qi_cas, rh_cas)
257    deallocate(th_cas, thl_cas, thv_cas, rv_cas)
258    deallocate(u_cas, v_cas, vitw_cas, omega_cas, tke_cas)
259
260    !forcing
261    deallocate(ht_cas, vt_cas, dt_cas, dtrad_cas)
262    deallocate(hq_cas, vq_cas, dq_cas)
263    deallocate(hth_cas, vth_cas, dth_cas)
264    deallocate(hr_cas, vr_cas, dr_cas)
265    deallocate(hu_cas, vu_cas, du_cas)
266    deallocate(hv_cas, vv_cas, dv_cas)
267    deallocate(ug_cas)
268    deallocate(vg_cas)
269    deallocate(lat_cas, sens_cas, tskin_cas, ts_cas, ps_cas, ustar_cas, tkes_cas, uw_cas, vw_cas, q1_cas, q2_cas)
270
271    !champs interpoles
272    deallocate(plev_prof_cas)
273    deallocate(t_prof_cas)
274    deallocate(theta_prof_cas)
275    deallocate(thl_prof_cas)
276    deallocate(thv_prof_cas)
277    deallocate(q_prof_cas)
278    deallocate(qv_prof_cas)
279    deallocate(ql_prof_cas)
280    deallocate(qi_prof_cas)
281    deallocate(rh_prof_cas)
282    deallocate(rv_prof_cas)
283    deallocate(u_prof_cas)
284    deallocate(v_prof_cas)
285    deallocate(vitw_prof_cas)
286    deallocate(omega_prof_cas)
287    deallocate(tke_prof_cas)
288    deallocate(ug_prof_cas)
289    deallocate(vg_prof_cas)
290    deallocate(temp_nudg_prof_cas, qv_nudg_prof_cas, u_nudg_prof_cas, v_nudg_prof_cas)
291    deallocate(invtau_temp_nudg_prof_cas, invtau_qv_nudg_prof_cas, invtau_u_nudg_prof_cas, invtau_v_nudg_prof_cas)
292    deallocate(ht_prof_cas)
293    deallocate(hq_prof_cas)
294    deallocate(hu_prof_cas)
295    deallocate(hv_prof_cas)
296    deallocate(vt_prof_cas)
297    deallocate(vq_prof_cas)
298    deallocate(vu_prof_cas)
299    deallocate(vv_prof_cas)
300    deallocate(dt_prof_cas)
301    deallocate(dtrad_prof_cas)
302    deallocate(dq_prof_cas)
303    deallocate(du_prof_cas)
304    deallocate(dv_prof_cas)
305    deallocate(t_prof_cas)
306    deallocate(u_prof_cas)
307    deallocate(v_prof_cas)
308    deallocate(uw_prof_cas)
309    deallocate(vw_prof_cas)
310    deallocate(q1_prof_cas)
311    deallocate(q2_prof_cas)
312
313  END SUBROUTINE deallocate2_1D_cases
314
315
316  !=====================================================================
317  SUBROUTINE read_SCM(nid, nlevel, ntime, &
318          ap, bp, zz, pp, zzh, pph, temp, theta, thv, thl, qv, ql, qi, rh, rv, u, v, vitw, omega, tke, ug, vg, &
319          temp_nudg, qv_nudg, u_nudg, v_nudg, &
320          invtau_temp_nudg, invtau_qv_nudg, invtau_u_nudg, invtau_v_nudg, &
321          du, hu, vu, dv, hv, vv, dt, ht, vt, dq, hq, vq, &
322          dth, hth, vth, dr, hr, vr, dtrad, sens, flat, ts, tskin, ps, ustar, tkes, uw, vw, q1, q2, &
323          orog_cas, albedo_cas, emiss_cas, q_skin_cas, mom_rough, &
324          heat_rough, o3_cas, rugos_cas, clay_cas, sand_cas)
325
326    !program reading forcing of the case study
327    USE lmdz_compar1d
328
329    IMPLICIT NONE
330
331    INTEGER ntime, nlevel, k, t
332
333    REAL ap(nlevel + 1), bp(nlevel + 1)
334    REAL zz(nlevel, ntime), zzh(nlevel + 1)
335    REAL pp(nlevel, ntime), pph(nlevel + 1)
336    !profils initiaux
337    REAL temp0(nlevel), qv0(nlevel), ql0(nlevel), qi0(nlevel), u0(nlevel), v0(nlevel), tke0(nlevel)
338    REAL pp0(nlevel)
339    REAL temp(nlevel, ntime), qv(nlevel, ntime), ql(nlevel, ntime), qi(nlevel, ntime), rh(nlevel, ntime)
340    REAL theta(nlevel, ntime), thv(nlevel, ntime), thl(nlevel, ntime), rv(nlevel, ntime)
341    REAL u(nlevel, ntime), v(nlevel, ntime), tkes(ntime)
342    REAL temp_nudg(nlevel, ntime), qv_nudg(nlevel, ntime), u_nudg(nlevel, ntime), v_nudg(nlevel, ntime)
343    REAL invtau_temp_nudg(nlevel, ntime), invtau_qv_nudg(nlevel, ntime), invtau_u_nudg(nlevel, ntime), invtau_v_nudg(nlevel, ntime)
344    REAL ug(nlevel, ntime), vg(nlevel, ntime)
345    REAL vitw(nlevel, ntime), omega(nlevel, ntime), tke(nlevel, ntime)
346    REAL du(nlevel, ntime), hu(nlevel, ntime), vu(nlevel, ntime)
347    REAL dv(nlevel, ntime), hv(nlevel, ntime), vv(nlevel, ntime)
348    REAL dt(nlevel, ntime), ht(nlevel, ntime), vt(nlevel, ntime)
349    REAL dtrad(nlevel, ntime)
350    REAL dq(nlevel, ntime), hq(nlevel, ntime), vq(nlevel, ntime)
351    REAL dth(nlevel, ntime), hth(nlevel, ntime), vth(nlevel, ntime), hthl(nlevel, ntime)
352    REAL dr(nlevel, ntime), hr(nlevel, ntime), vr(nlevel, ntime)
353    REAL flat(ntime), sens(ntime), ustar(ntime)
354    REAL uw(nlevel, ntime), vw(nlevel, ntime), q1(nlevel, ntime), q2(nlevel, ntime)
355    REAL ts(ntime), tskin(ntime), ps(ntime)
356    REAL orog_cas, albedo_cas, emiss_cas, q_skin_cas, mom_rough, heat_rough, o3_cas, rugos_cas, clay_cas, sand_cas
357    REAL apbp(nlevel + 1), resul(nlevel, ntime), resul1(nlevel), resul2(ntime), resul3
358
359    INTEGER nid, ierr, ierr1, ierr2, rid, i, int_test
360    INTEGER nbvar3d
361    parameter(nbvar3d = 78)
362    INTEGER var3didin(nbvar3d), missing_var(nbvar3d)
363    CHARACTER*13 name_var(1:nbvar3d)
364
365
366    !      data name_var/ &
367    !     ! coordonnees pression (n+1 niveaux) #4
368    !     & 'coor_par_a','coor_par_b','height_h','pressure_h',& ! #1-#4
369    !     ! coordonnees pression (n niveaux) #8
370    !     &'temp','qv','ql','qi','u','v','tke','pressure',& ! #5-#12
371    !     ! coordonnees pression + temps #42
372    !     &'w','omega','ug','vg','uadv','uadvh','uadvv','vadv','vadvh','vadvv','temp_adv','tadvh','tadvv',& !  #13 - #25
373    !     &'qv_adv','qvadvh','qvadvv','thadv','thadvh','thadvv','thladvh',                             & ! #26 - #32
374    !     & 'radv','radvh','radvv','radcool','q1','q2','ustress','vstress',                           & ! #33 - #40
375    !     & 'rh','temp_nudging','qv_nudging','u_nudging','v_nudging',                                       & ! #41-45
376    !     &'height_f','pressure_forc','tempt','theta','thv','thl','qvt','qlt','qit','rv','ut','vt',   & ! #46-58
377    !     ! coordonnees temps #12
378    !     &'tkes','sfc_sens_flx','sfc_lat_flx','ts','ps','ustar',&
379    !     &'orog','albedo','emiss','t_skin','q_skin','mom_rough','heat_rough',&
380    !     ! scalaires #4
381    !     &'o3','rugos','clay','sand'/
382
383    data name_var/ &
384            ! coordonnees pression (n+1 niveaux) #4
385            'coor_par_a', 'coor_par_b', 'zf', 'pressure_h', & ! #1-#4
386            ! coordonnees pression (n niveaux) #8
387            'ta', 'qv', 'ql', 'qi', 'ua', 'va', 'tke', 'pa', & ! #5-#12
388            ! coordonnees pression + temps #46
389            'wa', 'wap', 'ug', 'vg', 'tnua_adv', 'tnua_advh', 'tnua_advv', 'tnva_adv', 'tnva_advh', 'tnva_advv', 'tnta_adv', 'tnta_advh', & !  #13 - #25
390            'tnta_advv', 'tnqv_adv', 'tnqv_advh', 'tnqv_advv', 'thadv', 'thadvh', 'thadvv', 'thladvh', & ! #26 - #32
391            'radv', 'radvh', 'radvv', 'tnta_rad', 'q1', 'q2', 'ustress', 'vstress', & ! #33 - #40
392            'rh', 'ta_nud', 'qv_nud', 'ua_nud', 'va_nud', & ! #41-45
393            'zh_forc', 'pa_forc', 'tat', 'thetat', 'thetavt', 'thetalt', 'qvt', 'qlt', 'qit', 'rvt', 'uat', 'vat', & ! #46-57
394            'nudging_constant_ta', 'nudging_constant_qv', 'nudging_constant_ua', 'nudging_constant_va', & ! # 58-61
395            ! coordonnees temps #12
396            'tkes', 'hfss', 'hfls', 'ts_forc', 'tskin', 'ps_forc', 'ustar', &                     ! 62-68
397            ! scalaires
398            'orog', 'albedo', 'emiss', 'q_skin', 'z0', 'z0h', &                    ! 69-74
399            'O3', 'rugos', 'clay', 'sand'/                                                      ! 75-78
400
401
402    !-----------------------------------------------------------------------
403    ! First check that we are using a version > v2 of the 1D standard format
404    ! use the difference between 'temp' (old version) and 'ta' (new version)
405    !-----------------------------------------------------------------------
406
407    ierr = nf90_inq_varid(nid, 'ta', int_test)
408    IF(ierr/=nf90_noerr) THEN
409      PRINT*, '++++++++++++++++++++++++++++++'
410      PRINT*, 'variable ta missing in cas.nc '
411      PRINT*, 'You are probably using an obsolete version of the 1D cases'
412      PRINT*, 'please dowload the last version of the 1D archive from https://lmdz.lmd.jussieu.fr/pub/'
413      PRINT*, '++++++++++++++++++++++++++++++'
414      CALL abort_gcm ('mod_1D_cases_read_std', 'bad version of 1D directory', 0)
415    endif
416
417    !-----------------------------------------------------------------------
418    ! Checking availability of variable #i in the cas.nc file
419    !     missing_var=1 if the variable is missing
420    !-----------------------------------------------------------------------
421
422    DO i = 1, nbvar3d
423      missing_var(i) = 0.
424      ierr = nf90_inq_varid(nid, name_var(i), var3didin(i))
425      PRINT*, 'name_var(i)', name_var(i), var3didin(i)
426      IF(ierr/=nf90_noerr) THEN
427        print *, 'Variable manquante dans cas.nc:', i, name_var(i)
428        ierr = nf90_noerr
429        missing_var(i) = 1
430      else
431
432        !-----------------------------------------------------------------------
433        ! Activating keys depending on the presence of specific variables in cas.nc
434        !-----------------------------------------------------------------------
435        IF (1 == 1) THEN
436          ! A MODIFIER: il faudrait dire nudging_temp mais faut le declarer dans compar1d.h etc...
437          !           if ( name_var(i) == 'temp_nudging' .AND. nint(nudging_t)==0) stop 'Nudging inconsistency temp'
438          IF (name_var(i) == 'qv_nud' .AND. nint(nudging_qv)==0) stop 'Nudging inconsistency qv'
439          IF (name_var(i) == 'ua_nud' .AND. nint(nudging_u)==0) stop 'Nudging inconsistency u'
440          IF (name_var(i) == 'va_nud' .AND. nint(nudging_v)==0) stop 'Nudging inconsistency v'
441        ELSE
442          PRINT*, 'GUIDAGE : CONSISTENCY CHECK DEACTIVATED FOR TESTS of SANDU/REF'
443        ENDIF
444
445        !-----------------------------------------------------------------------
446        ! Reading variables 1D (N+1) vertical variables (nlevelp1,lat,lon)
447        !-----------------------------------------------------------------------
448        IF(i<=4) THEN
449          ierr = nf90_get_var(nid, var3didin(i), apbp)
450          print *, 'read_SCM(apbp), on a lu ', i, name_var(i)
451          IF(ierr/=nf90_noerr) THEN
452            print *, 'B Pb a la lecture de cas.nc: ', name_var(i)
453            stop "getvarup"
454          endif
455
456          !-----------------------------------------------------------------------
457          !  Reading 1D (N) vertical varialbes    (nlevel,lat,lon)
458          !-----------------------------------------------------------------------
459        else IF(i>4.AND.i<=12) THEN
460          ierr = nf90_get_var(nid, var3didin(i), resul1)
461          print *, 'read_SCM(resul1), on a lu ', i, name_var(i)
462          IF(ierr/=nf90_noerr) THEN
463            print *, 'C Pb a la lecture de cas.nc: ', name_var(i)
464            stop "getvarup"
465          endif
466          PRINT*, 'Lecture de la variable #i ', i, name_var(i), minval(resul1), maxval(resul1)
467
468          !-----------------------------------------------------------------------
469          !  Reading 2D tim-vertical variables  (time,nlevel,lat,lon)
470          !  TBD : seems to be the same as above.
471          !-----------------------------------------------------------------------
472        else IF(i>12.AND.i<=61) THEN
473          ierr = nf90_get_var(nid, var3didin(i), resul)
474          print *, 'read_SCM(resul), on a lu ', i, name_var(i)
475          IF(ierr/=nf90_noerr) THEN
476            print *, 'D Pb a la lecture de cas.nc: ', name_var(i)
477            stop "getvarup"
478          endif
479          PRINT*, 'Lecture de la variable #i ', i, name_var(i), minval(resul), maxval(resul)
480
481          !-----------------------------------------------------------------------
482          !  Reading 1D time variables (time,lat,lon)
483          !-----------------------------------------------------------------------
484        ELSE IF (i>62.AND.i<=75) THEN
485          ierr = nf90_get_var(nid, var3didin(i), resul2)
486          print *, 'read_SCM(resul2), on a lu ', i, name_var(i)
487          IF(ierr/=nf90_noerr) THEN
488            print *, 'E Pb a la lecture de cas.nc: ', name_var(i)
489            stop "getvarup"
490          endif
491          PRINT*, 'Lecture de la variable #i  ', i, name_var(i), minval(resul2), maxval(resul2)
492
493          !-----------------------------------------------------------------------
494          ! Reading scalar variables (lat,lon)
495          !-----------------------------------------------------------------------
496        else
497          ierr = nf90_get_var(nid, var3didin(i), resul3)
498          print *, 'read_SCM(resul3), on a lu ', i, name_var(i)
499          IF(ierr/=nf90_noerr) THEN
500            print *, 'F Pb a la lecture de cas.nc: ', name_var(i)
501            stop "getvarup"
502          endif
503          PRINT*, 'Lecture de la variable #i ', i, name_var(i), resul3
504        endif
505      endif
506
507      !-----------------------------------------------------------------------
508      ! Attributing variables
509      !-----------------------------------------------------------------------
510      select case(i)
511        !case(1) ; ap=apbp       ! donnees indexees en nlevel+1
512        ! case(2) ; bp=apbp
513      case(3) ; zzh = apbp
514      case(4) ; pph = apbp
515      case(5) ; temp0 = resul1    ! donnees initiales
516      case(6) ; qv0 = resul1
517      case(7) ; ql0 = resul1
518      case(8) ; qi0 = resul1
519      case(9) ; u0 = resul1
520      case(10) ; v0 = resul1
521      case(11) ; tke0 = resul1
522      case(12) ; pp0 = resul1
523      case(13) ; vitw = resul    ! donnees indexees en nlevel,time
524      case(14) ; omega = resul
525      case(15) ; ug = resul
526      case(16) ; vg = resul
527      case(17) ; du = resul
528      case(18) ; hu = resul
529      case(19) ; vu = resul
530      case(20) ; dv = resul
531      case(21) ; hv = resul
532      case(22) ; vv = resul
533      case(23) ; dt = resul
534      case(24) ; ht = resul
535      case(25) ; vt = resul
536      case(26) ; dq = resul
537      case(27) ; hq = resul
538      case(28) ; vq = resul
539      case(29) ; dth = resul
540      case(30) ; hth = resul
541      case(31) ; vth = resul
542      case(32) ; hthl = resul
543      case(33) ; dr = resul
544      case(34) ; hr = resul
545      case(35) ; vr = resul
546      case(36) ; dtrad = resul
547      case(37) ; q1 = resul
548      case(38) ; q2 = resul
549      case(39) ; uw = resul
550      case(40) ; vw = resul
551      case(41) ; rh = resul
552      case(42) ; temp_nudg = resul
553      case(43) ; qv_nudg = resul
554      case(44) ; u_nudg = resul
555      case(45) ; v_nudg = resul
556      case(46) ; zz = resul      ! donnees en time,nlevel pour profil initial
557      case(47) ; pp = resul
558      case(48) ; temp = resul
559      case(49) ; theta = resul
560      case(50) ; thv = resul
561      case(51) ; thl = resul
562      case(52) ; qv = resul
563      case(53) ; ql = resul
564      case(54) ; qi = resul
565      case(55) ; rv = resul
566      case(56) ; u = resul
567      case(57) ; v = resul
568      case(58) ; invtau_temp_nudg = resul
569      case(59) ; invtau_qv_nudg = resul
570      case(60) ; invtau_u_nudg = resul
571      case(61) ; invtau_v_nudg = resul
572      case(62) ; tkes = resul2   ! donnees indexees en time
573      case(63) ; sens = resul2
574      case(64) ; flat = resul2
575      case(65) ; ts = resul2
576      case(66) ; tskin = resul2
577      case(67) ; ps = resul2
578      case(68) ; ustar = resul2
579      case(69) ; orog_cas = resul3      ! constantes
580      case(70) ; albedo_cas = resul3
581      case(71) ; emiss_cas = resul3
582      case(72) ; q_skin_cas = resul3
583      case(73) ; mom_rough = resul3
584      case(74) ; heat_rough = resul3
585      case(75) ; o3_cas = resul3
586      case(76) ; rugos_cas = resul3
587      case(77) ; clay_cas = resul3
588      case(78) ; sand_cas = resul3
589      end select
590      resul = 0.
591      resul1 = 0.
592      resul2 = 0.
593      resul3 = 0.
594    enddo
595    PRINT*, 'Lecture de la variable APRES ,sens ', minval(sens), maxval(sens)
596    PRINT*, 'Lecture de la variable APRES ,flat ', minval(flat), maxval(flat)
597
598    !CR:ATTENTION EN ATTENTE DE REGLER LA QUESTION DU PAS DE TEMPS INITIAL
599    DO t = 1, ntime
600      DO k = 1, nlevel
601        temp(k, t) = temp0(k)
602        qv(k, t) = qv0(k)
603        ql(k, t) = ql0(k)
604        qi(k, t) = qi0(k)
605        u(k, t) = u0(k)
606        v(k, t) = v0(k)
607        tke(k, t) = tke0(k)
608      enddo
609    enddo
610    !!!! TRAVAIL : EN FONCTION DES DECISIONS SUR LA SPECIFICATION DE W
611    !!!omega=-vitw*pres*rg/(rd*temp)
612    !-----------------------------------------------------------------------
613
614  END SUBROUTINE read_SCM
615  !======================================================================
616
617  !======================================================================
618
619  !**********************************************************************************************
620
621  !**********************************************************************************************
622  SUBROUTINE interp_case_time_std(day, day1, annee_ref                           &
623          !    &         ,year_cas,day_cas,nt_cas,pdt_forc,nlev_cas                         &
624          , nt_cas, nlev_cas                                                   &
625          , ts_cas, tskin_cas, ps_cas, plev_cas, t_cas, theta_cas, thv_cas, thl_cas            &
626          , qv_cas, ql_cas, qi_cas, u_cas, v_cas                                  &
627          , ug_cas, vg_cas, temp_nudg_cas, qv_nudg_cas, u_nudg_cas, v_nudg_cas     &
628          , invtau_temp_nudg_cas, invtau_qv_nudg_cas, invtau_u_nudg_cas, invtau_v_nudg_cas     &
629          , vitw_cas, omega_cas, tke_cas, du_cas, hu_cas, vu_cas             &
630          , dv_cas, hv_cas, vv_cas, dt_cas, ht_cas, vt_cas, dtrad_cas               &
631          , dq_cas, hq_cas, vq_cas, dth_cas, hth_cas, vth_cas                      &
632          , lat_cas, sens_cas, ustar_cas                                        &
633          , uw_cas, vw_cas, q1_cas, q2_cas, tkes_cas                               &
634
635          , ts_prof_cas, tskin_prof_cas, ps_prof_cas, plev_prof_cas, t_prof_cas, theta_prof_cas               &
636          , thv_prof_cas, thl_prof_cas, qv_prof_cas, ql_prof_cas, qi_prof_cas     &
637          , u_prof_cas, v_prof_cas, ug_prof_cas, vg_prof_cas                     &
638          , temp_nudg_prof_cas, qv_nudg_prof_cas, u_nudg_prof_cas, v_nudg_prof_cas     &
639          , invtau_temp_nudg_prof_cas, invtau_qv_nudg_prof_cas, invtau_u_nudg_prof_cas, invtau_v_nudg_prof_cas     &
640          , vitw_prof_cas, omega_prof_cas, tke_prof_cas, du_prof_cas, hu_prof_cas, vu_prof_cas  &
641          , dv_prof_cas, hv_prof_cas, vv_prof_cas, dt_prof_cas                   &
642          , ht_prof_cas, vt_prof_cas, dtrad_prof_cas, dq_prof_cas                &
643          , hq_prof_cas, vq_prof_cas, dth_prof_cas, hth_prof_cas, vth_prof_cas    &
644          , lat_prof_cas, sens_prof_cas                                        &
645          , ustar_prof_cas, uw_prof_cas, vw_prof_cas, q1_prof_cas, q2_prof_cas, tkes_prof_cas)
646
647    USE lmdz_compar1d
648    USE lmdz_date_cas, ONLY: year_ini_cas, mth_ini_cas, day_deb, heure_ini_cas, pdt_cas, day_ju_ini_cas
649
650    IMPLICIT NONE
651
652    !---------------------------------------------------------------------------------------
653    ! Time interpolation of a 2D field to the timestep corresponding to day
654
655    ! day: current julian day (e.g. 717538.2)
656    ! day1: first day of the simulation
657    ! nt_cas: total nb of data in the forcing
658    ! pdt_cas: total time interval (in sec) between 2 forcing data
659    !---------------------------------------------------------------------------------------
660
661    ! inputs:
662    INTEGER annee_ref
663    INTEGER nt_cas, nlev_cas
664    REAL day, day1, day_cas
665    REAL ts_cas(nt_cas), tskin_cas(nt_cas), ps_cas(nt_cas)
666    REAL plev_cas(nlev_cas, nt_cas)
667    REAL t_cas(nlev_cas, nt_cas), theta_cas(nlev_cas, nt_cas)
668    REAL thv_cas(nlev_cas, nt_cas), thl_cas(nlev_cas, nt_cas)
669    REAL qv_cas(nlev_cas, nt_cas), ql_cas(nlev_cas, nt_cas), qi_cas(nlev_cas, nt_cas)
670    REAL u_cas(nlev_cas, nt_cas), v_cas(nlev_cas, nt_cas)
671    REAL ug_cas(nlev_cas, nt_cas), vg_cas(nlev_cas, nt_cas)
672    REAL temp_nudg_cas(nlev_cas, nt_cas), qv_nudg_cas(nlev_cas, nt_cas)
673    REAL u_nudg_cas(nlev_cas, nt_cas), v_nudg_cas(nlev_cas, nt_cas)
674
675    REAL invtau_temp_nudg_cas(nlev_cas, nt_cas), invtau_qv_nudg_cas(nlev_cas, nt_cas)
676    REAL invtau_u_nudg_cas(nlev_cas, nt_cas), invtau_v_nudg_cas(nlev_cas, nt_cas)
677
678    REAL vitw_cas(nlev_cas, nt_cas), omega_cas(nlev_cas, nt_cas), tke_cas(nlev_cas, nt_cas)
679    REAL du_cas(nlev_cas, nt_cas), hu_cas(nlev_cas, nt_cas), vu_cas(nlev_cas, nt_cas)
680    REAL dv_cas(nlev_cas, nt_cas), hv_cas(nlev_cas, nt_cas), vv_cas(nlev_cas, nt_cas)
681    REAL dt_cas(nlev_cas, nt_cas), ht_cas(nlev_cas, nt_cas), vt_cas(nlev_cas, nt_cas)
682    REAL dth_cas(nlev_cas, nt_cas), hth_cas(nlev_cas, nt_cas), vth_cas(nlev_cas, nt_cas)
683    REAL dtrad_cas(nlev_cas, nt_cas)
684    REAL dq_cas(nlev_cas, nt_cas), hq_cas(nlev_cas, nt_cas), vq_cas(nlev_cas, nt_cas)
685    REAL lat_cas(nt_cas), sens_cas(nt_cas), tkes_cas(nt_cas)
686    REAL ustar_cas(nt_cas), uw_cas(nlev_cas, nt_cas), vw_cas(nlev_cas, nt_cas)
687    REAL q1_cas(nlev_cas, nt_cas), q2_cas(nlev_cas, nt_cas)
688
689    ! outputs:
690    REAL plev_prof_cas(nlev_cas)
691    REAL t_prof_cas(nlev_cas), theta_prof_cas(nlev_cas), thl_prof_cas(nlev_cas), thv_prof_cas(nlev_cas)
692    REAL qv_prof_cas(nlev_cas), ql_prof_cas(nlev_cas), qi_prof_cas(nlev_cas)
693    REAL u_prof_cas(nlev_cas), v_prof_cas(nlev_cas)
694    REAL ug_prof_cas(nlev_cas), vg_prof_cas(nlev_cas)
695    REAL temp_nudg_prof_cas(nlev_cas), qv_nudg_prof_cas(nlev_cas)
696    REAL u_nudg_prof_cas(nlev_cas), v_nudg_prof_cas(nlev_cas)
697
698    REAL invtau_temp_nudg_prof_cas(nlev_cas), invtau_qv_nudg_prof_cas(nlev_cas)
699    REAL invtau_u_nudg_prof_cas(nlev_cas), invtau_v_nudg_prof_cas(nlev_cas)
700
701    REAL vitw_prof_cas(nlev_cas), omega_prof_cas(nlev_cas), tke_prof_cas(nlev_cas)
702    REAL du_prof_cas(nlev_cas), hu_prof_cas(nlev_cas), vu_prof_cas(nlev_cas)
703    REAL dv_prof_cas(nlev_cas), hv_prof_cas(nlev_cas), vv_prof_cas(nlev_cas)
704    REAL dt_prof_cas(nlev_cas), ht_prof_cas(nlev_cas), vt_prof_cas(nlev_cas)
705    REAL dth_prof_cas(nlev_cas), hth_prof_cas(nlev_cas), vth_prof_cas(nlev_cas)
706    REAL dtrad_prof_cas(nlev_cas)
707    REAL dq_prof_cas(nlev_cas), hq_prof_cas(nlev_cas), vq_prof_cas(nlev_cas)
708    REAL lat_prof_cas, sens_prof_cas, tkes_prof_cas, ts_prof_cas, tskin_prof_cas, ps_prof_cas, ustar_prof_cas
709    REAL uw_prof_cas(nlev_cas), vw_prof_cas(nlev_cas), q1_prof_cas(nlev_cas), q2_prof_cas(nlev_cas)
710    ! local:
711    INTEGER it_cas1, it_cas2, k
712    REAL timeit, time_cas1, time_cas2, frac
713
714    PRINT*, 'Check time', day1, day_ju_ini_cas, day_deb + 1, pdt_cas
715    !       do k=1,nlev_cas
716    !       PRINT*,'debut de interp2_case_time, plev_cas=',k,plev_cas(k,1)
717    !       enddo
718
719    ! On teste si la date du cas AMMA est correcte.
720    ! C est pour memoire car en fait les fichiers .def
721    ! sont censes etre corrects.
722    ! A supprimer a terme (MPL 20150623)
723    !     if ((forcing_type.EQ.10).AND.(1.EQ.0)) THEN
724    ! Check that initial day of the simulation consistent with AMMA case:
725    !      if (annee_ref.NE.2006) THEN
726    !       PRINT*,'Pour AMMA, annee_ref doit etre 2006'
727    !       PRINT*,'Changer annee_ref dans run.def'
728    !       stop
729    !      endif
730    !      if (annee_ref.EQ.2006 .AND. day1.lt.day_cas) THEN
731    !       PRINT*,'AMMA a debute le 10 juillet 2006',day1,day_cas
732    !       PRINT*,'Changer dayref dans run.def'
733    !       stop
734    !      endif
735    !      if (annee_ref.EQ.2006 .AND. day1.gt.day_cas+1) THEN
736    !       PRINT*,'AMMA a fini le 11 juillet'
737    !       PRINT*,'Changer dayref ou nday dans run.def'
738    !       stop
739    !      endif
740    !      endif
741
742    ! Determine timestep relative to the 1st day:
743    !       timeit=(day-day1)*86400.
744    !       if (annee_ref.EQ.1992) THEN
745    !        timeit=(day-day_cas)*86400.
746    !       else
747    !        timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
748    !       endif
749    timeit = (day - day_ju_ini_cas) * 86400
750    print *, 'day=', day
751    print *, 'day_ju_ini_cas=', day_ju_ini_cas
752    print *, 'pdt_cas=', pdt_cas
753    print *, 'timeit=', timeit
754    print *, 'nt_cas=', nt_cas
755
756    ! Determine the closest observation times:
757    !       it_cas1=INT(timeit/pdt_cas)+1
758    !       it_cas2=it_cas1 + 1
759    !       time_cas1=(it_cas1-1)*pdt_cas
760    !       time_cas2=(it_cas2-1)*pdt_cas
761
762    it_cas1 = INT(timeit / pdt_cas) + 1
763    IF (it_cas1 == nt_cas) THEN
764      it_cas2 = it_cas1
765    ELSE
766      it_cas2 = it_cas1 + 1
767    ENDIF
768    time_cas1 = (it_cas1 - 1) * pdt_cas
769    time_cas2 = (it_cas2 - 1) * pdt_cas
770    !     print *,'timeit,pdt_cas,nt_cas=',timeit,pdt_cas,nt_cas
771    !     print *,'it_cas1,it_cas2,time_cas1,time_cas2=',it_cas1,it_cas2,time_cas1,time_cas2
772
773    IF (it_cas1 > nt_cas) THEN
774      WRITE(*, *) 'PB-stop: day, day_ju_ini_cas,it_cas1, it_cas2, timeit: '            &
775              , day, day_ju_ini_cas, it_cas1, it_cas2, timeit
776      stop
777    endif
778
779    ! time interpolation:
780    IF (it_cas1 == it_cas2) THEN
781      frac = 0.
782    ELSE
783      frac = (time_cas2 - timeit) / (time_cas2 - time_cas1)
784      frac = max(frac, 0.0)
785    ENDIF
786
787    lat_prof_cas = lat_cas(it_cas2)                                   &
788            - frac * (lat_cas(it_cas2) - lat_cas(it_cas1))
789    sens_prof_cas = sens_cas(it_cas2)                                 &
790            - frac * (sens_cas(it_cas2) - sens_cas(it_cas1))
791    tkes_prof_cas = tkes_cas(it_cas2)                                   &
792            - frac * (tkes_cas(it_cas2) - tkes_cas(it_cas1))
793    ts_prof_cas = ts_cas(it_cas2)                                     &
794            - frac * (ts_cas(it_cas2) - ts_cas(it_cas1))
795    tskin_prof_cas = tskin_cas(it_cas2)                                     &
796            - frac * (tskin_cas(it_cas2) - tskin_cas(it_cas1))
797    ps_prof_cas = ps_cas(it_cas2)                                     &
798            - frac * (ps_cas(it_cas2) - ps_cas(it_cas1))
799    ustar_prof_cas = ustar_cas(it_cas2)                               &
800            - frac * (ustar_cas(it_cas2) - ustar_cas(it_cas1))
801
802    DO k = 1, nlev_cas
803      plev_prof_cas(k) = plev_cas(k, it_cas2)                           &
804              - frac * (plev_cas(k, it_cas2) - plev_cas(k, it_cas1))
805      t_prof_cas(k) = t_cas(k, it_cas2)                                 &
806              - frac * (t_cas(k, it_cas2) - t_cas(k, it_cas1))
807      !print *,'k,frac,plev_cas1,plev_cas2=',k,frac,plev_cas(k,it_cas1),plev_cas(k,it_cas2)
808      theta_prof_cas(k) = theta_cas(k, it_cas2)                         &
809              - frac * (theta_cas(k, it_cas2) - theta_cas(k, it_cas1))
810      thv_prof_cas(k) = thv_cas(k, it_cas2)                             &
811              - frac * (thv_cas(k, it_cas2) - thv_cas(k, it_cas1))
812      thl_prof_cas(k) = thl_cas(k, it_cas2)                             &
813              - frac * (thl_cas(k, it_cas2) - thl_cas(k, it_cas1))
814      qv_prof_cas(k) = qv_cas(k, it_cas2)                               &
815              - frac * (qv_cas(k, it_cas2) - qv_cas(k, it_cas1))
816      ql_prof_cas(k) = ql_cas(k, it_cas2)                               &
817              - frac * (ql_cas(k, it_cas2) - ql_cas(k, it_cas1))
818      qi_prof_cas(k) = qi_cas(k, it_cas2)                               &
819              - frac * (qi_cas(k, it_cas2) - qi_cas(k, it_cas1))
820      u_prof_cas(k) = u_cas(k, it_cas2)                                 &
821              - frac * (u_cas(k, it_cas2) - u_cas(k, it_cas1))
822      v_prof_cas(k) = v_cas(k, it_cas2)                                 &
823              - frac * (v_cas(k, it_cas2) - v_cas(k, it_cas1))
824      ug_prof_cas(k) = ug_cas(k, it_cas2)                               &
825              - frac * (ug_cas(k, it_cas2) - ug_cas(k, it_cas1))
826      vg_prof_cas(k) = vg_cas(k, it_cas2)                               &
827              - frac * (vg_cas(k, it_cas2) - vg_cas(k, it_cas1))
828      temp_nudg_prof_cas(k) = temp_nudg_cas(k, it_cas2)                    &
829              - frac * (temp_nudg_cas(k, it_cas2) - temp_nudg_cas(k, it_cas1))
830      qv_nudg_prof_cas(k) = qv_nudg_cas(k, it_cas2)                        &
831              - frac * (qv_nudg_cas(k, it_cas2) - qv_nudg_cas(k, it_cas1))
832      u_nudg_prof_cas(k) = u_nudg_cas(k, it_cas2)                          &
833              - frac * (u_nudg_cas(k, it_cas2) - u_nudg_cas(k, it_cas1))
834      v_nudg_prof_cas(k) = v_nudg_cas(k, it_cas2)                          &
835              - frac * (v_nudg_cas(k, it_cas2) - v_nudg_cas(k, it_cas1))
836      invtau_temp_nudg_prof_cas(k) = invtau_temp_nudg_cas(k, it_cas2)                    &
837              - frac * (invtau_temp_nudg_cas(k, it_cas2) - invtau_temp_nudg_cas(k, it_cas1))
838      invtau_qv_nudg_prof_cas(k) = invtau_qv_nudg_cas(k, it_cas2)                        &
839              - frac * (invtau_qv_nudg_cas(k, it_cas2) - invtau_qv_nudg_cas(k, it_cas1))
840      invtau_u_nudg_prof_cas(k) = invtau_u_nudg_cas(k, it_cas2)                          &
841              - frac * (invtau_u_nudg_cas(k, it_cas2) - invtau_u_nudg_cas(k, it_cas1))
842      invtau_v_nudg_prof_cas(k) = invtau_v_nudg_cas(k, it_cas2)                          &
843              - frac * (invtau_v_nudg_cas(k, it_cas2) - invtau_v_nudg_cas(k, it_cas1))
844      vitw_prof_cas(k) = vitw_cas(k, it_cas2)                           &
845              - frac * (vitw_cas(k, it_cas2) - vitw_cas(k, it_cas1))
846      omega_prof_cas(k) = omega_cas(k, it_cas2)                         &
847              - frac * (omega_cas(k, it_cas2) - omega_cas(k, it_cas1))
848      tke_prof_cas(k) = tke_cas(k, it_cas2)                         &
849              - frac * (tke_cas(k, it_cas2) - tke_cas(k, it_cas1))
850      du_prof_cas(k) = du_cas(k, it_cas2)                               &
851              - frac * (du_cas(k, it_cas2) - du_cas(k, it_cas1))
852      hu_prof_cas(k) = hu_cas(k, it_cas2)                               &
853              - frac * (hu_cas(k, it_cas2) - hu_cas(k, it_cas1))
854      vu_prof_cas(k) = vu_cas(k, it_cas2)                               &
855              - frac * (vu_cas(k, it_cas2) - vu_cas(k, it_cas1))
856      dv_prof_cas(k) = dv_cas(k, it_cas2)                               &
857              - frac * (dv_cas(k, it_cas2) - dv_cas(k, it_cas1))
858      hv_prof_cas(k) = hv_cas(k, it_cas2)                               &
859              - frac * (hv_cas(k, it_cas2) - hv_cas(k, it_cas1))
860      vv_prof_cas(k) = vv_cas(k, it_cas2)                               &
861              - frac * (vv_cas(k, it_cas2) - vv_cas(k, it_cas1))
862      dt_prof_cas(k) = dt_cas(k, it_cas2)                               &
863              - frac * (dt_cas(k, it_cas2) - dt_cas(k, it_cas1))
864      ht_prof_cas(k) = ht_cas(k, it_cas2)                               &
865              - frac * (ht_cas(k, it_cas2) - ht_cas(k, it_cas1))
866      vt_prof_cas(k) = vt_cas(k, it_cas2)                               &
867              - frac * (vt_cas(k, it_cas2) - vt_cas(k, it_cas1))
868      dth_prof_cas(k) = dth_cas(k, it_cas2)                             &
869              - frac * (dth_cas(k, it_cas2) - dth_cas(k, it_cas1))
870      hth_prof_cas(k) = hth_cas(k, it_cas2)                             &
871              - frac * (hth_cas(k, it_cas2) - hth_cas(k, it_cas1))
872      vth_prof_cas(k) = vth_cas(k, it_cas2)                             &
873              - frac * (vth_cas(k, it_cas2) - vth_cas(k, it_cas1))
874      dtrad_prof_cas(k) = dtrad_cas(k, it_cas2)                         &
875              - frac * (dtrad_cas(k, it_cas2) - dtrad_cas(k, it_cas1))
876      dq_prof_cas(k) = dq_cas(k, it_cas2)                               &
877              - frac * (dq_cas(k, it_cas2) - dq_cas(k, it_cas1))
878      hq_prof_cas(k) = hq_cas(k, it_cas2)                               &
879              - frac * (hq_cas(k, it_cas2) - hq_cas(k, it_cas1))
880      vq_prof_cas(k) = vq_cas(k, it_cas2)                               &
881              - frac * (vq_cas(k, it_cas2) - vq_cas(k, it_cas1))
882      uw_prof_cas(k) = uw_cas(k, it_cas2)                                &
883              - frac * (uw_cas(k, it_cas2) - uw_cas(k, it_cas1))
884      vw_prof_cas(k) = vw_cas(k, it_cas2)                                &
885              - frac * (vw_cas(k, it_cas2) - vw_cas(k, it_cas1))
886      q1_prof_cas(k) = q1_cas(k, it_cas2)                                &
887              - frac * (q1_cas(k, it_cas2) - q1_cas(k, it_cas1))
888      q2_prof_cas(k) = q2_cas(k, it_cas2)                                &
889              - frac * (q2_cas(k, it_cas2) - q2_cas(k, it_cas1))
890    enddo
891
892  END SUBROUTINE interp_case_time_std
893
894  !**********************************************************************************************
895  !=====================================================================
896  SUBROUTINE interp2_case_vertical_std(play, plev, nlev_cas, plev_prof_cas                           &
897          , t_prof_cas, th_prof_cas, thv_prof_cas, thl_prof_cas                                       &
898          , qv_prof_cas, ql_prof_cas, qi_prof_cas, u_prof_cas, v_prof_cas                              &
899          , ug_prof_cas, vg_prof_cas                                                                &
900          , temp_nudg_prof_cas, qv_nudg_prof_cas, u_nudg_prof_cas, v_nudg_prof_cas                    &
901          , invtau_temp_nudg_prof_cas, invtau_qv_nudg_prof_cas, invtau_u_nudg_prof_cas, invtau_v_nudg_prof_cas &
902          , vitw_prof_cas, omega_prof_cas, tke_prof_cas                                              &
903          , du_prof_cas, hu_prof_cas, vu_prof_cas, dv_prof_cas, hv_prof_cas, vv_prof_cas                &
904          , dt_prof_cas, ht_prof_cas, vt_prof_cas, dtrad_prof_cas, dq_prof_cas, hq_prof_cas, vq_prof_cas &
905          , dth_prof_cas, hth_prof_cas, vth_prof_cas                                                 &
906
907          , t_mod_cas, theta_mod_cas, thv_mod_cas, thl_mod_cas                                        &
908          , qv_mod_cas, ql_mod_cas, qi_mod_cas, u_mod_cas, v_mod_cas                                   &
909          , ug_mod_cas, vg_mod_cas                                                                  &
910          , temp_nudg_mod_cas, qv_nudg_mod_cas, u_nudg_mod_cas, v_nudg_mod_cas                        &
911          , invtau_temp_nudg_mod_cas, invtau_qv_nudg_mod_cas, invtau_u_nudg_mod_cas, invtau_v_nudg_mod_cas                        &
912          , w_mod_cas, omega_mod_cas, tke_mod_cas                                                    &
913          , du_mod_cas, hu_mod_cas, vu_mod_cas, dv_mod_cas, hv_mod_cas, vv_mod_cas                      &
914          , dt_mod_cas, ht_mod_cas, vt_mod_cas, dtrad_mod_cas, dq_mod_cas, hq_mod_cas, vq_mod_cas        &
915          , dth_mod_cas, hth_mod_cas, vth_mod_cas, mxcalc)
916
917    USE lmdz_yomcst
918
919    IMPLICIT NONE
920
921    INCLUDE "dimensions.h"
922
923    !-------------------------------------------------------------------------
924    ! Vertical interpolation of generic case forcing data onto mod_casel levels
925    !-------------------------------------------------------------------------
926
927    INTEGER nlevmax
928    parameter (nlevmax = 41)
929    INTEGER nlev_cas, mxcalc
930    !       real play(llm), plev_prof(nlevmax)
931    !       real t_prof(nlevmax),q_prof(nlevmax)
932    !       real u_prof(nlevmax),v_prof(nlevmax), w_prof(nlevmax)
933    !       real ht_prof(nlevmax),vt_prof(nlevmax)
934    !       real hq_prof(nlevmax),vq_prof(nlevmax)
935
936    REAL play(llm), plev(llm + 1), plev_prof_cas(nlev_cas)
937    REAL t_prof_cas(nlev_cas), th_prof_cas(nlev_cas), thv_prof_cas(nlev_cas), thl_prof_cas(nlev_cas)
938    REAL qv_prof_cas(nlev_cas), ql_prof_cas(nlev_cas), qi_prof_cas(nlev_cas)
939    REAL u_prof_cas(nlev_cas), v_prof_cas(nlev_cas)
940    REAL ug_prof_cas(nlev_cas), vg_prof_cas(nlev_cas), vitw_prof_cas(nlev_cas), omega_prof_cas(nlev_cas), tke_prof_cas(nlev_cas)
941    REAL temp_nudg_prof_cas(nlev_cas), qv_nudg_prof_cas(nlev_cas)
942    REAL u_nudg_prof_cas(nlev_cas), v_nudg_prof_cas(nlev_cas)
943    REAL invtau_temp_nudg_prof_cas(nlev_cas), invtau_qv_nudg_prof_cas(nlev_cas)
944    REAL invtau_u_nudg_prof_cas(nlev_cas), invtau_v_nudg_prof_cas(nlev_cas)
945
946    REAL du_prof_cas(nlev_cas), hu_prof_cas(nlev_cas), vu_prof_cas(nlev_cas)
947    REAL dv_prof_cas(nlev_cas), hv_prof_cas(nlev_cas), vv_prof_cas(nlev_cas)
948    REAL dt_prof_cas(nlev_cas), ht_prof_cas(nlev_cas), vt_prof_cas(nlev_cas), dtrad_prof_cas(nlev_cas)
949    REAL dth_prof_cas(nlev_cas), hth_prof_cas(nlev_cas), vth_prof_cas(nlev_cas)
950    REAL dq_prof_cas(nlev_cas), hq_prof_cas(nlev_cas), vq_prof_cas(nlev_cas)
951
952    REAL t_mod_cas(llm), theta_mod_cas(llm), thv_mod_cas(llm), thl_mod_cas(llm)
953    REAL qv_mod_cas(llm), ql_mod_cas(llm), qi_mod_cas(llm)
954    REAL u_mod_cas(llm), v_mod_cas(llm)
955    REAL ug_mod_cas(llm), vg_mod_cas(llm), w_mod_cas(llm), omega_mod_cas(llm), tke_mod_cas(llm + 1)
956    REAL temp_nudg_mod_cas(llm), qv_nudg_mod_cas(llm)
957    REAL u_nudg_mod_cas(llm), v_nudg_mod_cas(llm)
958    REAL invtau_temp_nudg_mod_cas(llm), invtau_qv_nudg_mod_cas(llm)
959    REAL invtau_u_nudg_mod_cas(llm), invtau_v_nudg_mod_cas(llm)
960    REAL du_mod_cas(llm), hu_mod_cas(llm), vu_mod_cas(llm)
961    REAL dv_mod_cas(llm), hv_mod_cas(llm), vv_mod_cas(llm)
962    REAL dt_mod_cas(llm), ht_mod_cas(llm), vt_mod_cas(llm), dtrad_mod_cas(llm)
963    REAL dth_mod_cas(llm), hth_mod_cas(llm), vth_mod_cas(llm)
964    REAL dq_mod_cas(llm), hq_mod_cas(llm), vq_mod_cas(llm)
965
966    INTEGER l, k, k1, k2
967    REAL frac, frac1, frac2, fact
968
969
970
971    ! for variables defined at the middle of layers
972
973    DO l = 1, llm
974
975      IF (play(l)>=plev_prof_cas(nlev_cas)) THEN
976        mxcalc = l
977        !        print *,'debut interp2, mxcalc=',mxcalc
978        k1 = 0
979        k2 = 0
980
981        IF (play(l)<=plev_prof_cas(1)) THEN
982          DO k = 1, nlev_cas - 1
983            IF (play(l)<=plev_prof_cas(k).AND. play(l)>plev_prof_cas(k + 1)) THEN
984              k1 = k
985              k2 = k + 1
986            endif
987          enddo
988
989          IF (k1==0 .OR. k2==0) THEN
990            WRITE(*, *) 'PB! k1, k2 = ', k1, k2
991            WRITE(*, *) 'l,play(l) = ', l, play(l) / 100
992            DO k = 1, nlev_cas - 1
993              WRITE(*, *) 'k,plev_prof_cas(k) = ', k, plev_prof_cas(k) / 100
994            enddo
995          endif
996
997          frac = (plev_prof_cas(k2) - play(l)) / (plev_prof_cas(k2) - plev_prof_cas(k1))
998
999          t_mod_cas(l) = t_prof_cas(k2) - frac * (t_prof_cas(k2) - t_prof_cas(k1))
1000          theta_mod_cas(l) = th_prof_cas(k2) - frac * (th_prof_cas(k2) - th_prof_cas(k1))
1001          IF(theta_mod_cas(l)/=0) t_mod_cas(l) = theta_mod_cas(l) * (play(l) / 100000.)**(RD / RCPD)
1002          thv_mod_cas(l) = thv_prof_cas(k2) - frac * (thv_prof_cas(k2) - thv_prof_cas(k1))
1003          thl_mod_cas(l) = thl_prof_cas(k2) - frac * (thl_prof_cas(k2) - thl_prof_cas(k1))
1004          qv_mod_cas(l) = qv_prof_cas(k2) - frac * (qv_prof_cas(k2) - qv_prof_cas(k1))
1005          ql_mod_cas(l) = ql_prof_cas(k2) - frac * (ql_prof_cas(k2) - ql_prof_cas(k1))
1006          qi_mod_cas(l) = qi_prof_cas(k2) - frac * (qi_prof_cas(k2) - qi_prof_cas(k1))
1007          u_mod_cas(l) = u_prof_cas(k2) - frac * (u_prof_cas(k2) - u_prof_cas(k1))
1008          v_mod_cas(l) = v_prof_cas(k2) - frac * (v_prof_cas(k2) - v_prof_cas(k1))
1009          ug_mod_cas(l) = ug_prof_cas(k2) - frac * (ug_prof_cas(k2) - ug_prof_cas(k1))
1010          vg_mod_cas(l) = vg_prof_cas(k2) - frac * (vg_prof_cas(k2) - vg_prof_cas(k1))
1011          temp_nudg_mod_cas(l) = temp_nudg_prof_cas(k2) - frac * (temp_nudg_prof_cas(k2) - temp_nudg_prof_cas(k1))
1012          qv_nudg_mod_cas(l) = qv_nudg_prof_cas(k2) - frac * (qv_nudg_prof_cas(k2) - qv_nudg_prof_cas(k1))
1013          u_nudg_mod_cas(l) = u_nudg_prof_cas(k2) - frac * (u_nudg_prof_cas(k2) - u_nudg_prof_cas(k1))
1014          v_nudg_mod_cas(l) = v_nudg_prof_cas(k2) - frac * (v_nudg_prof_cas(k2) - v_nudg_prof_cas(k1))
1015
1016          invtau_temp_nudg_mod_cas(l) = invtau_temp_nudg_prof_cas(k2) &
1017                  - frac * (invtau_temp_nudg_prof_cas(k2) - invtau_temp_nudg_prof_cas(k1))
1018          invtau_qv_nudg_mod_cas(l) = invtau_qv_nudg_prof_cas(k2) - frac * (invtau_qv_nudg_prof_cas(k2) - invtau_qv_nudg_prof_cas(k1))
1019          invtau_u_nudg_mod_cas(l) = invtau_u_nudg_prof_cas(k2) - frac * (invtau_u_nudg_prof_cas(k2) - invtau_u_nudg_prof_cas(k1))
1020          invtau_v_nudg_mod_cas(l) = invtau_v_nudg_prof_cas(k2) - frac * (invtau_v_nudg_prof_cas(k2) - invtau_v_nudg_prof_cas(k1))
1021
1022          w_mod_cas(l) = vitw_prof_cas(k2) - frac * (vitw_prof_cas(k2) - vitw_prof_cas(k1))
1023          omega_mod_cas(l) = omega_prof_cas(k2) - frac * (omega_prof_cas(k2) - omega_prof_cas(k1))
1024          du_mod_cas(l) = du_prof_cas(k2) - frac * (du_prof_cas(k2) - du_prof_cas(k1))
1025          hu_mod_cas(l) = hu_prof_cas(k2) - frac * (hu_prof_cas(k2) - hu_prof_cas(k1))
1026          vu_mod_cas(l) = vu_prof_cas(k2) - frac * (vu_prof_cas(k2) - vu_prof_cas(k1))
1027          dv_mod_cas(l) = dv_prof_cas(k2) - frac * (dv_prof_cas(k2) - dv_prof_cas(k1))
1028          hv_mod_cas(l) = hv_prof_cas(k2) - frac * (hv_prof_cas(k2) - hv_prof_cas(k1))
1029          vv_mod_cas(l) = vv_prof_cas(k2) - frac * (vv_prof_cas(k2) - vv_prof_cas(k1))
1030          dt_mod_cas(l) = dt_prof_cas(k2) - frac * (dt_prof_cas(k2) - dt_prof_cas(k1))
1031          ht_mod_cas(l) = ht_prof_cas(k2) - frac * (ht_prof_cas(k2) - ht_prof_cas(k1))
1032          vt_mod_cas(l) = vt_prof_cas(k2) - frac * (vt_prof_cas(k2) - vt_prof_cas(k1))
1033          dth_mod_cas(l) = dth_prof_cas(k2) - frac * (dth_prof_cas(k2) - dth_prof_cas(k1))
1034          hth_mod_cas(l) = hth_prof_cas(k2) - frac * (hth_prof_cas(k2) - hth_prof_cas(k1))
1035          vth_mod_cas(l) = vth_prof_cas(k2) - frac * (vth_prof_cas(k2) - vth_prof_cas(k1))
1036          dq_mod_cas(l) = dq_prof_cas(k2) - frac * (dq_prof_cas(k2) - dq_prof_cas(k1))
1037          hq_mod_cas(l) = hq_prof_cas(k2) - frac * (hq_prof_cas(k2) - hq_prof_cas(k1))
1038          vq_mod_cas(l) = vq_prof_cas(k2) - frac * (vq_prof_cas(k2) - vq_prof_cas(k1))
1039          dtrad_mod_cas(l) = dtrad_prof_cas(k2) - frac * (dtrad_prof_cas(k2) - dtrad_prof_cas(k1))
1040
1041        else !play>plev_prof_cas(1)
1042
1043          k1 = 1
1044          k2 = 2
1045          print *, 'interp2_vert, k1,k2=', plev_prof_cas(k1), plev_prof_cas(k2)
1046          frac1 = (play(l) - plev_prof_cas(k2)) / (plev_prof_cas(k1) - plev_prof_cas(k2))
1047          frac2 = (play(l) - plev_prof_cas(k1)) / (plev_prof_cas(k1) - plev_prof_cas(k2))
1048          t_mod_cas(l) = frac1 * t_prof_cas(k1) - frac2 * t_prof_cas(k2)
1049          theta_mod_cas(l) = frac1 * th_prof_cas(k1) - frac2 * th_prof_cas(k2)
1050          IF(theta_mod_cas(l)/=0) t_mod_cas(l) = theta_mod_cas(l) * (play(l) / 100000.)**(RD / RCPD)
1051          thv_mod_cas(l) = frac1 * thv_prof_cas(k1) - frac2 * thv_prof_cas(k2)
1052          thl_mod_cas(l) = frac1 * thl_prof_cas(k1) - frac2 * thl_prof_cas(k2)
1053          qv_mod_cas(l) = frac1 * qv_prof_cas(k1) - frac2 * qv_prof_cas(k2)
1054          ql_mod_cas(l) = frac1 * ql_prof_cas(k1) - frac2 * ql_prof_cas(k2)
1055          qi_mod_cas(l) = frac1 * qi_prof_cas(k1) - frac2 * qi_prof_cas(k2)
1056          u_mod_cas(l) = frac1 * u_prof_cas(k1) - frac2 * u_prof_cas(k2)
1057          v_mod_cas(l) = frac1 * v_prof_cas(k1) - frac2 * v_prof_cas(k2)
1058          ug_mod_cas(l) = frac1 * ug_prof_cas(k1) - frac2 * ug_prof_cas(k2)
1059          vg_mod_cas(l) = frac1 * vg_prof_cas(k1) - frac2 * vg_prof_cas(k2)
1060          temp_nudg_mod_cas(l) = frac1 * temp_nudg_prof_cas(k1) - frac2 * temp_nudg_prof_cas(k2)
1061          qv_nudg_mod_cas(l) = frac1 * qv_nudg_prof_cas(k1) - frac2 * qv_nudg_prof_cas(k2)
1062          u_nudg_mod_cas(l) = frac1 * u_nudg_prof_cas(k1) - frac2 * u_nudg_prof_cas(k2)
1063          v_nudg_mod_cas(l) = frac1 * v_nudg_prof_cas(k1) - frac2 * v_nudg_prof_cas(k2)
1064
1065          invtau_temp_nudg_mod_cas(l) = frac1 * invtau_temp_nudg_prof_cas(k1) - frac2 * invtau_temp_nudg_prof_cas(k2)
1066          invtau_qv_nudg_mod_cas(l) = frac1 * invtau_qv_nudg_prof_cas(k1) - frac2 * invtau_qv_nudg_prof_cas(k2)
1067          invtau_u_nudg_mod_cas(l) = frac1 * invtau_u_nudg_prof_cas(k1) - frac2 * invtau_u_nudg_prof_cas(k2)
1068          invtau_v_nudg_mod_cas(l) = frac1 * invtau_v_nudg_prof_cas(k1) - frac2 * invtau_v_nudg_prof_cas(k2)
1069
1070          w_mod_cas(l) = frac1 * vitw_prof_cas(k1) - frac2 * vitw_prof_cas(k2)
1071          omega_mod_cas(l) = frac1 * omega_prof_cas(k1) - frac2 * omega_prof_cas(k2)
1072          du_mod_cas(l) = frac1 * du_prof_cas(k1) - frac2 * du_prof_cas(k2)
1073          hu_mod_cas(l) = frac1 * hu_prof_cas(k1) - frac2 * hu_prof_cas(k2)
1074          vu_mod_cas(l) = frac1 * vu_prof_cas(k1) - frac2 * vu_prof_cas(k2)
1075          dv_mod_cas(l) = frac1 * dv_prof_cas(k1) - frac2 * dv_prof_cas(k2)
1076          hv_mod_cas(l) = frac1 * hv_prof_cas(k1) - frac2 * hv_prof_cas(k2)
1077          vv_mod_cas(l) = frac1 * vv_prof_cas(k1) - frac2 * vv_prof_cas(k2)
1078          dt_mod_cas(l) = frac1 * dt_prof_cas(k1) - frac2 * dt_prof_cas(k2)
1079          ht_mod_cas(l) = frac1 * ht_prof_cas(k1) - frac2 * ht_prof_cas(k2)
1080          vt_mod_cas(l) = frac1 * vt_prof_cas(k1) - frac2 * vt_prof_cas(k2)
1081          dth_mod_cas(l) = frac1 * dth_prof_cas(k1) - frac2 * dth_prof_cas(k2)
1082          hth_mod_cas(l) = frac1 * hth_prof_cas(k1) - frac2 * hth_prof_cas(k2)
1083          vth_mod_cas(l) = frac1 * vth_prof_cas(k1) - frac2 * vth_prof_cas(k2)
1084          dq_mod_cas(l) = frac1 * dq_prof_cas(k1) - frac2 * dq_prof_cas(k2)
1085          hq_mod_cas(l) = frac1 * hq_prof_cas(k1) - frac2 * hq_prof_cas(k2)
1086          vq_mod_cas(l) = frac1 * vq_prof_cas(k1) - frac2 * vq_prof_cas(k2)
1087          dtrad_mod_cas(l) = frac1 * dtrad_prof_cas(k1) - frac2 * dtrad_prof_cas(k2)
1088
1089        endif ! play.le.plev_prof_cas(1)
1090
1091      else ! above max altitude of forcing file
1092
1093        !jyg
1094        fact = 20. * (plev_prof_cas(nlev_cas) - play(l)) / plev_prof_cas(nlev_cas) !jyg
1095        fact = max(fact, 0.)                                           !jyg
1096        fact = exp(-fact)                                             !jyg
1097        t_mod_cas(l) = t_prof_cas(nlev_cas)                            !jyg
1098        theta_mod_cas(l) = th_prof_cas(nlev_cas)                       !jyg
1099        thv_mod_cas(l) = thv_prof_cas(nlev_cas)                        !jyg
1100        thl_mod_cas(l) = thl_prof_cas(nlev_cas)                        !jyg
1101        qv_mod_cas(l) = qv_prof_cas(nlev_cas) * fact                     !jyg
1102        ql_mod_cas(l) = ql_prof_cas(nlev_cas) * fact                     !jyg
1103        qi_mod_cas(l) = qi_prof_cas(nlev_cas) * fact                     !jyg
1104        u_mod_cas(l) = u_prof_cas(nlev_cas) * fact                       !jyg
1105        v_mod_cas(l) = v_prof_cas(nlev_cas) * fact                       !jyg
1106        ug_mod_cas(l) = ug_prof_cas(nlev_cas)                          !jyg
1107        vg_mod_cas(l) = vg_prof_cas(nlev_cas)                          !jyg
1108        temp_nudg_mod_cas(l) = temp_nudg_prof_cas(nlev_cas)            !jyg
1109        qv_nudg_mod_cas(l) = qv_nudg_prof_cas(nlev_cas)                !jyg
1110        u_nudg_mod_cas(l) = u_nudg_prof_cas(nlev_cas)                  !jyg
1111        v_nudg_mod_cas(l) = v_nudg_prof_cas(nlev_cas)                  !jyg
1112
1113        invtau_temp_nudg_mod_cas(l) = invtau_temp_nudg_prof_cas(nlev_cas)            !jyg
1114        invtau_qv_nudg_mod_cas(l) = invtau_qv_nudg_prof_cas(nlev_cas)                !jyg
1115        invtau_u_nudg_mod_cas(l) = invtau_u_nudg_prof_cas(nlev_cas)                  !jyg
1116        invtau_v_nudg_mod_cas(l) = invtau_v_nudg_prof_cas(nlev_cas)                  !jyg
1117
1118        thv_mod_cas(l) = thv_prof_cas(nlev_cas)                        !jyg
1119        w_mod_cas(l) = 0.0                                             !jyg
1120        omega_mod_cas(l) = 0.0                                         !jyg
1121        du_mod_cas(l) = du_prof_cas(nlev_cas) * fact
1122        hu_mod_cas(l) = hu_prof_cas(nlev_cas) * fact                     !jyg
1123        vu_mod_cas(l) = vu_prof_cas(nlev_cas) * fact                     !jyg
1124        dv_mod_cas(l) = dv_prof_cas(nlev_cas) * fact
1125        hv_mod_cas(l) = hv_prof_cas(nlev_cas) * fact                     !jyg
1126        vv_mod_cas(l) = vv_prof_cas(nlev_cas) * fact                     !jyg
1127        dt_mod_cas(l) = dt_prof_cas(nlev_cas)
1128        ht_mod_cas(l) = ht_prof_cas(nlev_cas)                          !jyg
1129        vt_mod_cas(l) = vt_prof_cas(nlev_cas)                          !jyg
1130        dth_mod_cas(l) = dth_prof_cas(nlev_cas)
1131        hth_mod_cas(l) = hth_prof_cas(nlev_cas)                        !jyg
1132        vth_mod_cas(l) = vth_prof_cas(nlev_cas)                        !jyg
1133        dq_mod_cas(l) = dq_prof_cas(nlev_cas) * fact
1134        hq_mod_cas(l) = hq_prof_cas(nlev_cas) * fact                     !jyg
1135        vq_mod_cas(l) = vq_prof_cas(nlev_cas) * fact                     !jyg
1136        dtrad_mod_cas(l) = dtrad_prof_cas(nlev_cas) * fact               !jyg
1137
1138      endif ! play
1139
1140    enddo ! l
1141
1142    ! for variables defined at layer interfaces (EV):
1143
1144    DO l = 1, llm + 1
1145
1146      IF (plev(l)>=plev_prof_cas(nlev_cas)) THEN
1147        mxcalc = l
1148        k1 = 0
1149        k2 = 0
1150
1151        IF (plev(l)<=plev_prof_cas(1)) THEN
1152          DO k = 1, nlev_cas - 1
1153            IF (plev(l)<=plev_prof_cas(k).AND. plev(l)>plev_prof_cas(k + 1)) THEN
1154              k1 = k
1155              k2 = k + 1
1156            endif
1157          enddo
1158
1159          IF (k1==0 .OR. k2==0) THEN
1160            WRITE(*, *) 'PB! k1, k2 = ', k1, k2
1161            WRITE(*, *) 'l,plev(l) = ', l, plev(l) / 100
1162            DO k = 1, nlev_cas - 1
1163              WRITE(*, *) 'k,plev_prof_cas(k) = ', k, plev_prof_cas(k) / 100
1164            enddo
1165          endif
1166
1167          frac = (plev_prof_cas(k2) - plev(l)) / (plev_prof_cas(k2) - plev_prof_cas(k1))
1168          tke_mod_cas(l) = tke_prof_cas(k2) - frac * (tke_prof_cas(k2) - tke_prof_cas(k1))
1169        else !play>plev_prof_cas(1)
1170          k1 = 1
1171          k2 = 2
1172          frac1 = (play(l) - plev_prof_cas(k2)) / (plev_prof_cas(k1) - plev_prof_cas(k2))
1173          frac2 = (play(l) - plev_prof_cas(k1)) / (plev_prof_cas(k1) - plev_prof_cas(k2))
1174          tke_mod_cas(l) = frac1 * tke_prof_cas(k1) - frac2 * tke_prof_cas(k2)
1175
1176        endif ! plev.le.plev_prof_cas(1)
1177
1178      else ! above max altitude of forcing file
1179
1180        tke_mod_cas(l) = 0.0
1181
1182      endif ! plev
1183
1184    enddo ! l
1185
1186  END SUBROUTINE  interp2_case_vertical_std
1187  !*****************************************************************************
1188
1189
1190END MODULE mod_1D_cases_read_std
Note: See TracBrowser for help on using the repository browser.