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

Last change on this file since 4274 was 4274, checked in by lguez, 20 months ago

Replace nf_get_var_type by nf90_get_var

The immediate motivation is a bug fix: nf_get_var_type was called
with scalar resul3 or lat, lon, alt, phis instead of array actual
argument for dummy array argument rvals or dvals. Correcting this, we
might as well take the opportunity to use nf90_get_var, so we no
longer need to test NC_DOUBLE and we have half as many calls.

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