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

Last change on this file since 4115 was 4104, checked in by evignon, 2 years ago

modification des routines pour lecture de la v1 des cas 1D
au format standard + corrections de bug pour le geopotentiel à la surface

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