Changeset 2191 for LMDZ5/trunk/libf


Ignore:
Timestamp:
Feb 3, 2015, 11:00:57 AM (10 years ago)
Author:
fhourdin
Message:

Ajout du cas 1D CINDY-DYNAMO, utilisant le nouveau format standard, amené à être étendu aux autres cas.
Addition of the CINDY-DYNAMO 1D case, using the new standard format for 1D cases, that will be extended to all 1D cases.
Catherine Rio

Location:
LMDZ5/trunk/libf/phylmd
Files:
1 added
7 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/phylmd/1DUTILS.h

    r2181 r2191  
    9999!             LS convergence imposed from  RICO (cst)
    100100!         = 6 ==> forcing_amma = .true.
     101!         = 10 ==> forcing_case = .true.
     102!             initial profiles from case.nc file
    101103!         = 40 ==> forcing_GCSSold = .true.
    102104!             initial profile from GCSS file
    … …  
    132134          CALL getin('turb_fcg',xTurb_fcg_gcssold)
    133135        ENDIF
     136
     137!Paramètres de forçage
     138!Config  Key  = tend_t
     139!Config  Desc = forcage ou non par advection de T
     140!Config  Def  = false
     141!Config  Help = forcage ou non par advection de T
     142       tend_t =0
     143       CALL getin('tend_t',tend_t)
     144
     145!Config  Key  = tend_q
     146!Config  Desc = forcage ou non par advection de q
     147!Config  Def  = false
     148!Config  Help = forcage ou non par advection de q
     149       tend_q =0
     150       CALL getin('tend_q',tend_q)
     151
     152!Config  Key  = tend_u
     153!Config  Desc = forcage ou non par advection de u
     154!Config  Def  = false
     155!Config  Help = forcage ou non par advection de u
     156       tend_u =0
     157       CALL getin('tend_u',tend_u)
     158
     159!Config  Key  = tend_v
     160!Config  Desc = forcage ou non par advection de v
     161!Config  Def  = false
     162!Config  Help = forcage ou non par advection de v
     163       tend_v =0
     164       CALL getin('tend_v',tend_v)
     165
     166!Config  Key  = tend_w
     167!Config  Desc = forcage ou non par vitesse verticale
     168!Config  Def  = false
     169!Config  Help = forcage ou non par vitesse verticale
     170       tend_w =0
     171       CALL getin('tend_w',tend_w)
     172
     173!Config  Key  = tend_rayo
     174!Config  Desc = forcage ou non par dtrad
     175!Config  Def  = false
     176!Config  Help = forcage ou non par dtrad
     177       tend_rayo =0
     178       CALL getin('tend_rayo',tend_rayo)
     179
     180
     181!Config  Key  = nudge_t
     182!Config  Desc = constante de nudging de T
     183!Config  Def  = false
     184!Config  Help = constante de nudging de T
     185       nudge_t =0.
     186       CALL getin('nudge_t',nudge_t)
     187
     188!Config  Key  = nudge_q
     189!Config  Desc = constante de nudging de q
     190!Config  Def  = false
     191!Config  Help = constante de nudging de q
     192       nudge_q =0.
     193       CALL getin('nudge_q',nudge_q)
     194
     195!Config  Key  = nudge_u
     196!Config  Desc = constante de nudging de u
     197!Config  Def  = false
     198!Config  Help = constante de nudging de u
     199       nudge_u =0.
     200       CALL getin('nudge_u',nudge_u)
     201
     202!Config  Key  = nudge_v
     203!Config  Desc = constante de nudging de v
     204!Config  Def  = false
     205!Config  Help = constante de nudging de v
     206       nudge_v =0.
     207       CALL getin('nudge_v',nudge_v)
     208
     209!Config  Key  = nudge_w
     210!Config  Desc = constante de nudging de w
     211!Config  Def  = false
     212!Config  Help = constante de nudging de w
     213       nudge_w =0.
     214       CALL getin('nudge_w',nudge_w)
     215
    134216
    135217!Config  Key  = iflag_nudge
    … …  
    24312513 
    24322514!=====================================================================
     2515       SUBROUTINE interp_case_vertical(play,nlev_cas,plev_prof_cas            &
     2516     &         ,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas                         &
     2517     &         ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas           &
     2518     &         ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas                            &
     2519     &         ,t_mod_cas,q_mod_cas,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas                              &
     2520     &         ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas               &
     2521     &         ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas,mxcalc)
     2522 
     2523       implicit none
     2524 
     2525#include "dimensions.h"
     2526
     2527!-------------------------------------------------------------------------
     2528! Vertical interpolation of TOGA-COARE forcing data onto mod_casel levels
     2529!-------------------------------------------------------------------------
     2530 
     2531       integer nlevmax
     2532       parameter (nlevmax=41)
     2533       integer nlev_cas,mxcalc
     2534!       real play(llm), plev_prof(nlevmax)
     2535!       real t_prof(nlevmax),q_prof(nlevmax)
     2536!       real u_prof(nlevmax),v_prof(nlevmax), w_prof(nlevmax)
     2537!       real ht_prof(nlevmax),vt_prof(nlevmax)
     2538!       real hq_prof(nlevmax),vq_prof(nlevmax)
     2539 
     2540       real play(llm), plev_prof_cas(nlev_cas)
     2541       real t_prof_cas(nlev_cas),q_prof_cas(nlev_cas)
     2542       real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas)
     2543       real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas), vitw_prof_cas(nlev_cas)
     2544       real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas)
     2545       real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas)
     2546       real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas),dtrad_prof_cas(nlev_cas)
     2547       real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas)
     2548 
     2549       real t_mod_cas(llm),q_mod_cas(llm)
     2550       real u_mod_cas(llm),v_mod_cas(llm)
     2551       real ug_mod_cas(llm),vg_mod_cas(llm), w_mod_cas(llm)
     2552       real du_mod_cas(llm),hu_mod_cas(llm),vu_mod_cas(llm)
     2553       real dv_mod_cas(llm),hv_mod_cas(llm),vv_mod_cas(llm)
     2554       real dt_mod_cas(llm),ht_mod_cas(llm),vt_mod_cas(llm),dtrad_mod_cas(llm)
     2555       real dq_mod_cas(llm),hq_mod_cas(llm),vq_mod_cas(llm)
     2556 
     2557       integer l,k,k1,k2
     2558       real frac,frac1,frac2,fact
     2559 
     2560       do l = 1, llm
     2561
     2562        if (play(l).ge.plev_prof_cas(nlev_cas)) then
     2563 
     2564        mxcalc=l
     2565         k1=0
     2566         k2=0
     2567
     2568         if (play(l).le.plev_prof_cas(1)) then
     2569
     2570         do k = 1, nlev_cas-1
     2571          if (play(l).le.plev_prof_cas(k).and. play(l).gt.plev_prof_cas(k+1)) then
     2572            k1=k
     2573            k2=k+1
     2574          endif
     2575         enddo
     2576
     2577         if (k1.eq.0 .or. k2.eq.0) then
     2578          write(*,*) 'PB! k1, k2 = ',k1,k2
     2579          write(*,*) 'l,play(l) = ',l,play(l)/100
     2580         do k = 1, nlev_cas-1
     2581          write(*,*) 'k,plev_prof_cas(k) = ',k,plev_prof_cas(k)/100
     2582         enddo
     2583         endif
     2584
     2585         frac = (plev_prof_cas(k2)-play(l))/(plev_prof_cas(k2)-plev_prof_cas(k1))
     2586         t_mod_cas(l)= t_prof_cas(k2) - frac*(t_prof_cas(k2)-t_prof_cas(k1))
     2587         q_mod_cas(l)= q_prof_cas(k2) - frac*(q_prof_cas(k2)-q_prof_cas(k1))
     2588         u_mod_cas(l)= u_prof_cas(k2) - frac*(u_prof_cas(k2)-u_prof_cas(k1))
     2589         v_mod_cas(l)= v_prof_cas(k2) - frac*(v_prof_cas(k2)-v_prof_cas(k1))
     2590         ug_mod_cas(l)= ug_prof_cas(k2) - frac*(ug_prof_cas(k2)-ug_prof_cas(k1))
     2591         vg_mod_cas(l)= vg_prof_cas(k2) - frac*(vg_prof_cas(k2)-vg_prof_cas(k1))
     2592         w_mod_cas(l)= vitw_prof_cas(k2) - frac*(vitw_prof_cas(k2)-vitw_prof_cas(k1))
     2593         du_mod_cas(l)= du_prof_cas(k2) - frac*(du_prof_cas(k2)-du_prof_cas(k1))
     2594         hu_mod_cas(l)= hu_prof_cas(k2) - frac*(hu_prof_cas(k2)-hu_prof_cas(k1))
     2595         vu_mod_cas(l)= vu_prof_cas(k2) - frac*(vu_prof_cas(k2)-vu_prof_cas(k1))
     2596         dv_mod_cas(l)= dv_prof_cas(k2) - frac*(dv_prof_cas(k2)-dv_prof_cas(k1))
     2597         hv_mod_cas(l)= hv_prof_cas(k2) - frac*(hv_prof_cas(k2)-hv_prof_cas(k1))
     2598         vv_mod_cas(l)= vv_prof_cas(k2) - frac*(vv_prof_cas(k2)-vv_prof_cas(k1))
     2599         dt_mod_cas(l)= dt_prof_cas(k2) - frac*(dt_prof_cas(k2)-dt_prof_cas(k1))
     2600         ht_mod_cas(l)= ht_prof_cas(k2) - frac*(ht_prof_cas(k2)-ht_prof_cas(k1))
     2601         vt_mod_cas(l)= vt_prof_cas(k2) - frac*(vt_prof_cas(k2)-vt_prof_cas(k1))
     2602         dq_mod_cas(l)= dq_prof_cas(k2) - frac*(dq_prof_cas(k2)-dq_prof_cas(k1))
     2603         hq_mod_cas(l)= hq_prof_cas(k2) - frac*(hq_prof_cas(k2)-hq_prof_cas(k1))
     2604         vq_mod_cas(l)= vq_prof_cas(k2) - frac*(vq_prof_cas(k2)-vq_prof_cas(k1))
     2605     
     2606         else !play>plev_prof_cas(1)
     2607
     2608         k1=1
     2609         k2=2
     2610         frac1 = (play(l)-plev_prof_cas(k2))/(plev_prof_cas(k1)-plev_prof_cas(k2))
     2611         frac2 = (play(l)-plev_prof_cas(k1))/(plev_prof_cas(k1)-plev_prof_cas(k2))
     2612         t_mod_cas(l)= frac1*t_prof_cas(k1) - frac2*t_prof_cas(k2)
     2613         q_mod_cas(l)= frac1*q_prof_cas(k1) - frac2*q_prof_cas(k2)
     2614         u_mod_cas(l)= frac1*u_prof_cas(k1) - frac2*u_prof_cas(k2)
     2615         v_mod_cas(l)= frac1*v_prof_cas(k1) - frac2*v_prof_cas(k2)
     2616         ug_mod_cas(l)= frac1*ug_prof_cas(k1) - frac2*ug_prof_cas(k2)
     2617         vg_mod_cas(l)= frac1*vg_prof_cas(k1) - frac2*vg_prof_cas(k2)
     2618         w_mod_cas(l)= frac1*vitw_prof_cas(k1) - frac2*vitw_prof_cas(k2)
     2619         du_mod_cas(l)= frac1*du_prof_cas(k1) - frac2*du_prof_cas(k2)
     2620         hu_mod_cas(l)= frac1*hu_prof_cas(k1) - frac2*hu_prof_cas(k2)
     2621         vu_mod_cas(l)= frac1*vu_prof_cas(k1) - frac2*vu_prof_cas(k2)
     2622         dv_mod_cas(l)= frac1*dv_prof_cas(k1) - frac2*dv_prof_cas(k2)
     2623         hv_mod_cas(l)= frac1*hv_prof_cas(k1) - frac2*hv_prof_cas(k2)
     2624         vv_mod_cas(l)= frac1*vv_prof_cas(k1) - frac2*vv_prof_cas(k2)
     2625         dt_mod_cas(l)= frac1*dt_prof_cas(k1) - frac2*dt_prof_cas(k2)
     2626         ht_mod_cas(l)= frac1*ht_prof_cas(k1) - frac2*ht_prof_cas(k2)
     2627         vt_mod_cas(l)= frac1*vt_prof_cas(k1) - frac2*vt_prof_cas(k2)
     2628         dq_mod_cas(l)= frac1*dq_prof_cas(k1) - frac2*dq_prof_cas(k2)
     2629         hq_mod_cas(l)= frac1*hq_prof_cas(k1) - frac2*hq_prof_cas(k2)
     2630         vq_mod_cas(l)= frac1*vq_prof_cas(k1) - frac2*vq_prof_cas(k2)
     2631
     2632         endif ! play.le.plev_prof_cas(1)
     2633
     2634        else ! above max altitude of forcing file
     2635 
     2636!jyg
     2637         fact=20.*(plev_prof_cas(nlev_cas)-play(l))/plev_prof_cas(nlev_cas) !jyg
     2638         fact = max(fact,0.)                                           !jyg
     2639         fact = exp(-fact)                                             !jyg
     2640         t_mod_cas(l)= t_prof_cas(nlev_cas)                                   !jyg
     2641         q_mod_cas(l)= q_prof_cas(nlev_cas)*fact                              !jyg
     2642         u_mod_cas(l)= u_prof_cas(nlev_cas)*fact                              !jyg
     2643         v_mod_cas(l)= v_prof_cas(nlev_cas)*fact                              !jyg
     2644         ug_mod_cas(l)= ug_prof_cas(nlev_cas)*fact                              !jyg
     2645         vg_mod_cas(l)= vg_prof_cas(nlev_cas)*fact                              !jyg
     2646         w_mod_cas(l)= 0.0                                                 !jyg
     2647         du_mod_cas(l)= du_prof_cas(nlev_cas)*fact
     2648         hu_mod_cas(l)= hu_prof_cas(nlev_cas)*fact                            !jyg
     2649         vu_mod_cas(l)= vu_prof_cas(nlev_cas)*fact                            !jyg
     2650         dv_mod_cas(l)= dv_prof_cas(nlev_cas)*fact
     2651         hv_mod_cas(l)= hv_prof_cas(nlev_cas)*fact                            !jyg
     2652         vv_mod_cas(l)= vv_prof_cas(nlev_cas)*fact                            !jyg
     2653         dt_mod_cas(l)= dt_prof_cas(nlev_cas)
     2654         ht_mod_cas(l)= ht_prof_cas(nlev_cas)                                 !jyg
     2655         vt_mod_cas(l)= vt_prof_cas(nlev_cas)                                 !jyg
     2656         dq_mod_cas(l)= dq_prof_cas(nlev_cas)*fact
     2657         hq_mod_cas(l)= hq_prof_cas(nlev_cas)*fact                            !jyg
     2658         vq_mod_cas(l)= vq_prof_cas(nlev_cas)*fact                            !jyg
     2659 
     2660        endif ! play
     2661 
     2662       enddo ! l
     2663
     2664!       do l = 1,llm
     2665!       print *,'t_mod_cas(l),q_mod_cas(l),ht_mod_cas(l),hq_mod_cas(l) ',
     2666!     $        l,t_mod_cas(l),q_mod_cas(l),ht_mod_cas(l),hq_mod_cas(l)
     2667!       enddo
     2668 
     2669          return
     2670          end
     2671!*****************************************************************************
     2672!=====================================================================
    24332673       SUBROUTINE interp_dice_vertical(play,nlev_dice,nt_dice,plev_prof   &
    24342674     &         ,th_prof,qv_prof,u_prof,v_prof,o3_prof                     &

  • LMDZ5/trunk/libf/phylmd/1D_decl_cases.h

    r2126 r2191  
    240240        real u_profa(nlev_astex),v_profa(nlev_astex),w_profa(nlev_astex)
    241241        real tke_profa(nlev_astex)
    242 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    243 
     242
     243!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     244!Declarations specifiques au cas standard
     245
     246        real w_mod_cas(llm), t_mod_cas(llm),q_mod_cas(llm)
     247        real ug_mod_cas(llm),vg_mod_cas(llm)
     248        real u_mod_cas(llm),v_mod_cas(llm)
     249        real ht_mod_cas(llm),vt_mod_cas(llm),dt_mod_cas(llm),dtrad_mod_cas(llm)
     250        real hq_mod_cas(llm),vq_mod_cas(llm),dq_mod_cas(llm)
     251        real hu_mod_cas(llm),vu_mod_cas(llm),du_mod_cas(llm)
     252        real hv_mod_cas(llm),vv_mod_cas(llm),dv_mod_cas(llm)
     253!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     254

  • LMDZ5/trunk/libf/phylmd/1D_interp_cases.h

    r2126 r2191  
    597597      enddo
    598598      endif ! forcing_astex
    599 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    600 
     599
     600!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     601!---------------------------------------------------------------------
     602! Interpolation forcing standard case
     603!---------------------------------------------------------------------
     604      if (forcing_case) then
     605
     606        print*,                                                             &
     607     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/pdt_cas=',     &
     608     &    daytime,day1,(daytime-day1)*86400.,                               &
     609     &    (daytime-day1)*86400/pdt_cas
     610
     611! time interpolation:
     612        CALL interp_case_time(daytime,day1,annee_ref                        &
     613     &       ,year_ini_cas,day_ju_ini_cas,nt_cas,pdt_cas,nlev_cas       &
     614     &       ,ts_cas,plev_cas,t_cas,q_cas,u_cas,v_cas,ug_cas,vg_cas    &
     615     &       ,vitw_cas,du_cas,hu_cas,vu_cas           &
     616     &       ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas             &
     617     &       ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas                 &
     618     &       ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas         &
     619     &       ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas     &
     620     &       ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas       &
     621     &       ,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas               &
     622     &       ,sens_prof_cas)
     623
     624             ts_cur = ts_prof_cas
     625             psurf=plev_prof_cas(1)
     626
     627! vertical interpolation:
     628      CALL interp_case_vertical(play,nlev_cas,plev_prof_cas            &
     629     &         ,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas                         &
     630     &         ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas           &
     631     &         ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas           &
     632     &         ,t_mod_cas,q_mod_cas,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas                              &
     633     &         ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas               &
     634     &         ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas,mxcalc)
     635
     636
     637!calcul de l'advection verticale a partir du omega
     638!Calcul des gradients verticaux
     639!initialisation
     640      d_t_z(:)=0.
     641      d_q_z(:)=0.
     642      d_u_z(:)=0.
     643      d_v_z(:)=0.
     644      d_t_dyn_z(:)=0.
     645      d_q_dyn_z(:)=0.
     646      d_u_dyn_z(:)=0.
     647      d_v_dyn_z(:)=0.
     648      DO l=2,llm-1
     649       d_t_z(l)=(temp(l+1)-temp(l-1))/(play(l+1)-play(l-1))
     650       d_q_z(l)=(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1))
     651       d_u_z(l)=(u(l+1)-u(l-1))/(play(l+1)-play(l-1))
     652       d_v_z(l)=(v(l+1)-v(l-1))/(play(l+1)-play(l-1))
     653      ENDDO
     654      d_t_z(1)=d_t_z(2)
     655      d_q_z(1)=d_q_z(2)
     656      d_u_z(1)=d_u_z(2)
     657      d_v_z(1)=d_v_z(2)
     658      d_t_z(llm)=d_t_z(llm-1)
     659      d_q_z(llm)=d_q_z(llm-1)
     660      d_u_z(llm)=d_u_z(llm-1)
     661      d_v_z(llm)=d_v_z(llm-1)
     662
     663!Calcul de l advection verticale
     664      d_t_dyn_z(:)=w_mod_cas(:)*d_t_z(:)
     665      d_q_dyn_z(:)=w_mod_cas(:)*d_q_z(:)
     666      d_u_dyn_z(:)=w_mod_cas(:)*d_u_z(:)
     667      d_v_dyn_z(:)=w_mod_cas(:)*d_v_z(:)
     668
     669!wind nudging
     670      if (nudge_u.gt.0.) then
     671        do l=1,llm
     672           u(l)=u(l)+timestep*(u_mod_cas(l)-u(l))/(nudge_u)
     673        enddo
     674      else
     675        do l=1,llm
     676        u(l) = u_mod_cas(l)
     677        enddo
     678      endif
     679
     680      if (nudge_v.gt.0.) then
     681        do l=1,llm
     682           v(l)=v(l)+timestep*(v_mod_cas(l)-v(l))/(nudge_v)
     683        enddo
     684      else
     685        do l=1,llm
     686        v(l) = v_mod_cas(l)
     687        enddo
     688      endif
     689
     690      if (nudge_w.gt.0.) then
     691        do l=1,llm
     692           w(l)=w(l)+timestep*(w_mod_cas(l)-w(l))/(nudge_w)
     693        enddo
     694      else
     695        do l=1,llm
     696        w(l) = w_mod_cas(l)
     697        enddo
     698      endif
     699
     700!nudging of q and temp
     701      if (nudge_t.gt.0.) then
     702        do l=1,llm
     703           temp(l)=temp(l)+timestep*(t_mod_cas(l)-temp(l))/(nudge_t)
     704        enddo
     705      endif
     706      if (nudge_q.gt.0.) then
     707        do l=1,llm
     708           q(l,1)=q(l,1)+timestep*(q_mod_cas(l)-q(l,1))/(nudge_q)
     709        enddo
     710      endif
     711
     712      do l = 1, llm
     713       omega(l) = w_mod_cas(l)
     714       omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
     715       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
     716
     717!calcul advection
     718        if ((tend_u.eq.1).and.(tend_w.eq.0)) then
     719           d_u_adv(l)=du_mod_cas(l)
     720        else if ((tend_u.eq.1).and.(tend_w.eq.1)) then
     721           d_u_adv(l)=hu_mod_cas(l)-d_u_dyn_z(l)
     722        endif
     723
     724        if ((tend_v.eq.1).and.(tend_w.eq.0)) then
     725           d_v_adv(l)=dv_mod_cas(l)
     726        else if ((tend_v.eq.1).and.(tend_w.eq.1)) then
     727           d_v_adv(l)=hv_mod_cas(l)-d_v_dyn_z(l)
     728        endif
     729
     730        if ((tend_t.eq.1).and.(tend_w.eq.0)) then
     731!           d_th_adv(l)=alpha*omega(l)/rcpd+dt_mod_cas(l)
     732           d_th_adv(l)=alpha*omega(l)/rcpd-dt_mod_cas(l)
     733        else if ((tend_t.eq.1).and.(tend_w.eq.1)) then
     734!           d_th_adv(l)=alpha*omega(l)/rcpd+ht_mod_cas(l)-d_t_dyn_z(l)
     735           d_th_adv(l)=alpha*omega(l)/rcpd-ht_mod_cas(l)-d_t_dyn_z(l)
     736        endif
     737
     738        if ((tend_q.eq.1).and.(tend_w.eq.0)) then
     739!           d_q_adv(l,1)=dq_mod_cas(l)
     740           d_q_adv(l,1)=-1*dq_mod_cas(l)
     741        else if ((tend_q.eq.1).and.(tend_w.eq.1)) then
     742!           d_q_adv(l,1)=hq_mod_cas(l)-d_q_dyn_z(l)
     743           d_q_adv(l,1)=-1*hq_mod_cas(l)-d_q_dyn_z(l)
     744        endif
     745         
     746        if (tend_rayo.eq.1) then
     747           dt_cooling(l) = dtrad_mod_cas(l)
     748        else
     749           dt_cooling(l) = 0.0
     750        endif
     751      enddo
     752
     753      endif ! forcing_case
     754
     755
     756!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     757

  • LMDZ5/trunk/libf/phylmd/1D_read_forc_cases.h

    r2126 r2191  
    720720
    721721      endif ! forcing_astex
    722 
     722!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     723!---------------------------------------------------------------------
     724! Forcing from standard case :
     725!---------------------------------------------------------------------
     726
     727      if (forcing_case) then
     728
     729         write(*,*),'avant call read_1D_cas'
     730         call read_1D_cas
     731         write(*,*) 'Forcing read'
     732
     733!Time interpolation for initial conditions using TOGA interpolation routine
     734         write(*,*) 'AVT 1ere INTERPOLATION: day,day1 = ',daytime,day1   
     735      CALL interp_case_time(day,day1,annee_ref                &
     736     &         ,year_ini_cas,day_ju_ini_cas,nt_cas,pdt_cas,nlev_cas       &
     737     &         ,ts_cas,plev_cas,t_cas,q_cas,u_cas,v_cas               &
     738     &         ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas           &
     739     &         ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas             &
     740     &         ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas                 &
     741     &         ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas         &
     742     &         ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas     &
     743     &         ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas       &
     744     &         ,dq_prof_cas,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas)
     745
     746! vertical interpolation using TOGA interpolation routine:
     747!      write(*,*)'avant interp vert', t_prof
     748      CALL interp_case_vertical(play,nlev_cas,plev_prof_cas            &
     749     &         ,t_prof_cas,q_prof_cas,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas    &
     750     &         ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas           &
     751     &         ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas           &
     752     &         ,t_mod_cas,q_mod_cas,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas           &
     753     &         ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas               &
     754     &         ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas,mxcalc)
     755!       write(*,*) 'Profil initial forcing case interpole',t_mod
     756
     757! initial and boundary conditions :
     758!      tsurf = ts_prof_cas
     759      ts_cur = ts_prof_cas
     760      psurf=plev_prof_cas(1)
     761      write(*,*) 'SST initiale: ',tsurf
     762      do l = 1, llm
     763       temp(l) = t_mod_cas(l)
     764       q(l,1) = q_mod_cas(l)
     765       q(l,2) = 0.0
     766       u(l) = u_mod_cas(l)
     767       v(l) = v_mod_cas(l)
     768       omega(l) = w_mod_cas(l)
     769       omega2(l)=omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
     770
     771       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
     772!on applique le forcage total au premier pas de temps
     773!attention: signe different de toga
     774       d_th_adv(l) = alpha*omega(l)/rcpd+(ht_mod_cas(l)+vt_mod_cas(l))
     775       d_q_adv(l,1) = (hq_mod_cas(l)+vq_mod_cas(l))
     776       d_q_adv(l,2) = 0.0
     777       d_u_adv(l) = (hu_mod_cas(l)+vu_mod_cas(l))
     778       d_u_adv(l) = (hv_mod_cas(l)+vv_mod_cas(l))
     779      enddo     
     780       
     781      endif !forcing_case
     782!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

  • LMDZ5/trunk/libf/phylmd/compar1d.h

    r2181 r2191  
    33!
    44      integer :: forcing_type
     5      integer :: tend_u,tend_v,tend_w,tend_t,tend_q,tend_rayo
     6      real :: nudge_u,nudge_v,nudge_w,nudge_t,nudge_q
    57      integer :: iflag_nudge
    68      real :: nat_surf
    … …  
    3234     & qsol,qsurf,psurf,zsurf,albedo,time,time_ini,xlat,xlon,airefi,    &
    3335     & wtsurf,wqsurf,restart_runoff,xagesno,qsolinp,zpicinp,            &
    34      & forcing_type,                                                    &
     36     & forcing_type,tend_u,tend_v,tend_w,tend_t,tend_q,tend_rayo,       &
     37     & nudge_u,nudge_v,nudge_w,nudge_t,nudge_q,                         &
    3538     & iflag_nudge,                                                     &
    3639     & restart,ok_old_disvert

  • LMDZ5/trunk/libf/phylmd/lmdz1d.F90

    r2181 r2191  
    2222      USE phyaqua_mod
    2323      USE mod_1D_cases_read
     24      USE mod_1D_amma_read
    2425
    2526      implicit none
    … …  
    118119        logical :: forcing_astex   = .false.
    119120        logical :: forcing_fire    = .false.
     121        logical :: forcing_case    = .false.
    120122        integer :: type_ts_forcing ! 0 = SST constant; 1 = SST read from a file
    121123!                                                            (cf read_tsurf1d.F)
    … …  
    171173      real :: dt_cooling(llm),d_th_adv(llm),d_t_nudge(llm)
    172174      real :: d_u_nudge(llm),d_v_nudge(llm)
     175      real :: du_adv(llm),dv_adv(llm)
    173176      real :: alpha
    174177      real :: ttt
    … …  
    285288!             Different stages: soil model alone, atm. model alone
    286289!             then both models coupled
     290!forcing_type = 10 ==> forcing_case = .true.
     291!             initial profiles and large scale forcings in cas.nc
     292!             LS convergence, omega and SST imposed from CINDY-DYNAMO files
    287293!forcing_type = 40 ==> forcing_GCSSold = .true.
    288294!             initial profile from GCSS file
    … …  
    321327      elseif (forcing_type .eq.7) THEN
    322328       forcing_dice = .true.
     329      elseif (forcing_type .eq.10) THEN
     330       forcing_case = .true.
    323331      elseif (forcing_type .eq.40) THEN
    324332       forcing_GCSSold = .true.
    … …  
    428436     & (year_ini_dice,mth_ini_dice,day_ini_dice,heure_ini_dice             &
    429437     & ,day_ju_ini_dice)
     438      ELSEIF (forcing_type .eq.10) THEN
     439! Convert the initial date to Julian day
     440      print*,'time cindy',year_ini_cas,mth_ini_cas,day_ini_cas
     441      call ymds2ju                                                         &
     442     & (year_ini_cas,mth_ini_cas,day_ini_cas,heure_ini_cas &
     443     & ,day_ju_ini_cas)
     444      print*,'time cindy 2',day_ju_ini_cas
    430445      ELSEIF (forcing_type .eq.59) THEN
    431446! Convert the initial date of Sandu case to Julian day
    … …  
    943958        u(1:mxcalc)=u(1:mxcalc) + timestep*(                                &
    944959     &              du_phys(1:mxcalc)                                       &
    945      &             +du_age(1:mxcalc)                                        &
     960     &             +du_age(1:mxcalc)+du_adv(1:mxcalc)                       &
    946961     &             +d_u_nudge(1:mxcalc) )           
    947962        v(1:mxcalc)=v(1:mxcalc) + timestep*(                                 &
    948963     &              dv_phys(1:mxcalc)                                       &
    949      &             +dv_age(1:mxcalc)                                        &
     964     &             +dv_age(1:mxcalc)+dv_adv(1:mxcalc)                       &
    950965     &             +d_v_nudge(1:mxcalc) )
    951966        q(1:mxcalc,:)=q(1:mxcalc,:)+timestep*(                              &

  • LMDZ5/trunk/libf/phylmd/mod_1D_cases_read.F90

    r2119 r2191  
    22
    33!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    4 !Declarations specifiques au cas AMMA
    5         character*80 :: fich_amma
    6 ! Option du cas AMMA ou on impose la discretisation verticale (Ap,Bp)
    7         integer nlev_amma, nt_amma
    8 
    9 
    10         integer year_ini_amma, day_ini_amma, mth_ini_amma
    11         real heure_ini_amma
    12         real day_ju_ini_amma   ! Julian day of amma first day
    13         parameter (year_ini_amma=2006)
    14         parameter (mth_ini_amma=7)
    15         parameter (day_ini_amma=10)  ! 10 = 10Juil2006
    16         parameter (heure_ini_amma=0.) !0h en secondes
    17         real dt_amma
    18         parameter (dt_amma=1800.)
    19 
    20 !profils initiaux:
    21         real, allocatable::  plev_amma(:)
    22 
    23         real, allocatable::  z_amma(:)
    24         real, allocatable::  th_amma(:),q_amma(:)
    25         real, allocatable::  u_amma(:)
    26         real, allocatable::  v_amma(:)
    27 
    28         real, allocatable::  th_ammai(:),q_ammai(:)
    29         real, allocatable::  u_ammai(:)
    30         real, allocatable::  v_ammai(:)
    31         real, allocatable::  vitw_ammai(:)
    32         real, allocatable::  ht_ammai(:)
    33         real, allocatable::  hq_ammai(:)
    34         real, allocatable::  vt_ammai(:)
    35         real, allocatable::  vq_ammai(:)
    36 
    37 !forcings
    38         real, allocatable::  ht_amma(:,:)
    39         real, allocatable::  hq_amma(:,:)
    40         real, allocatable::  vitw_amma(:,:)
    41         real, allocatable::  lat_amma(:),sens_amma(:)
     4!Declarations specifiques au cas standard
     5        character*80 :: fich_cas
     6! Discrétisation
     7        integer nlev_cas, nt_cas
     8
     9
     10        integer year_ini_cas, day_ini_cas, mth_ini_cas
     11        real heure_ini_cas
     12        real day_ju_ini_cas   ! Julian day of case first day
     13        parameter (year_ini_cas=2011)
     14        parameter (mth_ini_cas=10)
     15        parameter (day_ini_cas=1)  ! 10 = 10Juil2006
     16        parameter (heure_ini_cas=0.) !0h en secondes
     17        real pdt_cas
     18        parameter (pdt_cas=3.*3600)
     19
     20!CR ATTENTION TEST AMMA
     21!        parameter (year_ini_cas=2006)
     22!        parameter (mth_ini_cas=7)
     23!        parameter (day_ini_cas=10)  ! 10 = 10Juil2006
     24!        parameter (heure_ini_cas=0.) !0h en secondes
     25!        parameter (pdt_cas=1800.)
     26
     27!profils environnementaux
     28        real, allocatable::  plev_cas(:,:)
     29
     30        real, allocatable::  z_cas(:,:)
     31        real, allocatable::  t_cas(:,:),q_cas(:,:),rh_cas(:,:)
     32        real, allocatable::  th_cas(:,:),rv_cas(:,:)
     33        real, allocatable::  u_cas(:,:)
     34        real, allocatable::  v_cas(:,:)
     35
     36!forcing
     37        real, allocatable::  ht_cas(:,:),vt_cas(:,:),dt_cas(:,:),dtrad_cas(:,:)
     38        real, allocatable::  hth_cas(:,:),vth_cas(:,:),dth_cas(:,:)
     39        real, allocatable::  hq_cas(:,:),vq_cas(:,:),dq_cas(:,:)
     40        real, allocatable::  hr_cas(:,:),vr_cas(:,:),dr_cas(:,:)
     41        real, allocatable::  hu_cas(:,:),vu_cas(:,:),du_cas(:,:)
     42        real, allocatable::  hv_cas(:,:),vv_cas(:,:),dv_cas(:,:)
     43        real, allocatable::  vitw_cas(:,:)
     44        real, allocatable::  ug_cas(:,:),vg_cas(:,:)
     45        real, allocatable::  lat_cas(:),sens_cas(:),ts_cas(:)
    4246
    4347!champs interpoles
    44         real, allocatable::  vitw_profamma(:)
    45         real, allocatable::  ht_profamma(:)
    46         real, allocatable::  hq_profamma(:)
    47         real lat_profamma,sens_profamma
    48         real, allocatable::  vt_profamma(:)
    49         real, allocatable::  vq_profamma(:)
    50         real, allocatable::  th_profamma(:)
    51         real, allocatable::  q_profamma(:)
    52         real, allocatable::  u_profamma(:)
    53         real, allocatable::  v_profamma(:)
     48        real, allocatable::  plev_prof_cas(:)
     49        real, allocatable::  t_prof_cas(:)
     50        real, allocatable::  q_prof_cas(:)
     51        real, allocatable::  u_prof_cas(:)
     52        real, allocatable::  v_prof_cas(:)       
     53
     54        real, allocatable::  vitw_prof_cas(:)
     55        real, allocatable::  ug_prof_cas(:)
     56        real, allocatable::  vg_prof_cas(:)
     57        real, allocatable::  ht_prof_cas(:)
     58        real, allocatable::  hq_prof_cas(:)
     59        real, allocatable::  vt_prof_cas(:)
     60        real, allocatable::  vq_prof_cas(:)
     61        real, allocatable::  dt_prof_cas(:)
     62        real, allocatable::  dtrad_prof_cas(:)
     63        real, allocatable::  dq_prof_cas(:)
     64        real, allocatable::  hu_prof_cas(:)
     65        real, allocatable::  hv_prof_cas(:)
     66        real, allocatable::  vu_prof_cas(:)
     67        real, allocatable::  vv_prof_cas(:)
     68        real, allocatable::  du_prof_cas(:)
     69        real, allocatable::  dv_prof_cas(:)
     70
     71
     72        real lat_prof_cas,sens_prof_cas,ts_prof_cas
     73     
    5474
    5575
    5676CONTAINS
    5777
    58 SUBROUTINE read_1D_cases
     78SUBROUTINE read_1D_cas
    5979      implicit none
    6080
    … …  
    6282
    6383      INTEGER nid,rid,ierr
    64 
    65       fich_amma='amma.nc'
    66       print*,'fich_amma ',fich_amma
    67       ierr = NF_OPEN(fich_amma,NF_NOWRITE,nid)
    68       print*,'fich_amma,NF_NOWRITE,nid ',fich_amma,NF_NOWRITE,nid
     84      INTEGER ii,jj
     85
     86      fich_cas='setup/cas.nc'
     87      print*,'fich_cas ',fich_cas
     88      ierr = NF_OPEN(fich_cas,NF_NOWRITE,nid)
     89      print*,'fich_cas,NF_NOWRITE,nid ',fich_cas,NF_NOWRITE,nid
    6990      if (ierr.NE.NF_NOERR) then
    7091         write(*,*) 'ERROR: GROS Pb opening forcings nc file '
    … …  
    7394      endif
    7495!.......................................................................
     96      ierr=NF_INQ_DIMID(nid,'lat',rid)
     97      IF (ierr.NE.NF_NOERR) THEN
     98         print*, 'Oh probleme lecture dimension lat'
     99      ENDIF
     100      ierr=NF_INQ_DIMLEN(nid,rid,ii)
     101      print*,'OK nid,rid,lat',nid,rid,ii
     102!.......................................................................
     103      ierr=NF_INQ_DIMID(nid,'lon',rid)
     104      IF (ierr.NE.NF_NOERR) THEN
     105         print*, 'Oh probleme lecture dimension lon'
     106      ENDIF
     107      ierr=NF_INQ_DIMLEN(nid,rid,jj)
     108      print*,'OK nid,rid,lat',nid,rid,jj
     109!.......................................................................
    75110      ierr=NF_INQ_DIMID(nid,'lev',rid)
    76111      IF (ierr.NE.NF_NOERR) THEN
    77112         print*, 'Oh probleme lecture dimension zz'
    78113      ENDIF
    79       ierr=NF_INQ_DIMLEN(nid,rid,nlev_amma)
    80       print*,'OK nid,rid,nlev_amma',nid,rid,nlev_amma
     114      ierr=NF_INQ_DIMLEN(nid,rid,nlev_cas)
     115      print*,'OK nid,rid,nlev_cas',nid,rid,nlev_cas
    81116!.......................................................................
    82117      ierr=NF_INQ_DIMID(nid,'time',rid)
    83118      print*,'nid,rid',nid,rid
    84       nt_amma=0
     119      nt_cas=0
    85120      IF (ierr.NE.NF_NOERR) THEN
    86121        stop 'probleme lecture dimension sens'
    87122      ENDIF
    88       ierr=NF_INQ_DIMLEN(nid,rid,nt_amma)
    89       print*,'nid,rid,nlev_amma',nid,rid,nt_amma
     123      ierr=NF_INQ_DIMLEN(nid,rid,nt_cas)
     124      print*,'nid,rid,nlev_cas',nid,rid,nt_cas
    90125
    91126!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    92 !profils initiaux:
    93         allocate(plev_amma(nlev_amma))
    94        
    95         allocate(z_amma(nlev_amma))
    96         allocate(th_amma(nlev_amma),q_amma(nlev_amma))
    97         allocate(u_amma(nlev_amma))
    98         allocate(v_amma(nlev_amma))
    99 
    100 !forcings
    101         allocate(ht_amma(nlev_amma,nt_amma))
    102         allocate(hq_amma(nlev_amma,nt_amma))
    103         allocate(vitw_amma(nlev_amma,nt_amma))
    104         allocate(lat_amma(nt_amma),sens_amma(nt_amma))
    105 
    106 !profils initiaux:
    107         allocate(th_ammai(nlev_amma),q_ammai(nlev_amma))
    108         allocate(u_ammai(nlev_amma))
    109         allocate(v_ammai(nlev_amma))
    110         allocate(vitw_ammai(nlev_amma) )
    111         allocate(ht_ammai(nlev_amma))
    112         allocate(hq_ammai(nlev_amma))
    113         allocate(vt_ammai(nlev_amma))
    114         allocate(vq_ammai(nlev_amma))
     127!profils moyens:
     128        allocate(plev_cas(nlev_cas,nt_cas))       
     129        allocate(z_cas(nlev_cas,nt_cas))
     130        allocate(t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas),rh_cas(nlev_cas,nt_cas))
     131        allocate(th_cas(nlev_cas,nt_cas),rv_cas(nlev_cas,nt_cas))
     132        allocate(u_cas(nlev_cas,nt_cas))
     133        allocate(v_cas(nlev_cas,nt_cas))
     134
     135!forcing
     136        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))
     137        allocate(hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas),dq_cas(nlev_cas,nt_cas))
     138        allocate(hth_cas(nlev_cas,nt_cas),vth_cas(nlev_cas,nt_cas),dth_cas(nlev_cas,nt_cas))
     139        allocate(hr_cas(nlev_cas,nt_cas),vr_cas(nlev_cas,nt_cas),dr_cas(nlev_cas,nt_cas))
     140        allocate(hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas),du_cas(nlev_cas,nt_cas))
     141        allocate(hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas),dv_cas(nlev_cas,nt_cas))
     142        allocate(vitw_cas(nlev_cas,nt_cas))
     143        allocate(ug_cas(nlev_cas,nt_cas))
     144        allocate(vg_cas(nlev_cas,nt_cas))
     145        allocate(lat_cas(nt_cas),sens_cas(nt_cas),ts_cas(nt_cas))
     146
    115147
    116148!champs interpoles
    117         allocate(vitw_profamma(nlev_amma))
    118         allocate(ht_profamma(nlev_amma))
    119         allocate(hq_profamma(nlev_amma))
    120         allocate(vt_profamma(nlev_amma))
    121         allocate(vq_profamma(nlev_amma))
    122         allocate(th_profamma(nlev_amma))
    123         allocate(q_profamma(nlev_amma))
    124         allocate(u_profamma(nlev_amma))
    125         allocate(v_profamma(nlev_amma))
     149        allocate(plev_prof_cas(nlev_cas))
     150        allocate(t_prof_cas(nlev_cas))
     151        allocate(q_prof_cas(nlev_cas))
     152        allocate(u_prof_cas(nlev_cas))
     153        allocate(v_prof_cas(nlev_cas))
     154
     155        allocate(vitw_prof_cas(nlev_cas))
     156        allocate(ug_prof_cas(nlev_cas))
     157        allocate(vg_prof_cas(nlev_cas))
     158        allocate(ht_prof_cas(nlev_cas))
     159        allocate(hq_prof_cas(nlev_cas))
     160        allocate(hu_prof_cas(nlev_cas))
     161        allocate(hv_prof_cas(nlev_cas))
     162        allocate(vt_prof_cas(nlev_cas))
     163        allocate(vq_prof_cas(nlev_cas))
     164        allocate(vu_prof_cas(nlev_cas))
     165        allocate(vv_prof_cas(nlev_cas))
     166        allocate(dt_prof_cas(nlev_cas))
     167        allocate(dtrad_prof_cas(nlev_cas))
     168        allocate(dq_prof_cas(nlev_cas))
     169        allocate(du_prof_cas(nlev_cas))
     170        allocate(dv_prof_cas(nlev_cas))
    126171
    127172        print*,'Allocations OK'
    128         call read_amma(nid,nlev_amma,nt_amma                                  &
    129      &     ,z_amma,plev_amma,th_amma,q_amma,u_amma,v_amma,vitw_amma         &
    130      &     ,ht_amma,hq_amma,sens_amma,lat_amma)
    131 
    132 END SUBROUTINE read_1D_cases
     173        call read_cas(nid,nlev_cas,nt_cas                                  &
     174     &     ,z_cas,plev_cas,t_cas,q_cas,rh_cas,th_cas,rv_cas,u_cas,v_cas     &
     175     &     ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas,dv_cas,hv_cas,vv_cas &
     176     &     ,dt_cas,dtrad_cas,ht_cas,vt_cas,dq_cas,hq_cas,vq_cas             &
     177     &     ,dth_cas,hth_cas,vth_cas,dr_cas,hr_cas,vr_cas,sens_cas,lat_cas,ts_cas)
     178
     179
     180END SUBROUTINE read_1D_cas
    133181
    134182
    … …  
    136184!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    137185SUBROUTINE deallocate_1D_cases
    138 !profils initiaux:
    139         deallocate(plev_amma)
     186!profils environnementaux:
     187        deallocate(plev_cas)
    140188       
    141         deallocate(z_amma)
    142         deallocate(th_amma,q_amma)
    143         deallocate(u_amma)
    144         deallocate(v_amma)
    145 
    146         deallocate(th_ammai,q_ammai)
    147         deallocate(u_ammai)
    148         deallocate(v_ammai)
    149         deallocate(vitw_ammai )
    150         deallocate(ht_ammai)
    151         deallocate(hq_ammai)
    152         deallocate(vt_ammai)
    153         deallocate(vq_ammai)
     189        deallocate(z_cas)
     190        deallocate(t_cas,q_cas,rh_cas)
     191        deallocate(th_cas,rv_cas)
     192        deallocate(u_cas)
     193        deallocate(v_cas)
    154194       
    155 !forcings
    156         deallocate(ht_amma)
    157         deallocate(hq_amma)
    158         deallocate(vitw_amma)
    159         deallocate(lat_amma,sens_amma)
     195!forcing
     196        deallocate(ht_cas,vt_cas,dt_cas,dtrad_cas)
     197        deallocate(hq_cas,vq_cas,dq_cas)
     198        deallocate(hth_cas,vth_cas,dth_cas)
     199        deallocate(hr_cas,vr_cas,dr_cas)
     200        deallocate(hu_cas,vu_cas,du_cas)
     201        deallocate(hv_cas,vv_cas,dv_cas)
     202        deallocate(vitw_cas)
     203        deallocate(ug_cas)
     204        deallocate(vg_cas)
     205        deallocate(lat_cas,sens_cas,ts_cas)
    160206
    161207!champs interpoles
    162         deallocate(vitw_profamma)
    163         deallocate(ht_profamma)
    164         deallocate(hq_profamma)
    165         deallocate(vt_profamma)
    166         deallocate(vq_profamma)
    167         deallocate(th_profamma)
    168         deallocate(q_profamma)
    169         deallocate(u_profamma)
    170         deallocate(v_profamma)
     208        deallocate(plev_prof_cas)
     209        deallocate(t_prof_cas)
     210        deallocate(q_prof_cas)
     211        deallocate(u_prof_cas)
     212        deallocate(v_prof_cas)
     213
     214        deallocate(vitw_prof_cas)
     215        deallocate(ug_prof_cas)
     216        deallocate(vg_prof_cas)
     217        deallocate(ht_prof_cas)
     218        deallocate(hq_prof_cas)
     219        deallocate(hu_prof_cas)
     220        deallocate(hv_prof_cas)
     221        deallocate(vt_prof_cas)
     222        deallocate(vq_prof_cas)
     223        deallocate(vu_prof_cas)
     224        deallocate(vv_prof_cas)
     225        deallocate(dt_prof_cas)
     226        deallocate(dtrad_prof_cas)
     227        deallocate(dq_prof_cas)
     228        deallocate(du_prof_cas)
     229        deallocate(dv_prof_cas)
     230        deallocate(t_prof_cas)
     231        deallocate(q_prof_cas)
     232        deallocate(u_prof_cas)
     233        deallocate(v_prof_cas)
     234
    171235END SUBROUTINE deallocate_1D_cases
    172236
    … …  
    174238END MODULE mod_1D_cases_read
    175239!=====================================================================
    176       subroutine read_amma(nid,nlevel,ntime                          &
    177      &     ,zz,pp,temp,qv,u,v,dw                   &
    178      &     ,dt,dq,sens,flat)
    179 
    180 !program reading forcings of the AMMA case study
     240      subroutine read_cas(nid,nlevel,ntime                          &
     241     &     ,zz,pp,temp,qv,rh,theta,rv,u,v,ug,vg,w,                   &
     242     &     du,hu,vu,dv,hv,vv,dt,dtrad,ht,vt,dq,hq,vq,                     &
     243     &     dth,hth,vth,dr,hr,vr,sens,flat,ts)
     244
     245!program reading forcing of the case study
    181246      implicit none
    182247#include "netcdf.inc"
    … …  
    184249      integer ntime,nlevel
    185250
    186       real zz(nlevel)
    187       real temp(nlevel),pp(nlevel)
    188       real qv(nlevel),u(nlevel)
    189       real v(nlevel)
    190       real dw(nlevel,ntime)
    191       real dt(nlevel,ntime)
    192       real dq(nlevel,ntime)
    193       real flat(ntime),sens(ntime)
     251      real zz(nlevel,ntime)
     252      real pp(nlevel,ntime)
     253      real temp(nlevel,ntime),qv(nlevel,ntime),rh(nlevel,ntime)
     254      real theta(nlevel,ntime),rv(nlevel,ntime)
     255      real u(nlevel,ntime)
     256      real v(nlevel,ntime)
     257      real ug(nlevel,ntime)
     258      real vg(nlevel,ntime)
     259      real w(nlevel,ntime)
     260      real du(nlevel,ntime),hu(nlevel,ntime),vu(nlevel,ntime)
     261      real dv(nlevel,ntime),hv(nlevel,ntime),vv(nlevel,ntime)
     262      real dt(nlevel,ntime),ht(nlevel,ntime),vt(nlevel,ntime)
     263      real dtrad(nlevel,ntime)
     264      real dq(nlevel,ntime),hq(nlevel,ntime),vq(nlevel,ntime)
     265      real dth(nlevel,ntime),hth(nlevel,ntime),vth(nlevel,ntime)
     266      real dr(nlevel,ntime),hr(nlevel,ntime),vr(nlevel,ntime)
     267      real flat(ntime),sens(ntime),ts(ntime)
    194268
    195269
    196270      integer nid, ierr,rid
    197271      integer nbvar3d
    198       parameter(nbvar3d=30)
     272      parameter(nbvar3d=34)
    199273      integer var3didin(nbvar3d)
    200274
    … …  
    204278           stop 'lev'
    205279         endif
    206 
    207 
    208       ierr=NF_INQ_VARID(nid,"temp",var3didin(2))
     280     
     281      ierr=NF_INQ_VARID(nid,"pp",var3didin(2))
     282         if(ierr/=NF_NOERR) then
     283           write(*,*) NF_STRERROR(ierr)
     284           stop 'plev'
     285         endif
     286
     287
     288      ierr=NF_INQ_VARID(nid,"temp",var3didin(3))
    209289         if(ierr/=NF_NOERR) then
    210290           write(*,*) NF_STRERROR(ierr)
    … …  
    212292         endif
    213293
    214       ierr=NF_INQ_VARID(nid,"qv",var3didin(3))
     294      ierr=NF_INQ_VARID(nid,"qv",var3didin(4))
    215295         if(ierr/=NF_NOERR) then
    216296           write(*,*) NF_STRERROR(ierr)
    … …  
    218298         endif
    219299
    220       ierr=NF_INQ_VARID(nid,"u",var3didin(4))
     300      ierr=NF_INQ_VARID(nid,"rh",var3didin(5))
     301         if(ierr/=NF_NOERR) then
     302           write(*,*) NF_STRERROR(ierr)
     303           stop 'rh'
     304         endif
     305
     306      ierr=NF_INQ_VARID(nid,"theta",var3didin(6))
     307         if(ierr/=NF_NOERR) then
     308           write(*,*) NF_STRERROR(ierr)
     309           stop 'theta'
     310         endif
     311
     312      ierr=NF_INQ_VARID(nid,"rv",var3didin(7))
     313         if(ierr/=NF_NOERR) then
     314           write(*,*) NF_STRERROR(ierr)
     315           stop 'rv'
     316         endif
     317
     318
     319      ierr=NF_INQ_VARID(nid,"u",var3didin(8))
    221320         if(ierr/=NF_NOERR) then
    222321           write(*,*) NF_STRERROR(ierr)
    … …  
    224323         endif
    225324
    226       ierr=NF_INQ_VARID(nid,"v",var3didin(5))
     325      ierr=NF_INQ_VARID(nid,"v",var3didin(9))
    227326         if(ierr/=NF_NOERR) then
    228327           write(*,*) NF_STRERROR(ierr)
    … …  
    230329         endif
    231330
    232       ierr=NF_INQ_VARID(nid,"dw",var3didin(6))
    233          if(ierr/=NF_NOERR) then
    234            write(*,*) NF_STRERROR(ierr)
    235            stop 'dw'
    236          endif
    237 
    238       ierr=NF_INQ_VARID(nid,"dt",var3didin(7))
    239          if(ierr/=NF_NOERR) then
    240            write(*,*) NF_STRERROR(ierr)
    241            stop 'dt'
    242          endif
    243 
    244       ierr=NF_INQ_VARID(nid,"dq",var3didin(8))
    245          if(ierr/=NF_NOERR) then
    246            write(*,*) NF_STRERROR(ierr)
    247            stop 'dq'
     331       ierr=NF_INQ_VARID(nid,"ug",var3didin(10))
     332         if(ierr/=NF_NOERR) then
     333           write(*,*) NF_STRERROR(ierr)
     334           stop 'ug'
     335         endif
     336
     337      ierr=NF_INQ_VARID(nid,"vg",var3didin(11))
     338         if(ierr/=NF_NOERR) then
     339           write(*,*) NF_STRERROR(ierr)
     340           stop 'vg'
     341         endif
     342
     343      ierr=NF_INQ_VARID(nid,"w",var3didin(12))
     344         if(ierr/=NF_NOERR) then
     345           write(*,*) NF_STRERROR(ierr)
     346           stop 'w'
     347         endif
     348
     349      ierr=NF_INQ_VARID(nid,"advu",var3didin(13))
     350         if(ierr/=NF_NOERR) then
     351           write(*,*) NF_STRERROR(ierr)
     352           stop 'advu'
     353         endif
     354
     355      ierr=NF_INQ_VARID(nid,"hu",var3didin(14))
     356         if(ierr/=NF_NOERR) then
     357           write(*,*) NF_STRERROR(ierr)
     358           stop 'hu'
     359         endif
     360
     361       ierr=NF_INQ_VARID(nid,"vu",var3didin(15))
     362         if(ierr/=NF_NOERR) then
     363           write(*,*) NF_STRERROR(ierr)
     364           stop 'vu'
     365         endif
     366
     367       ierr=NF_INQ_VARID(nid,"advv",var3didin(16))
     368         if(ierr/=NF_NOERR) then
     369           write(*,*) NF_STRERROR(ierr)
     370           stop 'advv'
     371         endif
     372
     373      ierr=NF_INQ_VARID(nid,"hv",var3didin(17))
     374         if(ierr/=NF_NOERR) then
     375           write(*,*) NF_STRERROR(ierr)
     376           stop 'hv'
     377         endif
     378
     379       ierr=NF_INQ_VARID(nid,"vv",var3didin(18))
     380         if(ierr/=NF_NOERR) then
     381           write(*,*) NF_STRERROR(ierr)
     382           stop 'vv'
     383         endif
     384
     385      ierr=NF_INQ_VARID(nid,"advT",var3didin(19))
     386         if(ierr/=NF_NOERR) then
     387           write(*,*) NF_STRERROR(ierr)
     388           stop 'advT'
     389         endif
     390
     391      ierr=NF_INQ_VARID(nid,"hT",var3didin(20))
     392         if(ierr/=NF_NOERR) then
     393           write(*,*) NF_STRERROR(ierr)
     394           stop 'hT'
     395         endif
     396
     397      ierr=NF_INQ_VARID(nid,"vT",var3didin(21))
     398         if(ierr/=NF_NOERR) then
     399           write(*,*) NF_STRERROR(ierr)
     400           stop 'vT'
     401         endif
     402
     403      ierr=NF_INQ_VARID(nid,"advq",var3didin(22))
     404         if(ierr/=NF_NOERR) then
     405           write(*,*) NF_STRERROR(ierr)
     406           stop 'advq'
    248407         endif
    249408     
    250       ierr=NF_INQ_VARID(nid,"sens",var3didin(9))
     409      ierr=NF_INQ_VARID(nid,"hq",var3didin(23))
     410         if(ierr/=NF_NOERR) then
     411           write(*,*) NF_STRERROR(ierr)
     412           stop 'hq'
     413         endif
     414
     415      ierr=NF_INQ_VARID(nid,"vq",var3didin(24))
     416         if(ierr/=NF_NOERR) then
     417           write(*,*) NF_STRERROR(ierr)
     418           stop 'vq'
     419         endif
     420
     421      ierr=NF_INQ_VARID(nid,"advth",var3didin(25))
     422         if(ierr/=NF_NOERR) then
     423           write(*,*) NF_STRERROR(ierr)
     424           stop 'advth'
     425         endif
     426
     427      ierr=NF_INQ_VARID(nid,"hth",var3didin(26))
     428         if(ierr/=NF_NOERR) then
     429           write(*,*) NF_STRERROR(ierr)
     430           stop 'hth'
     431         endif
     432
     433      ierr=NF_INQ_VARID(nid,"vth",var3didin(27))
     434         if(ierr/=NF_NOERR) then
     435           write(*,*) NF_STRERROR(ierr)
     436           stop 'vth'
     437         endif
     438
     439      ierr=NF_INQ_VARID(nid,"advr",var3didin(28))
     440         if(ierr/=NF_NOERR) then
     441           write(*,*) NF_STRERROR(ierr)
     442           stop 'advr'
     443         endif
     444     
     445      ierr=NF_INQ_VARID(nid,"hr",var3didin(29))
     446         if(ierr/=NF_NOERR) then
     447           write(*,*) NF_STRERROR(ierr)
     448           stop 'hr'
     449         endif
     450
     451      ierr=NF_INQ_VARID(nid,"vr",var3didin(30))
     452         if(ierr/=NF_NOERR) then
     453           write(*,*) NF_STRERROR(ierr)
     454           stop 'vr'
     455         endif
     456
     457      ierr=NF_INQ_VARID(nid,"radT",var3didin(31))
     458         if(ierr/=NF_NOERR) then
     459           write(*,*) NF_STRERROR(ierr)
     460           stop 'radT'
     461         endif
     462
     463      ierr=NF_INQ_VARID(nid,"sens",var3didin(32))
    251464         if(ierr/=NF_NOERR) then
    252465           write(*,*) NF_STRERROR(ierr)
    … …  
    254467         endif
    255468
    256       ierr=NF_INQ_VARID(nid,"flat",var3didin(10))
     469      ierr=NF_INQ_VARID(nid,"flat",var3didin(33))
    257470         if(ierr/=NF_NOERR) then
    258471           write(*,*) NF_STRERROR(ierr)
    … …  
    260473         endif
    261474
    262       ierr=NF_INQ_VARID(nid,"pp",var3didin(11))
    263          if(ierr/=NF_NOERR) then
    264            write(*,*) NF_STRERROR(ierr)
    265       endif
    266 
    267 !dimensions lecture
    268 !      call catchaxis(nid,ntime,nlevel,time,z,ierr)
     475      ierr=NF_INQ_VARID(nid,"ts",var3didin(34))
     476         if(ierr/=NF_NOERR) then
     477           write(*,*) NF_STRERROR(ierr)
     478           stop 'ts'
     479         endif
    269480 
    270481#ifdef NC_DOUBLE
    … …  
    280491
    281492#ifdef NC_DOUBLE
    282          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(2),temp)
    283 #else
    284          ierr = NF_GET_VAR_REAL(nid,var3didin(2),temp)
    285 #endif
    286          if(ierr/=NF_NOERR) then
    287             write(*,*) NF_STRERROR(ierr)
    288             stop "getvarup"
    289          endif
    290 !          write(*,*)'lecture th ok',temp
    291 
    292 #ifdef NC_DOUBLE
    293          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(3),qv)
    294 #else
    295          ierr = NF_GET_VAR_REAL(nid,var3didin(3),qv)
     493         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(2),pp)
     494#else
     495         ierr = NF_GET_VAR_REAL(nid,var3didin(2),pp)
     496#endif
     497         if(ierr/=NF_NOERR) then
     498            write(*,*) NF_STRERROR(ierr)
     499            stop "getvarup"
     500         endif
     501!          write(*,*)'lecture pp ok',pp
     502
     503
     504#ifdef NC_DOUBLE
     505         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(3),temp)
     506#else
     507         ierr = NF_GET_VAR_REAL(nid,var3didin(3),temp)
     508#endif
     509         if(ierr/=NF_NOERR) then
     510            write(*,*) NF_STRERROR(ierr)
     511            stop "getvarup"
     512         endif
     513!          write(*,*)'lecture T ok',temp
     514
     515#ifdef NC_DOUBLE
     516         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(4),qv)
     517#else
     518         ierr = NF_GET_VAR_REAL(nid,var3didin(4),qv)
    296519#endif
    297520         if(ierr/=NF_NOERR) then
    … …  
    302525 
    303526#ifdef NC_DOUBLE
    304          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(4),u)
    305 #else
    306          ierr = NF_GET_VAR_REAL(nid,var3didin(4),u)
     527         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(5),rh)
     528#else
     529         ierr = NF_GET_VAR_REAL(nid,var3didin(5),rh)
     530#endif
     531         if(ierr/=NF_NOERR) then
     532            write(*,*) NF_STRERROR(ierr)
     533            stop "getvarup"
     534         endif
     535!          write(*,*)'lecture rh ok',rh
     536
     537#ifdef NC_DOUBLE
     538         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(6),theta)
     539#else
     540         ierr = NF_GET_VAR_REAL(nid,var3didin(6),theta)
     541#endif
     542         if(ierr/=NF_NOERR) then
     543            write(*,*) NF_STRERROR(ierr)
     544            stop "getvarup"
     545         endif
     546!          write(*,*)'lecture theta ok',theta
     547
     548#ifdef NC_DOUBLE
     549         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(7),rv)
     550#else
     551         ierr = NF_GET_VAR_REAL(nid,var3didin(7),rv)
     552#endif
     553         if(ierr/=NF_NOERR) then
     554            write(*,*) NF_STRERROR(ierr)
     555            stop "getvarup"
     556         endif
     557!          write(*,*)'lecture rv ok',rv
     558
     559#ifdef NC_DOUBLE
     560         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(8),u)
     561#else
     562         ierr = NF_GET_VAR_REAL(nid,var3didin(8),u)
    307563#endif
    308564         if(ierr/=NF_NOERR) then
    … …  
    313569
    314570#ifdef NC_DOUBLE
    315          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(5),v)
    316 #else
    317          ierr = NF_GET_VAR_REAL(nid,var3didin(5),v)
     571         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(9),v)
     572#else
     573         ierr = NF_GET_VAR_REAL(nid,var3didin(9),v)
    318574#endif
    319575         if(ierr/=NF_NOERR) then
    … …  
    324580
    325581#ifdef NC_DOUBLE
    326          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(6),dw)
    327 #else
    328          ierr = NF_GET_VAR_REAL(nid,var3didin(6),dw)
    329 #endif
    330          if(ierr/=NF_NOERR) then
    331             write(*,*) NF_STRERROR(ierr)
    332             stop "getvarup"
    333          endif
    334 !          write(*,*)'lecture w ok',dw
    335 
    336 #ifdef NC_DOUBLE
    337          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(7),dt)
    338 #else
    339          ierr = NF_GET_VAR_REAL(nid,var3didin(7),dt)
     582         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(10),ug)
     583#else
     584         ierr = NF_GET_VAR_REAL(nid,var3didin(10),ug)
     585#endif
     586         if(ierr/=NF_NOERR) then
     587            write(*,*) NF_STRERROR(ierr)
     588            stop "getvarup"
     589         endif
     590!          write(*,*)'lecture ug ok',ug
     591
     592#ifdef NC_DOUBLE
     593         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(11),vg)
     594#else
     595         ierr = NF_GET_VAR_REAL(nid,var3didin(11),vg)
     596#endif
     597         if(ierr/=NF_NOERR) then
     598            write(*,*) NF_STRERROR(ierr)
     599            stop "getvarup"
     600         endif
     601!          write(*,*)'lecture vg ok',vg
     602
     603#ifdef NC_DOUBLE
     604         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(12),w)
     605#else
     606         ierr = NF_GET_VAR_REAL(nid,var3didin(12),w)
     607#endif
     608         if(ierr/=NF_NOERR) then
     609            write(*,*) NF_STRERROR(ierr)
     610            stop "getvarup"
     611         endif
     612!          write(*,*)'lecture w ok',w
     613
     614#ifdef NC_DOUBLE
     615         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(13),du)
     616#else
     617         ierr = NF_GET_VAR_REAL(nid,var3didin(13),du)
     618#endif
     619         if(ierr/=NF_NOERR) then
     620            write(*,*) NF_STRERROR(ierr)
     621            stop "getvarup"
     622         endif
     623!          write(*,*)'lecture du ok',du
     624
     625#ifdef NC_DOUBLE
     626         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(14),hu)
     627#else
     628         ierr = NF_GET_VAR_REAL(nid,var3didin(14),hu)
     629#endif
     630         if(ierr/=NF_NOERR) then
     631            write(*,*) NF_STRERROR(ierr)
     632            stop "getvarup"
     633         endif
     634!          write(*,*)'lecture hu ok',hu
     635
     636#ifdef NC_DOUBLE
     637         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(15),vu)
     638#else
     639         ierr = NF_GET_VAR_REAL(nid,var3didin(15),vu)
     640#endif
     641         if(ierr/=NF_NOERR) then
     642            write(*,*) NF_STRERROR(ierr)
     643            stop "getvarup"
     644         endif
     645!          write(*,*)'lecture vu ok',vu
     646
     647#ifdef NC_DOUBLE
     648         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(16),dv)
     649#else
     650         ierr = NF_GET_VAR_REAL(nid,var3didin(16),dv)
     651#endif
     652         if(ierr/=NF_NOERR) then
     653            write(*,*) NF_STRERROR(ierr)
     654            stop "getvarup"
     655         endif
     656!          write(*,*)'lecture dv ok',dv
     657
     658#ifdef NC_DOUBLE
     659         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(17),hv)
     660#else
     661         ierr = NF_GET_VAR_REAL(nid,var3didin(17),hv)
     662#endif
     663         if(ierr/=NF_NOERR) then
     664            write(*,*) NF_STRERROR(ierr)
     665            stop "getvarup"
     666         endif
     667!          write(*,*)'lecture hv ok',hv
     668
     669#ifdef NC_DOUBLE
     670         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(18),vv)
     671#else
     672         ierr = NF_GET_VAR_REAL(nid,var3didin(18),vv)
     673#endif
     674         if(ierr/=NF_NOERR) then
     675            write(*,*) NF_STRERROR(ierr)
     676            stop "getvarup"
     677         endif
     678!          write(*,*)'lecture vv ok',vv
     679
     680#ifdef NC_DOUBLE
     681         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(19),dt)
     682#else
     683         ierr = NF_GET_VAR_REAL(nid,var3didin(19),dt)
    340684#endif
    341685         if(ierr/=NF_NOERR) then
    … …  
    346690
    347691#ifdef NC_DOUBLE
    348          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(8),dq)
    349 #else
    350          ierr = NF_GET_VAR_REAL(nid,var3didin(8),dq)
     692         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(20),ht)
     693#else
     694         ierr = NF_GET_VAR_REAL(nid,var3didin(20),ht)
     695#endif
     696         if(ierr/=NF_NOERR) then
     697            write(*,*) NF_STRERROR(ierr)
     698            stop "getvarup"
     699         endif
     700!          write(*,*)'lecture ht ok',ht
     701
     702#ifdef NC_DOUBLE
     703         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(21),vt)
     704#else
     705         ierr = NF_GET_VAR_REAL(nid,var3didin(21),vt)
     706#endif
     707         if(ierr/=NF_NOERR) then
     708            write(*,*) NF_STRERROR(ierr)
     709            stop "getvarup"
     710         endif
     711!          write(*,*)'lecture vt ok',vt
     712
     713#ifdef NC_DOUBLE
     714         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(22),dq)
     715#else
     716         ierr = NF_GET_VAR_REAL(nid,var3didin(22),dq)
    351717#endif
    352718         if(ierr/=NF_NOERR) then
    … …  
    357723
    358724#ifdef NC_DOUBLE
    359          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(9),sens)
    360 #else
    361          ierr = NF_GET_VAR_REAL(nid,var3didin(9),sens)
     725         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(23),hq)
     726#else
     727         ierr = NF_GET_VAR_REAL(nid,var3didin(23),hq)
     728#endif
     729         if(ierr/=NF_NOERR) then
     730            write(*,*) NF_STRERROR(ierr)
     731            stop "getvarup"
     732         endif
     733!          write(*,*)'lecture hq ok',hq
     734
     735#ifdef NC_DOUBLE
     736         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(24),vq)
     737#else
     738         ierr = NF_GET_VAR_REAL(nid,var3didin(24),vq)
     739#endif
     740         if(ierr/=NF_NOERR) then
     741            write(*,*) NF_STRERROR(ierr)
     742            stop "getvarup"
     743         endif
     744!          write(*,*)'lecture vq ok',vq
     745
     746#ifdef NC_DOUBLE
     747         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(25),dth)
     748#else
     749         ierr = NF_GET_VAR_REAL(nid,var3didin(25),dth)
     750#endif
     751         if(ierr/=NF_NOERR) then
     752            write(*,*) NF_STRERROR(ierr)
     753            stop "getvarup"
     754         endif
     755!          write(*,*)'lecture dth ok',dth
     756
     757#ifdef NC_DOUBLE
     758         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(26),hth)
     759#else
     760         ierr = NF_GET_VAR_REAL(nid,var3didin(26),hth)
     761#endif
     762         if(ierr/=NF_NOERR) then
     763            write(*,*) NF_STRERROR(ierr)
     764            stop "getvarup"
     765         endif
     766!          write(*,*)'lecture hth ok',hth
     767
     768#ifdef NC_DOUBLE
     769         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(27),vth)
     770#else
     771         ierr = NF_GET_VAR_REAL(nid,var3didin(27),vth)
     772#endif
     773         if(ierr/=NF_NOERR) then
     774            write(*,*) NF_STRERROR(ierr)
     775            stop "getvarup"
     776         endif
     777!          write(*,*)'lecture vth ok',vth
     778
     779#ifdef NC_DOUBLE
     780         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(28),dr)
     781#else
     782         ierr = NF_GET_VAR_REAL(nid,var3didin(28),dr)
     783#endif
     784         if(ierr/=NF_NOERR) then
     785            write(*,*) NF_STRERROR(ierr)
     786            stop "getvarup"
     787         endif
     788!          write(*,*)'lecture dr ok',dr
     789
     790#ifdef NC_DOUBLE
     791         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(29),hr)
     792#else
     793         ierr = NF_GET_VAR_REAL(nid,var3didin(29),hr)
     794#endif
     795         if(ierr/=NF_NOERR) then
     796            write(*,*) NF_STRERROR(ierr)
     797            stop "getvarup"
     798         endif
     799!          write(*,*)'lecture hr ok',hr
     800
     801#ifdef NC_DOUBLE
     802         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(30),vr)
     803#else
     804         ierr = NF_GET_VAR_REAL(nid,var3didin(30),vr)
     805#endif
     806         if(ierr/=NF_NOERR) then
     807            write(*,*) NF_STRERROR(ierr)
     808            stop "getvarup"
     809         endif
     810!          write(*,*)'lecture vr ok',vr
     811
     812#ifdef NC_DOUBLE
     813         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(31),dtrad)
     814#else
     815         ierr = NF_GET_VAR_REAL(nid,var3didin(31),dtrad)
     816#endif
     817         if(ierr/=NF_NOERR) then
     818            write(*,*) NF_STRERROR(ierr)
     819            stop "getvarup"
     820         endif
     821!          write(*,*)'lecture dtrad ok',dtrad
     822
     823#ifdef NC_DOUBLE
     824         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(32),sens)
     825#else
     826         ierr = NF_GET_VAR_REAL(nid,var3didin(32),sens)
    362827#endif
    363828         if(ierr/=NF_NOERR) then
    … …  
    368833
    369834#ifdef NC_DOUBLE
    370          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(10),flat)
    371 #else
    372          ierr = NF_GET_VAR_REAL(nid,var3didin(10),flat)
     835         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(33),flat)
     836#else
     837         ierr = NF_GET_VAR_REAL(nid,var3didin(33),flat)
    373838#endif
    374839         if(ierr/=NF_NOERR) then
    … …  
    379844
    380845#ifdef NC_DOUBLE
    381          ierr = NF_GET_VAR_DOUBLE(nid,var3didin(11),pp)
    382 #else
    383          ierr = NF_GET_VAR_REAL(nid,var3didin(11),pp)
    384 #endif
    385          if(ierr/=NF_NOERR) then
    386             write(*,*) NF_STRERROR(ierr)
    387             stop "getvarup"
    388          endif
    389 !          write(*,*)'lecture pp ok',pp
     846         ierr = NF_GET_VAR_DOUBLE(nid,var3didin(34),ts)
     847#else
     848         ierr = NF_GET_VAR_REAL(nid,var3didin(34),ts)
     849#endif
     850         if(ierr/=NF_NOERR) then
     851            write(*,*) NF_STRERROR(ierr)
     852            stop "getvarup"
     853         endif
     854!          write(*,*)'lecture ts ok',ts
    390855
    391856         return
    392          end subroutine read_amma
     857         end subroutine read_cas
    393858!======================================================================
    394         SUBROUTINE interp_amma_time(day,day1,annee_ref                     &
    395      &         ,year_ini_amma,day_ini_amma,nt_amma,dt_amma,nlev_amma       &
    396      &         ,vitw_amma,ht_amma,hq_amma,lat_amma,sens_amma               &
    397      &         ,vitw_prof,ht_prof,hq_prof,lat_prof,sens_prof)
     859        SUBROUTINE interp_case_time(day,day1,annee_ref                &
     860     &         ,year_ini_cas,day_ini_cas,nt_cas,pdt_cas,nlev_cas      &
     861     &         ,ts_cas,plev_cas,t_cas,q_cas,u_cas,v_cas               &
     862     &         ,ug_cas,vg_cas,vitw_cas,du_cas,hu_cas,vu_cas           &
     863     &         ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas   &
     864     &         ,dq_cas,hq_cas,vq_cas,lat_cas,sens_cas          &
     865     &         ,ts_prof_cas,plev_prof_cas,t_prof_cas,q_prof_cas       &
     866     &         ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas         &
     867     &         ,vitw_prof_cas,du_prof_cas,hu_prof_cas,vu_prof_cas     &
     868     &         ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas       &
     869     &         ,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas    &
     870     &         ,hq_prof_cas,vq_prof_cas,lat_prof_cas,sens_prof_cas)
     871         
     872
    398873        implicit none
    399874
    … …  
    403878! day: current julian day (e.g. 717538.2)
    404879! day1: first day of the simulation
    405 ! nt_amma: total nb of data in the forcing (e.g. 48 for AMMA)
    406 ! dt_amma: total time interval (in sec) between 2 forcing data (e.g. 30min for AMMA)
     880! nt_cas: total nb of data in the forcing
     881! pdt_cas: total time interval (in sec) between 2 forcing data
    407882!---------------------------------------------------------------------------------------
    408883
    … …  
    411886! inputs:
    412887        integer annee_ref
    413         integer nt_amma,nlev_amma
    414         integer year_ini_amma
    415         real day, day1,day_ini_amma,dt_amma
    416         real vitw_amma(nlev_amma,nt_amma)
    417         real ht_amma(nlev_amma,nt_amma)
    418         real hq_amma(nlev_amma,nt_amma)
    419         real lat_amma(nt_amma)
    420         real sens_amma(nt_amma)
     888        integer nt_cas,nlev_cas
     889        integer year_ini_cas
     890        real day_ini_cas
     891        real day, day1,pdt_cas
     892        real ts_cas(nt_cas)
     893        real plev_cas(nlev_cas,nt_cas)
     894        real t_cas(nlev_cas,nt_cas),q_cas(nlev_cas,nt_cas)
     895        real u_cas(nlev_cas,nt_cas),v_cas(nlev_cas,nt_cas)
     896        real ug_cas(nlev_cas,nt_cas),vg_cas(nlev_cas,nt_cas)
     897        real vitw_cas(nlev_cas,nt_cas)
     898        real du_cas(nlev_cas,nt_cas),hu_cas(nlev_cas,nt_cas),vu_cas(nlev_cas,nt_cas)
     899        real dv_cas(nlev_cas,nt_cas),hv_cas(nlev_cas,nt_cas),vv_cas(nlev_cas,nt_cas)
     900        real dt_cas(nlev_cas,nt_cas),ht_cas(nlev_cas,nt_cas),vt_cas(nlev_cas,nt_cas)
     901        real dtrad_cas(nlev_cas,nt_cas)
     902        real dq_cas(nlev_cas,nt_cas),hq_cas(nlev_cas,nt_cas),vq_cas(nlev_cas,nt_cas)
     903        real lat_cas(nt_cas)
     904        real sens_cas(nt_cas)
     905
    421906! outputs:
    422         real vitw_prof(nlev_amma)
    423         real ht_prof(nlev_amma)
    424         real hq_prof(nlev_amma)
    425         real lat_prof,sens_prof
     907        real plev_prof_cas(nlev_cas)
     908        real t_prof_cas(nlev_cas),q_prof_cas(nlev_cas)
     909        real u_prof_cas(nlev_cas),v_prof_cas(nlev_cas)
     910        real ug_prof_cas(nlev_cas),vg_prof_cas(nlev_cas)
     911        real vitw_prof_cas(nlev_cas)
     912        real du_prof_cas(nlev_cas),hu_prof_cas(nlev_cas),vu_prof_cas(nlev_cas)
     913        real dv_prof_cas(nlev_cas),hv_prof_cas(nlev_cas),vv_prof_cas(nlev_cas)
     914        real dt_prof_cas(nlev_cas),ht_prof_cas(nlev_cas),vt_prof_cas(nlev_cas)
     915        real dtrad_prof_cas(nlev_cas)
     916        real dq_prof_cas(nlev_cas),hq_prof_cas(nlev_cas),vq_prof_cas(nlev_cas)
     917        real lat_prof_cas,sens_prof_cas,ts_prof_cas
    426918! local:
    427         integer it_amma1, it_amma2,k
    428         real timeit,time_amma1,time_amma2,frac
    429 
    430 
    431         if (forcing_type.eq.6) then
     919        integer it_cas1, it_cas2,k
     920        real timeit,time_cas1,time_cas2,frac
     921
     922
     923        print*,'Check time',day1,day_ini_cas,day_ini_cas+90
     924
     925        if ((forcing_type.eq.10).and.(1.eq.1)) then
     926! Check that initial day of the simulation consistent with the case:
     927       if (annee_ref.ne.2011) then
     928        print*,'Pour CINDY, annee_ref doit etre 2011'
     929        print*,'Changer annee_ref dans run.def'
     930        stop
     931       endif
     932       if (annee_ref.eq.2011 .and. day1.lt.day_ini_cas) then
     933        print*,'CINDY a débuté le 1 octobre 2011',day1,day_ini_cas
     934        print*,'Changer dayref dans run.def'
     935        stop
     936       endif
     937       if (annee_ref.eq.2011 .and. day1.gt.day_ini_cas+90) then
     938        print*,'CINDY a fini le 31 decembre'
     939        print*,'Changer dayref ou nday dans run.def',day1,day_ini_cas+90
     940        stop
     941       endif
     942       endif
     943!CR:ATTENTION TEST AMMA
     944      if ((forcing_type.eq.10).and.(1.eq.0)) then
    432945! Check that initial day of the simulation consistent with AMMA case:
    433946       if (annee_ref.ne.2006) then
    … …  
    436949        stop
    437950       endif
    438        if (annee_ref.eq.2006 .and. day1.lt.day_ini_amma) then
    439         print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_amma
     951       if (annee_ref.eq.2006 .and. day1.lt.day_ini_cas) then
     952        print*,'AMMA a débuté le 10 juillet 2006',day1,day_ini_cas
    440953        print*,'Changer dayref dans run.def'
    441954        stop
    442955       endif
    443        if (annee_ref.eq.2006 .and. day1.gt.day_ini_amma+1) then
     956       if (annee_ref.eq.2006 .and. day1.gt.day_ini_cas+1) then
    444957        print*,'AMMA a fini le 11 juillet'
    445958        print*,'Changer dayref ou nday dans run.def'
    … …  
    448961       endif
    449962
    450 ! Determine timestep relative to the 1st day of AMMA:
     963! Determine timestep relative to the 1st day:
    451964!       timeit=(day-day1)*86400.
    452965!       if (annee_ref.eq.1992) then
    453 !        timeit=(day-day_ini_toga)*86400.
     966!        timeit=(day-day_ini_cas)*86400.
    454967!       else
    455968!        timeit=(day+61.-1.)*86400. ! 61 days between Nov01 and Dec31 1992
    456969!       endif
    457       timeit=(day-day_ini_amma)*86400
     970      timeit=(day-day_ini_cas)*86400
    458971
    459972! Determine the closest observation times:
    460 !       it_amma1=INT(timeit/dt_amma)+1
    461 !       it_amma2=it_amma1 + 1
    462 !       time_amma1=(it_amma1-1)*dt_amma
    463 !       time_amma2=(it_amma2-1)*dt_amma
    464 
    465        it_amma1=INT(timeit/dt_amma)+1
    466        IF (it_amma1 .EQ. nt_amma) THEN
    467        it_amma2=it_amma1
     973!       it_cas1=INT(timeit/pdt_cas)+1
     974!       it_cas2=it_cas1 + 1
     975!       time_cas1=(it_cas1-1)*pdt_cas
     976!       time_cas2=(it_cas2-1)*pdt_cas
     977
     978       it_cas1=INT(timeit/pdt_cas)+1
     979       IF (it_cas1 .EQ. nt_cas) THEN
     980       it_cas2=it_cas1
    468981       ELSE
    469        it_amma2=it_amma1 + 1
     982       it_cas2=it_cas1 + 1
    470983       ENDIF
    471        time_amma1=(it_amma1-1)*dt_amma
    472        time_amma2=(it_amma2-1)*dt_amma
    473 
    474        if (it_amma1 .gt. nt_amma) then
    475         write(*,*) 'PB-stop: day, it_amma1, it_amma2, timeit: '            &
    476      &        ,day,day_ini_amma,it_amma1,it_amma2,timeit/86400.
     984       time_cas1=(it_cas1-1)*pdt_cas
     985       time_cas2=(it_cas2-1)*pdt_cas
     986
     987       if (it_cas1 .gt. nt_cas) then
     988        write(*,*) 'PB-stop: day, it_cas1, it_cas2, timeit: '            &
     989     &        ,day,day_ini_cas,it_cas1,it_cas2,timeit/86400.
    477990        stop
    478991       endif
    479992
    480993! time interpolation:
    481        frac=(time_amma2-timeit)/(time_amma2-time_amma1)
     994       frac=(time_cas2-timeit)/(time_cas2-time_cas1)
    482995       frac=max(frac,0.0)
    483996
    484        lat_prof = lat_amma(it_amma2)                                       &
    485      &          -frac*(lat_amma(it_amma2)-lat_amma(it_amma1))
    486        sens_prof = sens_amma(it_amma2)                                     &
    487      &          -frac*(sens_amma(it_amma2)-sens_amma(it_amma1))
    488 
    489        do k=1,nlev_amma
    490         vitw_prof(k) = vitw_amma(k,it_amma2)                               &
    491      &          -frac*(vitw_amma(k,it_amma2)-vitw_amma(k,it_amma1))
    492         ht_prof(k) = ht_amma(k,it_amma2)                                   &
    493      &          -frac*(ht_amma(k,it_amma2)-ht_amma(k,it_amma1))
    494         hq_prof(k) = hq_amma(k,it_amma2)                                   &
    495      &          -frac*(hq_amma(k,it_amma2)-hq_amma(k,it_amma1))
     997       lat_prof_cas = lat_cas(it_cas2)                                       &
     998     &          -frac*(lat_cas(it_cas2)-lat_cas(it_cas1))
     999       sens_prof_cas = sens_cas(it_cas2)                                     &
     1000     &          -frac*(sens_cas(it_cas2)-sens_cas(it_cas1))
     1001       ts_prof_cas = ts_cas(it_cas2)                                     &
     1002     &          -frac*(ts_cas(it_cas2)-ts_cas(it_cas1))
     1003
     1004       do k=1,nlev_cas
     1005        plev_prof_cas(k) = plev_cas(k,it_cas2)                               &
     1006     &          -frac*(plev_cas(k,it_cas2)-plev_cas(k,it_cas1))
     1007        t_prof_cas(k) = t_cas(k,it_cas2)                               &
     1008     &          -frac*(t_cas(k,it_cas2)-t_cas(k,it_cas1))
     1009        q_prof_cas(k) = q_cas(k,it_cas2)                               &
     1010     &          -frac*(q_cas(k,it_cas2)-q_cas(k,it_cas1))
     1011        u_prof_cas(k) = u_cas(k,it_cas2)                               &
     1012     &          -frac*(u_cas(k,it_cas2)-u_cas(k,it_cas1))
     1013        v_prof_cas(k) = v_cas(k,it_cas2)                               &
     1014     &          -frac*(v_cas(k,it_cas2)-v_cas(k,it_cas1))
     1015        ug_prof_cas(k) = ug_cas(k,it_cas2)                               &
     1016     &          -frac*(ug_cas(k,it_cas2)-ug_cas(k,it_cas1))
     1017        vg_prof_cas(k) = vg_cas(k,it_cas2)                               &
     1018     &          -frac*(vg_cas(k,it_cas2)-vg_cas(k,it_cas1))
     1019        vitw_prof_cas(k) = vitw_cas(k,it_cas2)                               &
     1020     &          -frac*(vitw_cas(k,it_cas2)-vitw_cas(k,it_cas1))
     1021        du_prof_cas(k) = du_cas(k,it_cas2)                                   &
     1022     &          -frac*(du_cas(k,it_cas2)-du_cas(k,it_cas1))
     1023        hu_prof_cas(k) = hu_cas(k,it_cas2)                                   &
     1024     &          -frac*(hu_cas(k,it_cas2)-hu_cas(k,it_cas1))
     1025        vu_prof_cas(k) = vu_cas(k,it_cas2)                                   &
     1026     &          -frac*(vu_cas(k,it_cas2)-vu_cas(k,it_cas1))
     1027        dv_prof_cas(k) = dv_cas(k,it_cas2)                                   &
     1028     &          -frac*(dv_cas(k,it_cas2)-dv_cas(k,it_cas1))
     1029        hv_prof_cas(k) = hv_cas(k,it_cas2)                                   &
     1030     &          -frac*(hv_cas(k,it_cas2)-hv_cas(k,it_cas1))
     1031        vv_prof_cas(k) = vv_cas(k,it_cas2)                                   &
     1032     &          -frac*(vv_cas(k,it_cas2)-vv_cas(k,it_cas1))
     1033        dt_prof_cas(k) = dt_cas(k,it_cas2)                                   &
     1034     &          -frac*(dt_cas(k,it_cas2)-dt_cas(k,it_cas1))
     1035        ht_prof_cas(k) = ht_cas(k,it_cas2)                                   &
     1036     &          -frac*(ht_cas(k,it_cas2)-ht_cas(k,it_cas1))
     1037        vt_prof_cas(k) = vt_cas(k,it_cas2)                                   &
     1038     &          -frac*(vt_cas(k,it_cas2)-vt_cas(k,it_cas1))
     1039        dtrad_prof_cas(k) = dtrad_cas(k,it_cas2)                                   &
     1040     &          -frac*(dtrad_cas(k,it_cas2)-dtrad_cas(k,it_cas1))
     1041        dq_prof_cas(k) = dq_cas(k,it_cas2)                                   &
     1042     &          -frac*(dq_cas(k,it_cas2)-dq_cas(k,it_cas1))
     1043        hq_prof_cas(k) = hq_cas(k,it_cas2)                                   &
     1044     &          -frac*(hq_cas(k,it_cas2)-hq_cas(k,it_cas1))
     1045        vq_prof_cas(k) = vq_cas(k,it_cas2)                                   &
     1046     &          -frac*(vq_cas(k,it_cas2)-vq_cas(k,it_cas1))
    4961047        enddo
    4971048
    … …  
    4991050        END
    5001051
     1052!**********************************************************************************************
Note: See TracChangeset for help on using the changeset viewer.