Index: LMDZ5/trunk/libf/phylmd/1DUTILS.h =================================================================== --- LMDZ5/trunk/libf/phylmd/1DUTILS.h (revision 2190) +++ LMDZ5/trunk/libf/phylmd/1DUTILS.h (revision 2191) @@ -99,4 +99,6 @@ ! LS convergence imposed from RICO (cst) ! = 6 ==> forcing_amma = .true. +! = 10 ==> forcing_case = .true. +! initial profiles from case.nc file ! = 40 ==> forcing_GCSSold = .true. ! initial profile from GCSS file @@ -132,4 +134,84 @@ CALL getin('turb_fcg',xTurb_fcg_gcssold) ENDIF + +!Paramètres de forçage +!Config Key = tend_t +!Config Desc = forcage ou non par advection de T +!Config Def = false +!Config Help = forcage ou non par advection de T + tend_t =0 + CALL getin('tend_t',tend_t) + +!Config Key = tend_q +!Config Desc = forcage ou non par advection de q +!Config Def = false +!Config Help = forcage ou non par advection de q + tend_q =0 + CALL getin('tend_q',tend_q) + +!Config Key = tend_u +!Config Desc = forcage ou non par advection de u +!Config Def = false +!Config Help = forcage ou non par advection de u + tend_u =0 + CALL getin('tend_u',tend_u) + +!Config Key = tend_v +!Config Desc = forcage ou non par advection de v +!Config Def = false +!Config Help = forcage ou non par advection de v + tend_v =0 + CALL getin('tend_v',tend_v) + +!Config Key = tend_w +!Config Desc = forcage ou non par vitesse verticale +!Config Def = false +!Config Help = forcage ou non par vitesse verticale + tend_w =0 + CALL getin('tend_w',tend_w) + +!Config Key = tend_rayo +!Config Desc = forcage ou non par dtrad +!Config Def = false +!Config Help = forcage ou non par dtrad + tend_rayo =0 + CALL getin('tend_rayo',tend_rayo) + + +!Config Key = nudge_t +!Config Desc = constante de nudging de T +!Config Def = false +!Config Help = constante de nudging de T + nudge_t =0. + CALL getin('nudge_t',nudge_t) + +!Config Key = nudge_q +!Config Desc = constante de nudging de q +!Config Def = false +!Config Help = constante de nudging de q + nudge_q =0. + CALL getin('nudge_q',nudge_q) + +!Config Key = nudge_u +!Config Desc = constante de nudging de u +!Config Def = false +!Config Help = constante de nudging de u + nudge_u =0. + CALL getin('nudge_u',nudge_u) + +!Config Key = nudge_v +!Config Desc = constante de nudging de v +!Config Def = false +!Config Help = constante de nudging de v + nudge_v =0. + CALL getin('nudge_v',nudge_v) + +!Config Key = nudge_w +!Config Desc = constante de nudging de w +!Config Def = false +!Config Help = constante de nudging de w + nudge_w =0. + CALL getin('nudge_w',nudge_w) + !Config Key = iflag_nudge @@ -2431,4 +2513,162 @@ !===================================================================== + SUBROUTINE interp_case_vertical(play,nlev_cas,plev_prof_cas & + & ,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas & + & ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas & + & ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas & + & ,t_mod_cas,q_mod_cas,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas & + & ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas & + & ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas,mxcalc) + + implicit none + +#include "dimensions.h" + +!------------------------------------------------------------------------- +! Vertical interpolation of TOGA-COARE forcing data onto mod_casel levels +!------------------------------------------------------------------------- + + integer nlevmax + parameter (nlevmax=41) + integer nlev_cas,mxcalc +! real play(llm), plev_prof(nlevmax) +! real t_prof(nlevmax),q_prof(nlevmax) +! real u_prof(nlevmax),v_prof(nlevmax), w_prof(nlevmax) +! real ht_prof(nlevmax),vt_prof(nlevmax) +! real hq_prof(nlevmax),vq_prof(nlevmax) + + real play(llm), plev_prof_cas(nlev_cas) + real t_prof_cas(nlev_cas),q_prof_cas(nlev_cas) + real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas) + real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas), vitw_prof_cas(nlev_cas) + real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas) + real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas) + real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas),dtrad_prof_cas(nlev_cas) + real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas) + + real t_mod_cas(llm),q_mod_cas(llm) + real u_mod_cas(llm),v_mod_cas(llm) + real ug_mod_cas(llm),vg_mod_cas(llm), w_mod_cas(llm) + real du_mod_cas(llm),hu_mod_cas(llm),vu_mod_cas(llm) + real dv_mod_cas(llm),hv_mod_cas(llm),vv_mod_cas(llm) + real dt_mod_cas(llm),ht_mod_cas(llm),vt_mod_cas(llm),dtrad_mod_cas(llm) + real dq_mod_cas(llm),hq_mod_cas(llm),vq_mod_cas(llm) + + integer l,k,k1,k2 + real frac,frac1,frac2,fact + + do l = 1, llm + + if (play(l).ge.plev_prof_cas(nlev_cas)) then + + mxcalc=l + k1=0 + k2=0 + + if (play(l).le.plev_prof_cas(1)) then + + do k = 1, nlev_cas-1 + if (play(l).le.plev_prof_cas(k).and. play(l).gt.plev_prof_cas(k+1)) then + k1=k + k2=k+1 + endif + enddo + + if (k1.eq.0 .or. k2.eq.0) then + write(*,*) 'PB! k1, k2 = ',k1,k2 + write(*,*) 'l,play(l) = ',l,play(l)/100 + do k = 1, nlev_cas-1 + write(*,*) 'k,plev_prof_cas(k) = ',k,plev_prof_cas(k)/100 + enddo + endif + + frac = (plev_prof_cas(k2)-play(l))/(plev_prof_cas(k2)-plev_prof_cas(k1)) + t_mod_cas(l)= t_prof_cas(k2) - frac*(t_prof_cas(k2)-t_prof_cas(k1)) + q_mod_cas(l)= q_prof_cas(k2) - frac*(q_prof_cas(k2)-q_prof_cas(k1)) + u_mod_cas(l)= u_prof_cas(k2) - frac*(u_prof_cas(k2)-u_prof_cas(k1)) + v_mod_cas(l)= v_prof_cas(k2) - frac*(v_prof_cas(k2)-v_prof_cas(k1)) + ug_mod_cas(l)= ug_prof_cas(k2) - frac*(ug_prof_cas(k2)-ug_prof_cas(k1)) + vg_mod_cas(l)= vg_prof_cas(k2) - frac*(vg_prof_cas(k2)-vg_prof_cas(k1)) + w_mod_cas(l)= vitw_prof_cas(k2) - frac*(vitw_prof_cas(k2)-vitw_prof_cas(k1)) + du_mod_cas(l)= du_prof_cas(k2) - frac*(du_prof_cas(k2)-du_prof_cas(k1)) + hu_mod_cas(l)= hu_prof_cas(k2) - frac*(hu_prof_cas(k2)-hu_prof_cas(k1)) + vu_mod_cas(l)= vu_prof_cas(k2) - frac*(vu_prof_cas(k2)-vu_prof_cas(k1)) + dv_mod_cas(l)= dv_prof_cas(k2) - frac*(dv_prof_cas(k2)-dv_prof_cas(k1)) + hv_mod_cas(l)= hv_prof_cas(k2) - frac*(hv_prof_cas(k2)-hv_prof_cas(k1)) + vv_mod_cas(l)= vv_prof_cas(k2) - frac*(vv_prof_cas(k2)-vv_prof_cas(k1)) + dt_mod_cas(l)= dt_prof_cas(k2) - frac*(dt_prof_cas(k2)-dt_prof_cas(k1)) + ht_mod_cas(l)= ht_prof_cas(k2) - frac*(ht_prof_cas(k2)-ht_prof_cas(k1)) + vt_mod_cas(l)= vt_prof_cas(k2) - frac*(vt_prof_cas(k2)-vt_prof_cas(k1)) + dq_mod_cas(l)= dq_prof_cas(k2) - frac*(dq_prof_cas(k2)-dq_prof_cas(k1)) + hq_mod_cas(l)= hq_prof_cas(k2) - frac*(hq_prof_cas(k2)-hq_prof_cas(k1)) + vq_mod_cas(l)= vq_prof_cas(k2) - frac*(vq_prof_cas(k2)-vq_prof_cas(k1)) + + else !play>plev_prof_cas(1) + + k1=1 + k2=2 + frac1 = (play(l)-plev_prof_cas(k2))/(plev_prof_cas(k1)-plev_prof_cas(k2)) + frac2 = (play(l)-plev_prof_cas(k1))/(plev_prof_cas(k1)-plev_prof_cas(k2)) + t_mod_cas(l)= frac1*t_prof_cas(k1) - frac2*t_prof_cas(k2) + q_mod_cas(l)= frac1*q_prof_cas(k1) - frac2*q_prof_cas(k2) + u_mod_cas(l)= frac1*u_prof_cas(k1) - frac2*u_prof_cas(k2) + v_mod_cas(l)= frac1*v_prof_cas(k1) - frac2*v_prof_cas(k2) + ug_mod_cas(l)= frac1*ug_prof_cas(k1) - frac2*ug_prof_cas(k2) + vg_mod_cas(l)= frac1*vg_prof_cas(k1) - frac2*vg_prof_cas(k2) + w_mod_cas(l)= frac1*vitw_prof_cas(k1) - frac2*vitw_prof_cas(k2) + du_mod_cas(l)= frac1*du_prof_cas(k1) - frac2*du_prof_cas(k2) + hu_mod_cas(l)= frac1*hu_prof_cas(k1) - frac2*hu_prof_cas(k2) + vu_mod_cas(l)= frac1*vu_prof_cas(k1) - frac2*vu_prof_cas(k2) + dv_mod_cas(l)= frac1*dv_prof_cas(k1) - frac2*dv_prof_cas(k2) + hv_mod_cas(l)= frac1*hv_prof_cas(k1) - frac2*hv_prof_cas(k2) + vv_mod_cas(l)= frac1*vv_prof_cas(k1) - frac2*vv_prof_cas(k2) + dt_mod_cas(l)= frac1*dt_prof_cas(k1) - frac2*dt_prof_cas(k2) + ht_mod_cas(l)= frac1*ht_prof_cas(k1) - frac2*ht_prof_cas(k2) + vt_mod_cas(l)= frac1*vt_prof_cas(k1) - frac2*vt_prof_cas(k2) + dq_mod_cas(l)= frac1*dq_prof_cas(k1) - frac2*dq_prof_cas(k2) + hq_mod_cas(l)= frac1*hq_prof_cas(k1) - frac2*hq_prof_cas(k2) + vq_mod_cas(l)= frac1*vq_prof_cas(k1) - frac2*vq_prof_cas(k2) + + endif ! play.le.plev_prof_cas(1) + + else ! above max altitude of forcing file + +!jyg + fact=20.*(plev_prof_cas(nlev_cas)-play(l))/plev_prof_cas(nlev_cas) !jyg + fact = max(fact,0.) !jyg + fact = exp(-fact) !jyg + t_mod_cas(l)= t_prof_cas(nlev_cas) !jyg + q_mod_cas(l)= q_prof_cas(nlev_cas)*fact !jyg + u_mod_cas(l)= u_prof_cas(nlev_cas)*fact !jyg + v_mod_cas(l)= v_prof_cas(nlev_cas)*fact !jyg + ug_mod_cas(l)= ug_prof_cas(nlev_cas)*fact !jyg + vg_mod_cas(l)= vg_prof_cas(nlev_cas)*fact !jyg + w_mod_cas(l)= 0.0 !jyg + du_mod_cas(l)= du_prof_cas(nlev_cas)*fact + hu_mod_cas(l)= hu_prof_cas(nlev_cas)*fact !jyg + vu_mod_cas(l)= vu_prof_cas(nlev_cas)*fact !jyg + dv_mod_cas(l)= dv_prof_cas(nlev_cas)*fact + hv_mod_cas(l)= hv_prof_cas(nlev_cas)*fact !jyg + vv_mod_cas(l)= vv_prof_cas(nlev_cas)*fact !jyg + dt_mod_cas(l)= dt_prof_cas(nlev_cas) + ht_mod_cas(l)= ht_prof_cas(nlev_cas) !jyg + vt_mod_cas(l)= vt_prof_cas(nlev_cas) !jyg + dq_mod_cas(l)= dq_prof_cas(nlev_cas)*fact + hq_mod_cas(l)= hq_prof_cas(nlev_cas)*fact !jyg + vq_mod_cas(l)= vq_prof_cas(nlev_cas)*fact !jyg + + endif ! play + + enddo ! l + +! do l = 1,llm +! print *,'t_mod_cas(l),q_mod_cas(l),ht_mod_cas(l),hq_mod_cas(l) ', +! $ l,t_mod_cas(l),q_mod_cas(l),ht_mod_cas(l),hq_mod_cas(l) +! enddo + + return + end +!***************************************************************************** +!===================================================================== SUBROUTINE interp_dice_vertical(play,nlev_dice,nt_dice,plev_prof & & ,th_prof,qv_prof,u_prof,v_prof,o3_prof & Index: LMDZ5/trunk/libf/phylmd/1D_decl_cases.h =================================================================== --- LMDZ5/trunk/libf/phylmd/1D_decl_cases.h (revision 2190) +++ LMDZ5/trunk/libf/phylmd/1D_decl_cases.h (revision 2191) @@ -240,4 +240,15 @@ real u_profa(nlev_astex),v_profa(nlev_astex),w_profa(nlev_astex) real tke_profa(nlev_astex) -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +!Declarations specifiques au cas standard + + real w_mod_cas(llm), t_mod_cas(llm),q_mod_cas(llm) + real ug_mod_cas(llm),vg_mod_cas(llm) + real u_mod_cas(llm),v_mod_cas(llm) + real ht_mod_cas(llm),vt_mod_cas(llm),dt_mod_cas(llm),dtrad_mod_cas(llm) + real hq_mod_cas(llm),vq_mod_cas(llm),dq_mod_cas(llm) + real hu_mod_cas(llm),vu_mod_cas(llm),du_mod_cas(llm) + real hv_mod_cas(llm),vv_mod_cas(llm),dv_mod_cas(llm) +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + Index: LMDZ5/trunk/libf/phylmd/1D_interp_cases.h =================================================================== --- LMDZ5/trunk/libf/phylmd/1D_interp_cases.h (revision 2190) +++ LMDZ5/trunk/libf/phylmd/1D_interp_cases.h (revision 2191) @@ -597,4 +597,161 @@ enddo endif ! forcing_astex -!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +!--------------------------------------------------------------------- +! Interpolation forcing standard case +!--------------------------------------------------------------------- + if (forcing_case) then + + print*, & + & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/pdt_cas=', & + & daytime,day1,(daytime-day1)*86400., & + & (daytime-day1)*86400/pdt_cas + +! time interpolation: + CALL interp_case_time(daytime,day1,annee_ref & + & ,year_ini_cas,day_ju_ini_cas,nt_cas,pdt_cas,nlev_cas & + & ,ts_cas,plev_cas,t_cas,q_cas,u_cas,v_cas,ug_cas,vg_cas & + & ,vitw_cas,du_cas,hu_cas,vu_cas & + & ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas & + & ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas & + & ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas & + & ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas & + & ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas & + & ,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas & + & ,sens_prof_cas) + + ts_cur = ts_prof_cas + psurf=plev_prof_cas(1) + +! vertical interpolation: + CALL interp_case_vertical(play,nlev_cas,plev_prof_cas & + & ,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas & + & ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas & + & ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas & + & ,t_mod_cas,q_mod_cas,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas & + & ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas & + & ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas,mxcalc) + + +!calcul de l'advection verticale a partir du omega +!Calcul des gradients verticaux +!initialisation + d_t_z(:)=0. + d_q_z(:)=0. + d_u_z(:)=0. + d_v_z(:)=0. + d_t_dyn_z(:)=0. + d_q_dyn_z(:)=0. + d_u_dyn_z(:)=0. + d_v_dyn_z(:)=0. + DO l=2,llm-1 + d_t_z(l)=(temp(l+1)-temp(l-1))/(play(l+1)-play(l-1)) + d_q_z(l)=(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1)) + d_u_z(l)=(u(l+1)-u(l-1))/(play(l+1)-play(l-1)) + d_v_z(l)=(v(l+1)-v(l-1))/(play(l+1)-play(l-1)) + ENDDO + d_t_z(1)=d_t_z(2) + d_q_z(1)=d_q_z(2) + d_u_z(1)=d_u_z(2) + d_v_z(1)=d_v_z(2) + d_t_z(llm)=d_t_z(llm-1) + d_q_z(llm)=d_q_z(llm-1) + d_u_z(llm)=d_u_z(llm-1) + d_v_z(llm)=d_v_z(llm-1) + +!Calcul de l advection verticale + d_t_dyn_z(:)=w_mod_cas(:)*d_t_z(:) + d_q_dyn_z(:)=w_mod_cas(:)*d_q_z(:) + d_u_dyn_z(:)=w_mod_cas(:)*d_u_z(:) + d_v_dyn_z(:)=w_mod_cas(:)*d_v_z(:) + +!wind nudging + if (nudge_u.gt.0.) then + do l=1,llm + u(l)=u(l)+timestep*(u_mod_cas(l)-u(l))/(nudge_u) + enddo + else + do l=1,llm + u(l) = u_mod_cas(l) + enddo + endif + + if (nudge_v.gt.0.) then + do l=1,llm + v(l)=v(l)+timestep*(v_mod_cas(l)-v(l))/(nudge_v) + enddo + else + do l=1,llm + v(l) = v_mod_cas(l) + enddo + endif + + if (nudge_w.gt.0.) then + do l=1,llm + w(l)=w(l)+timestep*(w_mod_cas(l)-w(l))/(nudge_w) + enddo + else + do l=1,llm + w(l) = w_mod_cas(l) + enddo + endif + +!nudging of q and temp + if (nudge_t.gt.0.) then + do l=1,llm + temp(l)=temp(l)+timestep*(t_mod_cas(l)-temp(l))/(nudge_t) + enddo + endif + if (nudge_q.gt.0.) then + do l=1,llm + q(l,1)=q(l,1)+timestep*(q_mod_cas(l)-q(l,1))/(nudge_q) + enddo + endif + + do l = 1, llm + omega(l) = w_mod_cas(l) + omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq + alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l) + +!calcul advection + if ((tend_u.eq.1).and.(tend_w.eq.0)) then + d_u_adv(l)=du_mod_cas(l) + else if ((tend_u.eq.1).and.(tend_w.eq.1)) then + d_u_adv(l)=hu_mod_cas(l)-d_u_dyn_z(l) + endif + + if ((tend_v.eq.1).and.(tend_w.eq.0)) then + d_v_adv(l)=dv_mod_cas(l) + else if ((tend_v.eq.1).and.(tend_w.eq.1)) then + d_v_adv(l)=hv_mod_cas(l)-d_v_dyn_z(l) + endif + + if ((tend_t.eq.1).and.(tend_w.eq.0)) then +! d_th_adv(l)=alpha*omega(l)/rcpd+dt_mod_cas(l) + d_th_adv(l)=alpha*omega(l)/rcpd-dt_mod_cas(l) + else if ((tend_t.eq.1).and.(tend_w.eq.1)) then +! d_th_adv(l)=alpha*omega(l)/rcpd+ht_mod_cas(l)-d_t_dyn_z(l) + d_th_adv(l)=alpha*omega(l)/rcpd-ht_mod_cas(l)-d_t_dyn_z(l) + endif + + if ((tend_q.eq.1).and.(tend_w.eq.0)) then +! d_q_adv(l,1)=dq_mod_cas(l) + d_q_adv(l,1)=-1*dq_mod_cas(l) + else if ((tend_q.eq.1).and.(tend_w.eq.1)) then +! d_q_adv(l,1)=hq_mod_cas(l)-d_q_dyn_z(l) + d_q_adv(l,1)=-1*hq_mod_cas(l)-d_q_dyn_z(l) + endif + + if (tend_rayo.eq.1) then + dt_cooling(l) = dtrad_mod_cas(l) + else + dt_cooling(l) = 0.0 + endif + enddo + + endif ! forcing_case + + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + Index: LMDZ5/trunk/libf/phylmd/1D_read_forc_cases.h =================================================================== --- LMDZ5/trunk/libf/phylmd/1D_read_forc_cases.h (revision 2190) +++ LMDZ5/trunk/libf/phylmd/1D_read_forc_cases.h (revision 2191) @@ -720,3 +720,63 @@ endif ! forcing_astex - +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +!--------------------------------------------------------------------- +! Forcing from standard case : +!--------------------------------------------------------------------- + + if (forcing_case) then + + write(*,*),'avant call read_1D_cas' + call read_1D_cas + write(*,*) 'Forcing read' + +!Time interpolation for initial conditions using TOGA interpolation routine + write(*,*) 'AVT 1ere INTERPOLATION: day,day1 = ',daytime,day1 + CALL interp_case_time(day,day1,annee_ref & + & ,year_ini_cas,day_ju_ini_cas,nt_cas,pdt_cas,nlev_cas & + & ,ts_cas,plev_cas,t_cas,q_cas,u_cas,v_cas & + & ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas & + & ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas & + & ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas & + & ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas & + & ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas & + & ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas & + & ,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas) + +! vertical interpolation using TOGA interpolation routine: +! write(*,*)'avant interp vert', t_prof + CALL interp_case_vertical(play,nlev_cas,plev_prof_cas & + & ,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas & + & ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas & + & ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas & + & ,t_mod_cas,q_mod_cas,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas & + & ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas & + & ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas,mxcalc) +! write(*,*) 'Profil initial forcing case interpole',t_mod + +! initial and boundary conditions : +! tsurf = ts_prof_cas + ts_cur = ts_prof_cas + psurf=plev_prof_cas(1) + write(*,*) 'SST initiale: ',tsurf + do l = 1, llm + temp(l) = t_mod_cas(l) + q(l,1) = q_mod_cas(l) + q(l,2) = 0.0 + u(l) = u_mod_cas(l) + v(l) = v_mod_cas(l) + omega(l) = w_mod_cas(l) + omega2(l)=omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq + + alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l) +!on applique le forcage total au premier pas de temps +!attention: signe different de toga + d_th_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l)) + d_q_adv(l,1) = (hq_mod_cas(l)+vq_mod_cas(l)) + d_q_adv(l,2) = 0.0 + d_u_adv(l) = (hu_mod_cas(l)+vu_mod_cas(l)) + d_u_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l)) + enddo + + endif !forcing_case +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! Index: LMDZ5/trunk/libf/phylmd/compar1d.h =================================================================== --- LMDZ5/trunk/libf/phylmd/compar1d.h (revision 2190) +++ LMDZ5/trunk/libf/phylmd/compar1d.h (revision 2191) @@ -3,4 +3,6 @@ ! integer :: forcing_type + integer :: tend_u,tend_v,tend_w,tend_t,tend_q,tend_rayo + real :: nudge_u,nudge_v,nudge_w,nudge_t,nudge_q integer :: iflag_nudge real :: nat_surf @@ -32,5 +34,6 @@ & qsol,qsurf,psurf,zsurf,albedo,time,time_ini,xlat,xlon,airefi, & & wtsurf,wqsurf,restart_runoff,xagesno,qsolinp,zpicinp, & - & forcing_type, & + & forcing_type,tend_u,tend_v,tend_w,tend_t,tend_q,tend_rayo, & + & nudge_u,nudge_v,nudge_w,nudge_t,nudge_q, & & iflag_nudge, & & restart,ok_old_disvert Index: LMDZ5/trunk/libf/phylmd/lmdz1d.F90 =================================================================== --- LMDZ5/trunk/libf/phylmd/lmdz1d.F90 (revision 2190) +++ LMDZ5/trunk/libf/phylmd/lmdz1d.F90 (revision 2191) @@ -22,4 +22,5 @@ USE phyaqua_mod USE mod_1D_cases_read + USE mod_1D_amma_read implicit none @@ -118,4 +119,5 @@ logical :: forcing_astex = .false. logical :: forcing_fire = .false. + logical :: forcing_case = .false. integer :: type_ts_forcing ! 0 = SST constant; 1 = SST read from a file ! (cf read_tsurf1d.F) @@ -171,4 +173,5 @@ real :: dt_cooling(llm),d_th_adv(llm),d_t_nudge(llm) real :: d_u_nudge(llm),d_v_nudge(llm) + real :: du_adv(llm),dv_adv(llm) real :: alpha real :: ttt @@ -285,4 +288,7 @@ ! Different stages: soil model alone, atm. model alone ! then both models coupled +!forcing_type = 10 ==> forcing_case = .true. +! initial profiles and large scale forcings in cas.nc +! LS convergence, omega and SST imposed from CINDY-DYNAMO files !forcing_type = 40 ==> forcing_GCSSold = .true. ! initial profile from GCSS file @@ -321,4 +327,6 @@ elseif (forcing_type .eq.7) THEN forcing_dice = .true. + elseif (forcing_type .eq.10) THEN + forcing_case = .true. elseif (forcing_type .eq.40) THEN forcing_GCSSold = .true. @@ -428,4 +436,11 @@ & (year_ini_dice,mth_ini_dice,day_ini_dice,heure_ini_dice & & ,day_ju_ini_dice) + ELSEIF (forcing_type .eq.10) THEN +! Convert the initial date to Julian day + print*,'time cindy',year_ini_cas,mth_ini_cas,day_ini_cas + call ymds2ju & + & (year_ini_cas,mth_ini_cas,day_ini_cas,heure_ini_cas & + & ,day_ju_ini_cas) + print*,'time cindy 2',day_ju_ini_cas ELSEIF (forcing_type .eq.59) THEN ! Convert the initial date of Sandu case to Julian day @@ -943,9 +958,9 @@ u(1:mxcalc)=u(1:mxcalc) + timestep*( & & du_phys(1:mxcalc) & - & +du_age(1:mxcalc) & + & +du_age(1:mxcalc)+du_adv(1:mxcalc) & & +d_u_nudge(1:mxcalc) ) v(1:mxcalc)=v(1:mxcalc) + timestep*( & & dv_phys(1:mxcalc) & - & +dv_age(1:mxcalc) & + & +dv_age(1:mxcalc)+dv_adv(1:mxcalc) & & +d_v_nudge(1:mxcalc) ) q(1:mxcalc,:)=q(1:mxcalc,:)+timestep*( & Index: LMDZ5/trunk/libf/phylmd/mod_1D_amma_read.F90 =================================================================== --- LMDZ5/trunk/libf/phylmd/mod_1D_amma_read.F90 (revision 2191) +++ LMDZ5/trunk/libf/phylmd/mod_1D_amma_read.F90 (revision 2191) @@ -0,0 +1,500 @@ +MODULE mod_1D_amma_read + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +!Declarations specifiques au cas AMMA + character*80 :: fich_amma +! Option du cas AMMA ou on impose la discretisation verticale (Ap,Bp) + integer nlev_amma, nt_amma + + + integer year_ini_amma, day_ini_amma, mth_ini_amma + real heure_ini_amma + real day_ju_ini_amma ! Julian day of amma first day + parameter (year_ini_amma=2006) + parameter (mth_ini_amma=7) + parameter (day_ini_amma=10) ! 10 = 10Juil2006 + parameter (heure_ini_amma=0.) !0h en secondes + real dt_amma + parameter (dt_amma=1800.) + +!profils initiaux: + real, allocatable:: plev_amma(:) + + real, allocatable:: z_amma(:) + real, allocatable:: th_amma(:),q_amma(:) + real, allocatable:: u_amma(:) + real, allocatable:: v_amma(:) + + real, allocatable:: th_ammai(:),q_ammai(:) + real, allocatable:: u_ammai(:) + real, allocatable:: v_ammai(:) + real, allocatable:: vitw_ammai(:) + real, allocatable:: ht_ammai(:) + real, allocatable:: hq_ammai(:) + real, allocatable:: vt_ammai(:) + real, allocatable:: vq_ammai(:) + +!forcings + real, allocatable:: ht_amma(:,:) + real, allocatable:: hq_amma(:,:) + real, allocatable:: vitw_amma(:,:) + real, allocatable:: lat_amma(:),sens_amma(:) + +!champs interpoles + real, allocatable:: vitw_profamma(:) + real, allocatable:: ht_profamma(:) + real, allocatable:: hq_profamma(:) + real lat_profamma,sens_profamma + real, allocatable:: vt_profamma(:) + real, allocatable:: vq_profamma(:) + real, allocatable:: th_profamma(:) + real, allocatable:: q_profamma(:) + real, allocatable:: u_profamma(:) + real, allocatable:: v_profamma(:) + + +CONTAINS + +SUBROUTINE read_1D_cases + implicit none + +#include "netcdf.inc" + + INTEGER nid,rid,ierr + + fich_amma='amma.nc' + print*,'fich_amma ',fich_amma + ierr = NF_OPEN(fich_amma,NF_NOWRITE,nid) + print*,'fich_amma,NF_NOWRITE,nid ',fich_amma,NF_NOWRITE,nid + if (ierr.NE.NF_NOERR) then + write(*,*) 'ERROR: GROS Pb opening forcings nc file ' + write(*,*) NF_STRERROR(ierr) + stop "" + endif +!....................................................................... + ierr=NF_INQ_DIMID(nid,'lev',rid) + IF (ierr.NE.NF_NOERR) THEN + print*, 'Oh probleme lecture dimension zz' + ENDIF + ierr=NF_INQ_DIMLEN(nid,rid,nlev_amma) + print*,'OK nid,rid,nlev_amma',nid,rid,nlev_amma +!....................................................................... + ierr=NF_INQ_DIMID(nid,'time',rid) + print*,'nid,rid',nid,rid + nt_amma=0 + IF (ierr.NE.NF_NOERR) THEN + stop 'probleme lecture dimension sens' + ENDIF + ierr=NF_INQ_DIMLEN(nid,rid,nt_amma) + print*,'nid,rid,nlev_amma',nid,rid,nt_amma + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +!profils initiaux: + allocate(plev_amma(nlev_amma)) + + allocate(z_amma(nlev_amma)) + allocate(th_amma(nlev_amma),q_amma(nlev_amma)) + allocate(u_amma(nlev_amma)) + allocate(v_amma(nlev_amma)) + +!forcings + allocate(ht_amma(nlev_amma,nt_amma)) + allocate(hq_amma(nlev_amma,nt_amma)) + allocate(vitw_amma(nlev_amma,nt_amma)) + allocate(lat_amma(nt_amma),sens_amma(nt_amma)) + +!profils initiaux: + allocate(th_ammai(nlev_amma),q_ammai(nlev_amma)) + allocate(u_ammai(nlev_amma)) + allocate(v_ammai(nlev_amma)) + allocate(vitw_ammai(nlev_amma) ) + allocate(ht_ammai(nlev_amma)) + allocate(hq_ammai(nlev_amma)) + allocate(vt_ammai(nlev_amma)) + allocate(vq_ammai(nlev_amma)) + +!champs interpoles + allocate(vitw_profamma(nlev_amma)) + allocate(ht_profamma(nlev_amma)) + allocate(hq_profamma(nlev_amma)) + allocate(vt_profamma(nlev_amma)) + allocate(vq_profamma(nlev_amma)) + allocate(th_profamma(nlev_amma)) + allocate(q_profamma(nlev_amma)) + allocate(u_profamma(nlev_amma)) + allocate(v_profamma(nlev_amma)) + + print*,'Allocations OK' + call read_amma(nid,nlev_amma,nt_amma & + & ,z_amma,plev_amma,th_amma,q_amma,u_amma,v_amma,vitw_amma & + & ,ht_amma,hq_amma,sens_amma,lat_amma) + +END SUBROUTINE read_1D_cases + + + +!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! +SUBROUTINE deallocate_1D_cases +!profils initiaux: + deallocate(plev_amma) + + deallocate(z_amma) + deallocate(th_amma,q_amma) + deallocate(u_amma) + deallocate(v_amma) + + deallocate(th_ammai,q_ammai) + deallocate(u_ammai) + deallocate(v_ammai) + deallocate(vitw_ammai ) + deallocate(ht_ammai) + deallocate(hq_ammai) + deallocate(vt_ammai) + deallocate(vq_ammai) + +!forcings + deallocate(ht_amma) + deallocate(hq_amma) + deallocate(vitw_amma) + deallocate(lat_amma,sens_amma) + +!champs interpoles + deallocate(vitw_profamma) + deallocate(ht_profamma) + deallocate(hq_profamma) + deallocate(vt_profamma) + deallocate(vq_profamma) + deallocate(th_profamma) + deallocate(q_profamma) + deallocate(u_profamma) + deallocate(v_profamma) +END SUBROUTINE deallocate_1D_cases + + +END MODULE mod_1D_amma_read +!===================================================================== + subroutine read_amma(nid,nlevel,ntime & + & ,zz,pp,temp,qv,u,v,dw & + & ,dt,dq,sens,flat) + +!program reading forcings of the AMMA case study + implicit none +#include "netcdf.inc" + + integer ntime,nlevel + + real zz(nlevel) + real temp(nlevel),pp(nlevel) + real qv(nlevel),u(nlevel) + real v(nlevel) + real dw(nlevel,ntime) + real dt(nlevel,ntime) + real dq(nlevel,ntime) + real flat(ntime),sens(ntime) + + + integer nid, ierr,rid + integer nbvar3d + parameter(nbvar3d=30) + integer var3didin(nbvar3d) + + ierr=NF_INQ_VARID(nid,"zz",var3didin(1)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'lev' + endif + + + ierr=NF_INQ_VARID(nid,"temp",var3didin(2)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'temp' + endif + + ierr=NF_INQ_VARID(nid,"qv",var3didin(3)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'qv' + endif + + ierr=NF_INQ_VARID(nid,"u",var3didin(4)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'u' + endif + + ierr=NF_INQ_VARID(nid,"v",var3didin(5)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'v' + endif + + ierr=NF_INQ_VARID(nid,"dw",var3didin(6)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'dw' + endif + + ierr=NF_INQ_VARID(nid,"dt",var3didin(7)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'dt' + endif + + ierr=NF_INQ_VARID(nid,"dq",var3didin(8)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'dq' + endif + + ierr=NF_INQ_VARID(nid,"sens",var3didin(9)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'sens' + endif + + ierr=NF_INQ_VARID(nid,"flat",var3didin(10)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'flat' + endif + + ierr=NF_INQ_VARID(nid,"pp",var3didin(11)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + endif + +!dimensions lecture +! call catchaxis(nid,ntime,nlevel,time,z,ierr) + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(1),zz) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(1),zz) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture z ok',zz + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(2),temp) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(2),temp) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture th ok',temp + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(3),qv) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(3),qv) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture qv ok',qv + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(4),u) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(4),u) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture u ok',u + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(5),v) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(5),v) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture v ok',v + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(6),dw) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(6),dw) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture w ok',dw + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(7),dt) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(7),dt) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture dt ok',dt + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(8),dq) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(8),dq) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture dq ok',dq + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(9),sens) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(9),sens) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture sens ok',sens + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(10),flat) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(10),flat) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture flat ok',flat + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(11),pp) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(11),pp) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture pp ok',pp + + return + end subroutine read_amma +!====================================================================== + SUBROUTINE interp_amma_time(day,day1,annee_ref & + & ,year_ini_amma,day_ini_amma,nt_amma,dt_amma,nlev_amma & + & ,vitw_amma,ht_amma,hq_amma,lat_amma,sens_amma & + & ,vitw_prof,ht_prof,hq_prof,lat_prof,sens_prof) + implicit none + +!--------------------------------------------------------------------------------------- +! Time interpolation of a 2D field to the timestep corresponding to day +! +! day: current julian day (e.g. 717538.2) +! day1: first day of the simulation +! nt_amma: total nb of data in the forcing (e.g. 48 for AMMA) +! dt_amma: total time interval (in sec) between 2 forcing data (e.g. 30min for AMMA) +!--------------------------------------------------------------------------------------- + +#include "compar1d.h" + +! inputs: + integer annee_ref + integer nt_amma,nlev_amma + integer year_ini_amma + real day, day1,day_ini_amma,dt_amma + real vitw_amma(nlev_amma,nt_amma) + real ht_amma(nlev_amma,nt_amma) + real hq_amma(nlev_amma,nt_amma) + real lat_amma(nt_amma) + real sens_amma(nt_amma) +! outputs: + real vitw_prof(nlev_amma) + real ht_prof(nlev_amma) + real hq_prof(nlev_amma) + real lat_prof,sens_prof +! local: + integer it_amma1, it_amma2,k + real timeit,time_amma1,time_amma2,frac + + + if (forcing_type.eq.6) then +! Check that initial day of the simulation consistent with AMMA case: + if (annee_ref.ne.2006) then + print*,'Pour AMMA, annee_ref doit etre 2006' + print*,'Changer annee_ref dans run.def' + stop + endif + if (annee_ref.eq.2006 .and. day1.lt.day_ini_amma) then + print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_amma + print*,'Changer dayref dans run.def' + stop + endif + if (annee_ref.eq.2006 .and. day1.gt.day_ini_amma+1) then + print*,'AMMA a fini le 11 juillet' + print*,'Changer dayref ou nday dans run.def' + stop + endif + endif + +! Determine timestep relative to the 1st day of AMMA: +! timeit=(day-day1)*86400. +! if (annee_ref.eq.1992) then +! timeit=(day-day_ini_toga)*86400. +! else +! timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992 +! endif + timeit=(day-day_ini_amma)*86400 + +! Determine the closest observation times: +! it_amma1=INT(timeit/dt_amma)+1 +! it_amma2=it_amma1 + 1 +! time_amma1=(it_amma1-1)*dt_amma +! time_amma2=(it_amma2-1)*dt_amma + + it_amma1=INT(timeit/dt_amma)+1 + IF (it_amma1 .EQ. nt_amma) THEN + it_amma2=it_amma1 + ELSE + it_amma2=it_amma1 + 1 + ENDIF + time_amma1=(it_amma1-1)*dt_amma + time_amma2=(it_amma2-1)*dt_amma + + if (it_amma1 .gt. nt_amma) then + write(*,*) 'PB-stop: day, it_amma1, it_amma2, timeit: ' & + & ,day,day_ini_amma,it_amma1,it_amma2,timeit/86400. + stop + endif + +! time interpolation: + frac=(time_amma2-timeit)/(time_amma2-time_amma1) + frac=max(frac,0.0) + + lat_prof = lat_amma(it_amma2) & + & -frac*(lat_amma(it_amma2)-lat_amma(it_amma1)) + sens_prof = sens_amma(it_amma2) & + & -frac*(sens_amma(it_amma2)-sens_amma(it_amma1)) + + do k=1,nlev_amma + vitw_prof(k) = vitw_amma(k,it_amma2) & + & -frac*(vitw_amma(k,it_amma2)-vitw_amma(k,it_amma1)) + ht_prof(k) = ht_amma(k,it_amma2) & + & -frac*(ht_amma(k,it_amma2)-ht_amma(k,it_amma1)) + hq_prof(k) = hq_amma(k,it_amma2) & + & -frac*(hq_amma(k,it_amma2)-hq_amma(k,it_amma1)) + enddo + + return + END + Index: LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90 =================================================================== --- LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90 (revision 2190) +++ LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90 (revision 2191) @@ -2,59 +2,79 @@ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -!Declarations specifiques au cas AMMA - character*80 :: fich_amma -! Option du cas AMMA ou on impose la discretisation verticale (Ap,Bp) - integer nlev_amma, nt_amma - - - integer year_ini_amma, day_ini_amma, mth_ini_amma - real heure_ini_amma - real day_ju_ini_amma ! Julian day of amma first day - parameter (year_ini_amma=2006) - parameter (mth_ini_amma=7) - parameter (day_ini_amma=10) ! 10 = 10Juil2006 - parameter (heure_ini_amma=0.) !0h en secondes - real dt_amma - parameter (dt_amma=1800.) - -!profils initiaux: - real, allocatable:: plev_amma(:) - - real, allocatable:: z_amma(:) - real, allocatable:: th_amma(:),q_amma(:) - real, allocatable:: u_amma(:) - real, allocatable:: v_amma(:) - - real, allocatable:: th_ammai(:),q_ammai(:) - real, allocatable:: u_ammai(:) - real, allocatable:: v_ammai(:) - real, allocatable:: vitw_ammai(:) - real, allocatable:: ht_ammai(:) - real, allocatable:: hq_ammai(:) - real, allocatable:: vt_ammai(:) - real, allocatable:: vq_ammai(:) - -!forcings - real, allocatable:: ht_amma(:,:) - real, allocatable:: hq_amma(:,:) - real, allocatable:: vitw_amma(:,:) - real, allocatable:: lat_amma(:),sens_amma(:) +!Declarations specifiques au cas standard + character*80 :: fich_cas +! Discrétisation + integer nlev_cas, nt_cas + + + integer year_ini_cas, day_ini_cas, mth_ini_cas + real heure_ini_cas + real day_ju_ini_cas ! Julian day of case first day + parameter (year_ini_cas=2011) + parameter (mth_ini_cas=10) + parameter (day_ini_cas=1) ! 10 = 10Juil2006 + parameter (heure_ini_cas=0.) !0h en secondes + real pdt_cas + parameter (pdt_cas=3.*3600) + +!CR ATTENTION TEST AMMA +! parameter (year_ini_cas=2006) +! parameter (mth_ini_cas=7) +! parameter (day_ini_cas=10) ! 10 = 10Juil2006 +! parameter (heure_ini_cas=0.) !0h en secondes +! parameter (pdt_cas=1800.) + +!profils environnementaux + real, allocatable:: plev_cas(:,:) + + real, allocatable:: z_cas(:,:) + real, allocatable:: t_cas(:,:),q_cas(:,:),rh_cas(:,:) + real, allocatable:: th_cas(:,:),rv_cas(:,:) + real, allocatable:: u_cas(:,:) + real, allocatable:: v_cas(:,:) + +!forcing + real, allocatable:: ht_cas(:,:),vt_cas(:,:),dt_cas(:,:),dtrad_cas(:,:) + real, allocatable:: hth_cas(:,:),vth_cas(:,:),dth_cas(:,:) + real, allocatable:: hq_cas(:,:),vq_cas(:,:),dq_cas(:,:) + real, allocatable:: hr_cas(:,:),vr_cas(:,:),dr_cas(:,:) + real, allocatable:: hu_cas(:,:),vu_cas(:,:),du_cas(:,:) + real, allocatable:: hv_cas(:,:),vv_cas(:,:),dv_cas(:,:) + real, allocatable:: vitw_cas(:,:) + real, allocatable:: ug_cas(:,:),vg_cas(:,:) + real, allocatable:: lat_cas(:),sens_cas(:),ts_cas(:) !champs interpoles - real, allocatable:: vitw_profamma(:) - real, allocatable:: ht_profamma(:) - real, allocatable:: hq_profamma(:) - real lat_profamma,sens_profamma - real, allocatable:: vt_profamma(:) - real, allocatable:: vq_profamma(:) - real, allocatable:: th_profamma(:) - real, allocatable:: q_profamma(:) - real, allocatable:: u_profamma(:) - real, allocatable:: v_profamma(:) + real, allocatable:: plev_prof_cas(:) + real, allocatable:: t_prof_cas(:) + real, allocatable:: q_prof_cas(:) + real, allocatable:: u_prof_cas(:) + real, allocatable:: v_prof_cas(:) + + real, allocatable:: vitw_prof_cas(:) + real, allocatable:: ug_prof_cas(:) + real, allocatable:: vg_prof_cas(:) + real, allocatable:: ht_prof_cas(:) + real, allocatable:: hq_prof_cas(:) + real, allocatable:: vt_prof_cas(:) + real, allocatable:: vq_prof_cas(:) + real, allocatable:: dt_prof_cas(:) + real, allocatable:: dtrad_prof_cas(:) + real, allocatable:: dq_prof_cas(:) + real, allocatable:: hu_prof_cas(:) + real, allocatable:: hv_prof_cas(:) + real, allocatable:: vu_prof_cas(:) + real, allocatable:: vv_prof_cas(:) + real, allocatable:: du_prof_cas(:) + real, allocatable:: dv_prof_cas(:) + + + real lat_prof_cas,sens_prof_cas,ts_prof_cas + CONTAINS -SUBROUTINE read_1D_cases +SUBROUTINE read_1D_cas implicit none @@ -62,9 +82,10 @@ INTEGER nid,rid,ierr - - fich_amma='amma.nc' - print*,'fich_amma ',fich_amma - ierr = NF_OPEN(fich_amma,NF_NOWRITE,nid) - print*,'fich_amma,NF_NOWRITE,nid ',fich_amma,NF_NOWRITE,nid + INTEGER ii,jj + + fich_cas='setup/cas.nc' + print*,'fich_cas ',fich_cas + ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid) + print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid if (ierr.NE.NF_NOERR) then write(*,*) 'ERROR: GROS Pb opening forcings nc file ' @@ -73,62 +94,89 @@ endif !....................................................................... + ierr=NF_INQ_DIMID(nid,'lat',rid) + IF (ierr.NE.NF_NOERR) THEN + print*, 'Oh probleme lecture dimension lat' + ENDIF + ierr=NF_INQ_DIMLEN(nid,rid,ii) + print*,'OK nid,rid,lat',nid,rid,ii +!....................................................................... + ierr=NF_INQ_DIMID(nid,'lon',rid) + IF (ierr.NE.NF_NOERR) THEN + print*, 'Oh probleme lecture dimension lon' + ENDIF + ierr=NF_INQ_DIMLEN(nid,rid,jj) + print*,'OK nid,rid,lat',nid,rid,jj +!....................................................................... ierr=NF_INQ_DIMID(nid,'lev',rid) IF (ierr.NE.NF_NOERR) THEN print*, 'Oh probleme lecture dimension zz' ENDIF - ierr=NF_INQ_DIMLEN(nid,rid,nlev_amma) - print*,'OK nid,rid,nlev_amma',nid,rid,nlev_amma + ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas) + print*,'OK nid,rid,nlev_cas',nid,rid,nlev_cas !....................................................................... ierr=NF_INQ_DIMID(nid,'time',rid) print*,'nid,rid',nid,rid - nt_amma=0 + nt_cas=0 IF (ierr.NE.NF_NOERR) THEN stop 'probleme lecture dimension sens' ENDIF - ierr=NF_INQ_DIMLEN(nid,rid,nt_amma) - print*,'nid,rid,nlev_amma',nid,rid,nt_amma + ierr=NF_INQ_DIMLEN(nid,rid,nt_cas) + print*,'nid,rid,nlev_cas',nid,rid,nt_cas !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! -!profils initiaux: - allocate(plev_amma(nlev_amma)) - - allocate(z_amma(nlev_amma)) - allocate(th_amma(nlev_amma),q_amma(nlev_amma)) - allocate(u_amma(nlev_amma)) - allocate(v_amma(nlev_amma)) - -!forcings - allocate(ht_amma(nlev_amma,nt_amma)) - allocate(hq_amma(nlev_amma,nt_amma)) - allocate(vitw_amma(nlev_amma,nt_amma)) - allocate(lat_amma(nt_amma),sens_amma(nt_amma)) - -!profils initiaux: - allocate(th_ammai(nlev_amma),q_ammai(nlev_amma)) - allocate(u_ammai(nlev_amma)) - allocate(v_ammai(nlev_amma)) - allocate(vitw_ammai(nlev_amma) ) - allocate(ht_ammai(nlev_amma)) - allocate(hq_ammai(nlev_amma)) - allocate(vt_ammai(nlev_amma)) - allocate(vq_ammai(nlev_amma)) +!profils moyens: + allocate(plev_cas(nlev_cas,nt_cas)) + allocate(z_cas(nlev_cas,nt_cas)) + allocate(t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas),rh_cas(nlev_cas,nt_cas)) + allocate(th_cas(nlev_cas,nt_cas),rv_cas(nlev_cas,nt_cas)) + allocate(u_cas(nlev_cas,nt_cas)) + allocate(v_cas(nlev_cas,nt_cas)) + +!forcing + 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)) + allocate(hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas),dq_cas(nlev_cas,nt_cas)) + allocate(hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas),dth_cas(nlev_cas,nt_cas)) + allocate(hr_cas(nlev_cas,nt_cas),vr_cas(nlev_cas,nt_cas),dr_cas(nlev_cas,nt_cas)) + allocate(hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas),du_cas(nlev_cas,nt_cas)) + allocate(hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas),dv_cas(nlev_cas,nt_cas)) + allocate(vitw_cas(nlev_cas,nt_cas)) + allocate(ug_cas(nlev_cas,nt_cas)) + allocate(vg_cas(nlev_cas,nt_cas)) + allocate(lat_cas(nt_cas),sens_cas(nt_cas),ts_cas(nt_cas)) + !champs interpoles - allocate(vitw_profamma(nlev_amma)) - allocate(ht_profamma(nlev_amma)) - allocate(hq_profamma(nlev_amma)) - allocate(vt_profamma(nlev_amma)) - allocate(vq_profamma(nlev_amma)) - allocate(th_profamma(nlev_amma)) - allocate(q_profamma(nlev_amma)) - allocate(u_profamma(nlev_amma)) - allocate(v_profamma(nlev_amma)) + allocate(plev_prof_cas(nlev_cas)) + allocate(t_prof_cas(nlev_cas)) + allocate(q_prof_cas(nlev_cas)) + allocate(u_prof_cas(nlev_cas)) + allocate(v_prof_cas(nlev_cas)) + + allocate(vitw_prof_cas(nlev_cas)) + allocate(ug_prof_cas(nlev_cas)) + allocate(vg_prof_cas(nlev_cas)) + allocate(ht_prof_cas(nlev_cas)) + allocate(hq_prof_cas(nlev_cas)) + allocate(hu_prof_cas(nlev_cas)) + allocate(hv_prof_cas(nlev_cas)) + allocate(vt_prof_cas(nlev_cas)) + allocate(vq_prof_cas(nlev_cas)) + allocate(vu_prof_cas(nlev_cas)) + allocate(vv_prof_cas(nlev_cas)) + allocate(dt_prof_cas(nlev_cas)) + allocate(dtrad_prof_cas(nlev_cas)) + allocate(dq_prof_cas(nlev_cas)) + allocate(du_prof_cas(nlev_cas)) + allocate(dv_prof_cas(nlev_cas)) print*,'Allocations OK' - call read_amma(nid,nlev_amma,nt_amma & - & ,z_amma,plev_amma,th_amma,q_amma,u_amma,v_amma,vitw_amma & - & ,ht_amma,hq_amma,sens_amma,lat_amma) - -END SUBROUTINE read_1D_cases + call read_cas(nid,nlev_cas,nt_cas & + & ,z_cas,plev_cas,t_cas,q_cas,rh_cas,th_cas,rv_cas,u_cas,v_cas & + & ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas,dv_cas,hv_cas,vv_cas & + & ,dt_cas,dtrad_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas & + & ,dth_cas,hth_cas,vth_cas,dr_cas,hr_cas,vr_cas,sens_cas,lat_cas,ts_cas) + + +END SUBROUTINE read_1D_cas @@ -136,37 +184,53 @@ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! SUBROUTINE deallocate_1D_cases -!profils initiaux: - deallocate(plev_amma) +!profils environnementaux: + deallocate(plev_cas) - deallocate(z_amma) - deallocate(th_amma,q_amma) - deallocate(u_amma) - deallocate(v_amma) - - deallocate(th_ammai,q_ammai) - deallocate(u_ammai) - deallocate(v_ammai) - deallocate(vitw_ammai ) - deallocate(ht_ammai) - deallocate(hq_ammai) - deallocate(vt_ammai) - deallocate(vq_ammai) + deallocate(z_cas) + deallocate(t_cas,q_cas,rh_cas) + deallocate(th_cas,rv_cas) + deallocate(u_cas) + deallocate(v_cas) -!forcings - deallocate(ht_amma) - deallocate(hq_amma) - deallocate(vitw_amma) - deallocate(lat_amma,sens_amma) +!forcing + deallocate(ht_cas,vt_cas,dt_cas,dtrad_cas) + deallocate(hq_cas,vq_cas,dq_cas) + deallocate(hth_cas,vth_cas,dth_cas) + deallocate(hr_cas,vr_cas,dr_cas) + deallocate(hu_cas,vu_cas,du_cas) + deallocate(hv_cas,vv_cas,dv_cas) + deallocate(vitw_cas) + deallocate(ug_cas) + deallocate(vg_cas) + deallocate(lat_cas,sens_cas,ts_cas) !champs interpoles - deallocate(vitw_profamma) - deallocate(ht_profamma) - deallocate(hq_profamma) - deallocate(vt_profamma) - deallocate(vq_profamma) - deallocate(th_profamma) - deallocate(q_profamma) - deallocate(u_profamma) - deallocate(v_profamma) + deallocate(plev_prof_cas) + deallocate(t_prof_cas) + deallocate(q_prof_cas) + deallocate(u_prof_cas) + deallocate(v_prof_cas) + + deallocate(vitw_prof_cas) + deallocate(ug_prof_cas) + deallocate(vg_prof_cas) + deallocate(ht_prof_cas) + deallocate(hq_prof_cas) + deallocate(hu_prof_cas) + deallocate(hv_prof_cas) + deallocate(vt_prof_cas) + deallocate(vq_prof_cas) + deallocate(vu_prof_cas) + deallocate(vv_prof_cas) + deallocate(dt_prof_cas) + deallocate(dtrad_prof_cas) + deallocate(dq_prof_cas) + deallocate(du_prof_cas) + deallocate(dv_prof_cas) + deallocate(t_prof_cas) + deallocate(q_prof_cas) + deallocate(u_prof_cas) + deallocate(v_prof_cas) + END SUBROUTINE deallocate_1D_cases @@ -174,9 +238,10 @@ END MODULE mod_1D_cases_read !===================================================================== - subroutine read_amma(nid,nlevel,ntime & - & ,zz,pp,temp,qv,u,v,dw & - & ,dt,dq,sens,flat) - -!program reading forcings of the AMMA case study + subroutine read_cas(nid,nlevel,ntime & + & ,zz,pp,temp,qv,rh,theta,rv,u,v,ug,vg,w, & + & du,hu,vu,dv,hv,vv,dt,dtrad,ht,vt,dq,hq,vq, & + & dth,hth,vth,dr,hr,vr,sens,flat,ts) + +!program reading forcing of the case study implicit none #include "netcdf.inc" @@ -184,17 +249,26 @@ integer ntime,nlevel - real zz(nlevel) - real temp(nlevel),pp(nlevel) - real qv(nlevel),u(nlevel) - real v(nlevel) - real dw(nlevel,ntime) - real dt(nlevel,ntime) - real dq(nlevel,ntime) - real flat(ntime),sens(ntime) + real zz(nlevel,ntime) + real pp(nlevel,ntime) + real temp(nlevel,ntime),qv(nlevel,ntime),rh(nlevel,ntime) + real theta(nlevel,ntime),rv(nlevel,ntime) + real u(nlevel,ntime) + real v(nlevel,ntime) + real ug(nlevel,ntime) + real vg(nlevel,ntime) + real w(nlevel,ntime) + real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime) + real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime) + real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime) + real dtrad(nlevel,ntime) + real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime) + real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime) + real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime) + real flat(ntime),sens(ntime),ts(ntime) integer nid, ierr,rid integer nbvar3d - parameter(nbvar3d=30) + parameter(nbvar3d=34) integer var3didin(nbvar3d) @@ -204,7 +278,13 @@ stop 'lev' endif - - - ierr=NF_INQ_VARID(nid,"temp",var3didin(2)) + + ierr=NF_INQ_VARID(nid,"pp",var3didin(2)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'plev' + endif + + + ierr=NF_INQ_VARID(nid,"temp",var3didin(3)) if(ierr/=NF_NOERR) then write(*,*) NF_STRERROR(ierr) @@ -212,5 +292,5 @@ endif - ierr=NF_INQ_VARID(nid,"qv",var3didin(3)) + ierr=NF_INQ_VARID(nid,"qv",var3didin(4)) if(ierr/=NF_NOERR) then write(*,*) NF_STRERROR(ierr) @@ -218,5 +298,24 @@ endif - ierr=NF_INQ_VARID(nid,"u",var3didin(4)) + ierr=NF_INQ_VARID(nid,"rh",var3didin(5)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'rh' + endif + + ierr=NF_INQ_VARID(nid,"theta",var3didin(6)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'theta' + endif + + ierr=NF_INQ_VARID(nid,"rv",var3didin(7)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'rv' + endif + + + ierr=NF_INQ_VARID(nid,"u",var3didin(8)) if(ierr/=NF_NOERR) then write(*,*) NF_STRERROR(ierr) @@ -224,5 +323,5 @@ endif - ierr=NF_INQ_VARID(nid,"v",var3didin(5)) + ierr=NF_INQ_VARID(nid,"v",var3didin(9)) if(ierr/=NF_NOERR) then write(*,*) NF_STRERROR(ierr) @@ -230,23 +329,137 @@ endif - ierr=NF_INQ_VARID(nid,"dw",var3didin(6)) - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - stop 'dw' - endif - - ierr=NF_INQ_VARID(nid,"dt",var3didin(7)) - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - stop 'dt' - endif - - ierr=NF_INQ_VARID(nid,"dq",var3didin(8)) - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - stop 'dq' + ierr=NF_INQ_VARID(nid,"ug",var3didin(10)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'ug' + endif + + ierr=NF_INQ_VARID(nid,"vg",var3didin(11)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vg' + endif + + ierr=NF_INQ_VARID(nid,"w",var3didin(12)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'w' + endif + + ierr=NF_INQ_VARID(nid,"advu",var3didin(13)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'advu' + endif + + ierr=NF_INQ_VARID(nid,"hu",var3didin(14)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'hu' + endif + + ierr=NF_INQ_VARID(nid,"vu",var3didin(15)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vu' + endif + + ierr=NF_INQ_VARID(nid,"advv",var3didin(16)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'advv' + endif + + ierr=NF_INQ_VARID(nid,"hv",var3didin(17)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'hv' + endif + + ierr=NF_INQ_VARID(nid,"vv",var3didin(18)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vv' + endif + + ierr=NF_INQ_VARID(nid,"advT",var3didin(19)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'advT' + endif + + ierr=NF_INQ_VARID(nid,"hT",var3didin(20)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'hT' + endif + + ierr=NF_INQ_VARID(nid,"vT",var3didin(21)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vT' + endif + + ierr=NF_INQ_VARID(nid,"advq",var3didin(22)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'advq' endif - ierr=NF_INQ_VARID(nid,"sens",var3didin(9)) + ierr=NF_INQ_VARID(nid,"hq",var3didin(23)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'hq' + endif + + ierr=NF_INQ_VARID(nid,"vq",var3didin(24)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vq' + endif + + ierr=NF_INQ_VARID(nid,"advth",var3didin(25)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'advth' + endif + + ierr=NF_INQ_VARID(nid,"hth",var3didin(26)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'hth' + endif + + ierr=NF_INQ_VARID(nid,"vth",var3didin(27)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vth' + endif + + ierr=NF_INQ_VARID(nid,"advr",var3didin(28)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'advr' + endif + + ierr=NF_INQ_VARID(nid,"hr",var3didin(29)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'hr' + endif + + ierr=NF_INQ_VARID(nid,"vr",var3didin(30)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'vr' + endif + + ierr=NF_INQ_VARID(nid,"radT",var3didin(31)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'radT' + endif + + ierr=NF_INQ_VARID(nid,"sens",var3didin(32)) if(ierr/=NF_NOERR) then write(*,*) NF_STRERROR(ierr) @@ -254,5 +467,5 @@ endif - ierr=NF_INQ_VARID(nid,"flat",var3didin(10)) + ierr=NF_INQ_VARID(nid,"flat",var3didin(33)) if(ierr/=NF_NOERR) then write(*,*) NF_STRERROR(ierr) @@ -260,11 +473,9 @@ endif - ierr=NF_INQ_VARID(nid,"pp",var3didin(11)) - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - endif - -!dimensions lecture -! call catchaxis(nid,ntime,nlevel,time,z,ierr) + ierr=NF_INQ_VARID(nid,"ts",var3didin(34)) + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop 'ts' + endif #ifdef NC_DOUBLE @@ -280,18 +491,30 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(2),temp) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(2),temp) -#endif - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - stop "getvarup" - endif -! write(*,*)'lecture th ok',temp - -#ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(3),qv) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(3),qv) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(2),pp) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(2),pp) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture pp ok',pp + + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(3),temp) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(3),temp) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture T ok',temp + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(4),qv) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(4),qv) #endif if(ierr/=NF_NOERR) then @@ -302,7 +525,40 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(4),u) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(4),u) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(5),rh) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(5),rh) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture rh ok',rh + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(6),theta) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(6),theta) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture theta ok',theta + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(7),rv) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(7),rv) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture rv ok',rv + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(8),u) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(8),u) #endif if(ierr/=NF_NOERR) then @@ -313,7 +569,7 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(5),v) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(5),v) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(9),v) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(9),v) #endif if(ierr/=NF_NOERR) then @@ -324,18 +580,106 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(6),dw) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(6),dw) -#endif - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - stop "getvarup" - endif -! write(*,*)'lecture w ok',dw - -#ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(7),dt) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(7),dt) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(10),ug) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(10),ug) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture ug ok',ug + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(11),vg) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(11),vg) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vg ok',vg + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(12),w) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(12),w) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture w ok',w + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(13),du) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(13),du) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture du ok',du + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(14),hu) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(14),hu) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture hu ok',hu + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(15),vu) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(15),vu) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vu ok',vu + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(16),dv) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(16),dv) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture dv ok',dv + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(17),hv) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(17),hv) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture hv ok',hv + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(18),vv) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(18),vv) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vv ok',vv + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(19),dt) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(19),dt) #endif if(ierr/=NF_NOERR) then @@ -346,7 +690,29 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(8),dq) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(8),dq) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(20),ht) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(20),ht) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture ht ok',ht + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(21),vt) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(21),vt) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vt ok',vt + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(22),dq) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(22),dq) #endif if(ierr/=NF_NOERR) then @@ -357,7 +723,106 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(9),sens) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(9),sens) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(23),hq) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(23),hq) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture hq ok',hq + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(24),vq) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(24),vq) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vq ok',vq + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(25),dth) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(25),dth) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture dth ok',dth + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(26),hth) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(26),hth) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture hth ok',hth + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(27),vth) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(27),vth) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vth ok',vth + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(28),dr) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(28),dr) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture dr ok',dr + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(29),hr) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(29),hr) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture hr ok',hr + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(30),vr) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(30),vr) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture vr ok',vr + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(31),dtrad) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(31),dtrad) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture dtrad ok',dtrad + +#ifdef NC_DOUBLE + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(32),sens) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(32),sens) #endif if(ierr/=NF_NOERR) then @@ -368,7 +833,7 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(10),flat) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(10),flat) + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(33),flat) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(33),flat) #endif if(ierr/=NF_NOERR) then @@ -379,21 +844,31 @@ #ifdef NC_DOUBLE - ierr = NF_GET_VAR_DOUBLE(nid,var3didin(11),pp) -#else - ierr = NF_GET_VAR_REAL(nid,var3didin(11),pp) -#endif - if(ierr/=NF_NOERR) then - write(*,*) NF_STRERROR(ierr) - stop "getvarup" - endif -! write(*,*)'lecture pp ok',pp + ierr = NF_GET_VAR_DOUBLE(nid,var3didin(34),ts) +#else + ierr = NF_GET_VAR_REAL(nid,var3didin(34),ts) +#endif + if(ierr/=NF_NOERR) then + write(*,*) NF_STRERROR(ierr) + stop "getvarup" + endif +! write(*,*)'lecture ts ok',ts return - end subroutine read_amma + end subroutine read_cas !====================================================================== - SUBROUTINE interp_amma_time(day,day1,annee_ref & - & ,year_ini_amma,day_ini_amma,nt_amma,dt_amma,nlev_amma & - & ,vitw_amma,ht_amma,hq_amma,lat_amma,sens_amma & - & ,vitw_prof,ht_prof,hq_prof,lat_prof,sens_prof) + SUBROUTINE interp_case_time(day,day1,annee_ref & + & ,year_ini_cas,day_ini_cas,nt_cas,pdt_cas,nlev_cas & + & ,ts_cas,plev_cas,t_cas,q_cas,u_cas,v_cas & + & ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas & + & ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas & + & ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas & + & ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas & + & ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas & + & ,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas & + & ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas & + & ,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas & + & ,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas) + + implicit none @@ -403,6 +878,6 @@ ! day: current julian day (e.g. 717538.2) ! day1: first day of the simulation -! nt_amma: total nb of data in the forcing (e.g. 48 for AMMA) -! dt_amma: total time interval (in sec) between 2 forcing data (e.g. 30min for AMMA) +! nt_cas: total nb of data in the forcing +! pdt_cas: total time interval (in sec) between 2 forcing data !--------------------------------------------------------------------------------------- @@ -411,23 +886,61 @@ ! inputs: integer annee_ref - integer nt_amma,nlev_amma - integer year_ini_amma - real day, day1,day_ini_amma,dt_amma - real vitw_amma(nlev_amma,nt_amma) - real ht_amma(nlev_amma,nt_amma) - real hq_amma(nlev_amma,nt_amma) - real lat_amma(nt_amma) - real sens_amma(nt_amma) + integer nt_cas,nlev_cas + integer year_ini_cas + real day_ini_cas + real day, day1,pdt_cas + real ts_cas(nt_cas) + real plev_cas(nlev_cas,nt_cas) + real t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas) + real u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas) + real ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas) + real vitw_cas(nlev_cas,nt_cas) + real du_cas(nlev_cas,nt_cas),hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas) + real dv_cas(nlev_cas,nt_cas),hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas) + real dt_cas(nlev_cas,nt_cas),ht_cas(nlev_cas,nt_cas),vt_cas(nlev_cas,nt_cas) + real dtrad_cas(nlev_cas,nt_cas) + real dq_cas(nlev_cas,nt_cas),hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas) + real lat_cas(nt_cas) + real sens_cas(nt_cas) + ! outputs: - real vitw_prof(nlev_amma) - real ht_prof(nlev_amma) - real hq_prof(nlev_amma) - real lat_prof,sens_prof + real plev_prof_cas(nlev_cas) + real t_prof_cas(nlev_cas),q_prof_cas(nlev_cas) + real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas) + real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas) + real vitw_prof_cas(nlev_cas) + real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas) + real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas) + real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas) + real dtrad_prof_cas(nlev_cas) + real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas) + real lat_prof_cas,sens_prof_cas,ts_prof_cas ! local: - integer it_amma1, it_amma2,k - real timeit,time_amma1,time_amma2,frac - - - if (forcing_type.eq.6) then + integer it_cas1, it_cas2,k + real timeit,time_cas1,time_cas2,frac + + + print*,'Check time',day1,day_ini_cas,day_ini_cas+90 + + if ((forcing_type.eq.10).and.(1.eq.1)) then +! Check that initial day of the simulation consistent with the case: + if (annee_ref.ne.2011) then + print*,'Pour CINDY, annee_ref doit etre 2011' + print*,'Changer annee_ref dans run.def' + stop + endif + if (annee_ref.eq.2011 .and. day1.lt.day_ini_cas) then + print*,'CINDY a débuté le 1 octobre 2011',day1,day_ini_cas + print*,'Changer dayref dans run.def' + stop + endif + if (annee_ref.eq.2011 .and. day1.gt.day_ini_cas+90) then + print*,'CINDY a fini le 31 decembre' + print*,'Changer dayref ou nday dans run.def',day1,day_ini_cas+90 + stop + endif + endif +!CR:ATTENTION TEST AMMA + if ((forcing_type.eq.10).and.(1.eq.0)) then ! Check that initial day of the simulation consistent with AMMA case: if (annee_ref.ne.2006) then @@ -436,10 +949,10 @@ stop endif - if (annee_ref.eq.2006 .and. day1.lt.day_ini_amma) then - print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_amma + if (annee_ref.eq.2006 .and. day1.lt.day_ini_cas) then + print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_cas print*,'Changer dayref dans run.def' stop endif - if (annee_ref.eq.2006 .and. day1.gt.day_ini_amma+1) then + if (annee_ref.eq.2006 .and. day1.gt.day_ini_cas+1) then print*,'AMMA a fini le 11 juillet' print*,'Changer dayref ou nday dans run.def' @@ -448,50 +961,88 @@ endif -! Determine timestep relative to the 1st day of AMMA: +! Determine timestep relative to the 1st day: ! timeit=(day-day1)*86400. ! if (annee_ref.eq.1992) then -! timeit=(day-day_ini_toga)*86400. +! timeit=(day-day_ini_cas)*86400. ! else ! timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992 ! endif - timeit=(day-day_ini_amma)*86400 + timeit=(day-day_ini_cas)*86400 ! Determine the closest observation times: -! it_amma1=INT(timeit/dt_amma)+1 -! it_amma2=it_amma1 + 1 -! time_amma1=(it_amma1-1)*dt_amma -! time_amma2=(it_amma2-1)*dt_amma - - it_amma1=INT(timeit/dt_amma)+1 - IF (it_amma1 .EQ. nt_amma) THEN - it_amma2=it_amma1 +! it_cas1=INT(timeit/pdt_cas)+1 +! it_cas2=it_cas1 + 1 +! time_cas1=(it_cas1-1)*pdt_cas +! time_cas2=(it_cas2-1)*pdt_cas + + it_cas1=INT(timeit/pdt_cas)+1 + IF (it_cas1 .EQ. nt_cas) THEN + it_cas2=it_cas1 ELSE - it_amma2=it_amma1 + 1 + it_cas2=it_cas1 + 1 ENDIF - time_amma1=(it_amma1-1)*dt_amma - time_amma2=(it_amma2-1)*dt_amma - - if (it_amma1 .gt. nt_amma) then - write(*,*) 'PB-stop: day, it_amma1, it_amma2, timeit: ' & - & ,day,day_ini_amma,it_amma1,it_amma2,timeit/86400. + time_cas1=(it_cas1-1)*pdt_cas + time_cas2=(it_cas2-1)*pdt_cas + + if (it_cas1 .gt. nt_cas) then + write(*,*) 'PB-stop: day, it_cas1, it_cas2, timeit: ' & + & ,day,day_ini_cas,it_cas1,it_cas2,timeit/86400. stop endif ! time interpolation: - frac=(time_amma2-timeit)/(time_amma2-time_amma1) + frac=(time_cas2-timeit)/(time_cas2-time_cas1) frac=max(frac,0.0) - lat_prof = lat_amma(it_amma2) & - & -frac*(lat_amma(it_amma2)-lat_amma(it_amma1)) - sens_prof = sens_amma(it_amma2) & - & -frac*(sens_amma(it_amma2)-sens_amma(it_amma1)) - - do k=1,nlev_amma - vitw_prof(k) = vitw_amma(k,it_amma2) & - & -frac*(vitw_amma(k,it_amma2)-vitw_amma(k,it_amma1)) - ht_prof(k) = ht_amma(k,it_amma2) & - & -frac*(ht_amma(k,it_amma2)-ht_amma(k,it_amma1)) - hq_prof(k) = hq_amma(k,it_amma2) & - & -frac*(hq_amma(k,it_amma2)-hq_amma(k,it_amma1)) + lat_prof_cas = lat_cas(it_cas2) & + & -frac*(lat_cas(it_cas2)-lat_cas(it_cas1)) + sens_prof_cas = sens_cas(it_cas2) & + & -frac*(sens_cas(it_cas2)-sens_cas(it_cas1)) + ts_prof_cas = ts_cas(it_cas2) & + & -frac*(ts_cas(it_cas2)-ts_cas(it_cas1)) + + do k=1,nlev_cas + plev_prof_cas(k) = plev_cas(k,it_cas2) & + & -frac*(plev_cas(k,it_cas2)-plev_cas(k,it_cas1)) + t_prof_cas(k) = t_cas(k,it_cas2) & + & -frac*(t_cas(k,it_cas2)-t_cas(k,it_cas1)) + q_prof_cas(k) = q_cas(k,it_cas2) & + & -frac*(q_cas(k,it_cas2)-q_cas(k,it_cas1)) + u_prof_cas(k) = u_cas(k,it_cas2) & + & -frac*(u_cas(k,it_cas2)-u_cas(k,it_cas1)) + v_prof_cas(k) = v_cas(k,it_cas2) & + & -frac*(v_cas(k,it_cas2)-v_cas(k,it_cas1)) + ug_prof_cas(k) = ug_cas(k,it_cas2) & + & -frac*(ug_cas(k,it_cas2)-ug_cas(k,it_cas1)) + vg_prof_cas(k) = vg_cas(k,it_cas2) & + & -frac*(vg_cas(k,it_cas2)-vg_cas(k,it_cas1)) + vitw_prof_cas(k) = vitw_cas(k,it_cas2) & + & -frac*(vitw_cas(k,it_cas2)-vitw_cas(k,it_cas1)) + du_prof_cas(k) = du_cas(k,it_cas2) & + & -frac*(du_cas(k,it_cas2)-du_cas(k,it_cas1)) + hu_prof_cas(k) = hu_cas(k,it_cas2) & + & -frac*(hu_cas(k,it_cas2)-hu_cas(k,it_cas1)) + vu_prof_cas(k) = vu_cas(k,it_cas2) & + & -frac*(vu_cas(k,it_cas2)-vu_cas(k,it_cas1)) + dv_prof_cas(k) = dv_cas(k,it_cas2) & + & -frac*(dv_cas(k,it_cas2)-dv_cas(k,it_cas1)) + hv_prof_cas(k) = hv_cas(k,it_cas2) & + & -frac*(hv_cas(k,it_cas2)-hv_cas(k,it_cas1)) + vv_prof_cas(k) = vv_cas(k,it_cas2) & + & -frac*(vv_cas(k,it_cas2)-vv_cas(k,it_cas1)) + dt_prof_cas(k) = dt_cas(k,it_cas2) & + & -frac*(dt_cas(k,it_cas2)-dt_cas(k,it_cas1)) + ht_prof_cas(k) = ht_cas(k,it_cas2) & + & -frac*(ht_cas(k,it_cas2)-ht_cas(k,it_cas1)) + vt_prof_cas(k) = vt_cas(k,it_cas2) & + & -frac*(vt_cas(k,it_cas2)-vt_cas(k,it_cas1)) + dtrad_prof_cas(k) = dtrad_cas(k,it_cas2) & + & -frac*(dtrad_cas(k,it_cas2)-dtrad_cas(k,it_cas1)) + dq_prof_cas(k) = dq_cas(k,it_cas2) & + & -frac*(dq_cas(k,it_cas2)-dq_cas(k,it_cas1)) + hq_prof_cas(k) = hq_cas(k,it_cas2) & + & -frac*(hq_cas(k,it_cas2)-hq_cas(k,it_cas1)) + vq_prof_cas(k) = vq_cas(k,it_cas2) & + & -frac*(vq_cas(k,it_cas2)-vq_cas(k,it_cas1)) enddo @@ -499,2 +1050,3 @@ END +!**********************************************************************************************