Ignore:
Timestamp:
Jun 3, 2011, 7:28:17 PM (13 years ago)
Author:
musat
Message:

Ajouts CFMIP2/CMIP5

  • 6eme fichier de sortie "stations" histstn.nc qui necessite 2 fichiers: PARAM/npCFMIP_param.data contenant le nombre de points (120 pour simulations AMIP, 73 pour aqua) PARAM/pointlocations.txt contenat le numero, les coordonnees (lon,lat) et le nom de chaque station
  • flag LOGICAL dans tous les appels histwrite_phy pour pouvoir sortir le fichier histstn.nc

NB: 1) les flags de type phys_ que l'on met dans le physiq.def_L* pour ajouter plus de sorties

necessitent dorenavant 6 valeurs, la 6eme correspondant au fichier histstn.nc

2) par defaut le fichier histstn.nc ne sort pas; pour le sortir ajouter les lignes suivantes

dans physiq.def_L*

### Type de fichier : global (n) ou stations (y)
phys_out_filestations = n n n n n y

  • introduction de 2 jeux de flags pour les taux des GES; taux actuels avec suffixes _act, taux futurs avec "_per" avec 2 appels au rayonnement si taux "_per" different des taux "_act" (utiles pour diags. CFMIP 4CO2)
  • flags "betaCRF" pour calculs CRF pour experiences sensibilite proprietes optiques eau liquide nuageuse avec initialisations par defaut; sinon besoin de fichier beta_crf.data

IM

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/branches/LMDZ4_AR5/libf/phylmd/physiq.F

    r1533 r1534  
    4242      use radlwsw_m, only: radlwsw
    4343
     44!IM stations CFMIP
     45      USE CFMIP_point_locations
    4446      IMPLICIT none
    4547c======================================================================
    … …  
    676678cAA
    677679      REAL coefh(klon,klev)     ! coef d'echange pour phytrac, valable pour 2<=k<=klev
     680      REAL coefm(klon,klev)     ! coef d'echange pour U, V
    678681      REAL u1(klon)             ! vents dans la premiere couche U
    679682      REAL v1(klon)             ! vents dans la premiere couche V
    … …  
    986989      REAL zx_tmp_fi2d(klon)      ! variable temporaire grille physique
    987990      REAL zx_tmp_fi3d(klon,klev) ! variable temporaire pour champs 3D
     991      REAL zx_tmp_fi3d1(klon,klev+1) !variable temporaire pour champs 3D (kelvp1)
    988992c#ifdef histNMC
    989993cym   A voir plus tard !!!!
    … …  
    10161020      REAL ue_lay(klon,klev) ! transport zonal de l'energie a chaque niveau vert.
    10171021      REAL uq_lay(klon,klev) ! transport zonal de l'eau a chaque niveau vert.
    1018 c
    1019 cIM 280405 END
    10201022c
    10211023      INTEGER nhori, nvert, nvert1, nvert3
    … …  
    11591161      REAL grain(1), gtsol(1), gt2m(1), gprw(1)
    11601162
     1163cIM stations CFMIP
     1164      INTEGER, SAVE :: nCFMIP
     1165c$OMP THREADPRIVATE(nCFMIP)
     1166      INTEGER, PARAMETER :: npCFMIP=120
     1167      INTEGER, ALLOCATABLE, SAVE :: tabCFMIP(:)
     1168      REAL, ALLOCATABLE, SAVE :: lonCFMIP(:), latCFMIP(:)
     1169c$OMP THREADPRIVATE(tabCFMIP, lonCFMIP, latCFMIP)
     1170      INTEGER, ALLOCATABLE, SAVE :: tabijGCM(:)
     1171      REAL, ALLOCATABLE, SAVE :: lonGCM(:), latGCM(:)
     1172c$OMP THREADPRIVATE(tabijGCM, lonGCM, latGCM)
     1173      INTEGER, ALLOCATABLE, SAVE :: iGCM(:), jGCM(:)
     1174c$OMP THREADPRIVATE(iGCM, jGCM)
     1175      logical, dimension(nfiles)            :: phys_out_filestations
     1176      logical, parameter :: lNMC=.FALSE.
     1177
     1178cIM betaCRF
     1179      REAL, SAVE :: pfree, beta_pbl, beta_free
     1180c$OMP THREADPRIVATE(pfree, beta_pbl, beta_free)
     1181      REAL, SAVE :: lon1_beta,  lon2_beta, lat1_beta, lat2_beta
     1182c$OMP THREADPRIVATE(lon1_beta,  lon2_beta, lat1_beta, lat2_beta)
     1183      LOGICAL, SAVE :: mskocean_beta
     1184c$OMP THREADPRIVATE(mskocean_beta)
     1185      REAL, dimension(klon, klev) :: beta       ! facteur sur cldtaurad et cldemirad pour evaluer les retros liees aux CRF
     1186      REAL, dimension(klon, klev) :: cldtaurad  ! epaisseur optique pour radlwsw,COSP
     1187      REAL, dimension(klon, klev) :: cldemirad  ! emissivite pour radlwsw,COSP
     1188
    11611189cIM for NMC files
    11621190      missing_val=nf90_fill_real
    … …  
    14431471
    14441472c================================================================================
    1445 
     1473cIM stations CFMIP
     1474      nCFMIP=npCFMIP
     1475      OPEN(98,file='npCFMIP_param.data',status='old',
     1476     $          form='formatted',err=999)
     1477      READ(98,*,end=998) nCFMIP
     1478998   CONTINUE
     1479      CLOSE(98)
     1480999   CONTINUE
     1481      IF(nCFMIP.GT.npCFMIP) THEN
     1482       print*,'nCFMIP > npCFMIP : augmenter npCFMIP et recompiler'
     1483       CALL abort
     1484      else
     1485       print*,'physiq npCFMIP=',npCFMIP,'nCFMIP=',nCFMIP
     1486      ENDIF
     1487c
     1488      ALLOCATE(tabCFMIP(nCFMIP))
     1489      ALLOCATE(lonCFMIP(nCFMIP), latCFMIP(nCFMIP))
     1490      ALLOCATE(tabijGCM(nCFMIP))
     1491      ALLOCATE(lonGCM(nCFMIP), latGCM(nCFMIP))
     1492      ALLOCATE(iGCM(nCFMIP), jGCM(nCFMIP))
     1493c
     1494c lecture des nCFMIP stations CFMIP, de leur numero
     1495c et des coordonnees geographiques lonCFMIP, latCFMIP
     1496c
     1497         CALL read_CFMIP_point_locations(nCFMIP, tabCFMIP,
     1498     $lonCFMIP, latCFMIP)
     1499c
     1500c identification des
     1501c 1) coordonnees lonGCM, latGCM des points CFMIP dans la grille de LMDZ
     1502c 2) indices points tabijGCM de la grille physique 1d sur klon points
     1503c 3) indices iGCM, jGCM de la grille physique 2d
     1504c
     1505         CALL LMDZ_CFMIP_point_locations(nCFMIP, lonCFMIP, latCFMIP,
     1506     $tabijGCM, lonGCM, latGCM, iGCM, jGCM)
     1507c
    14461508         ENDIF !debut
    1447 
     1509 
    14481510           DO i=1,klon
    14491511             rugoro(i) = f_rugoro * MAX(1.0e-05, zstd(i)*zsig(i)/2.0)
    … …  
    14831545     .                   lmt_pas
    14841546c
    1485 cIM 030306 END
    1486 
    14871547      capemaxcels = 't_max(X)'
    14881548      t2mincels = 't_min(X)'
    … …  
    15011561
    15021562c$OMP MASTER
    1503        call phys_output_open(jjmp1,nlevSTD,clevSTD,nbteta,
    1504      &                        ctetaSTD,dtime,ok_veget,
    1505      &                        type_ocean,iflag_pbl,ok_mensuel,ok_journe,
    1506      &                        ok_hf,ok_instan,ok_LES,ok_ade,ok_aie,
    1507      &                        read_climoz, new_aod, aerosol_couple)
     1563       call phys_output_open(rlon,rlat,nCFMIP,tabijGCM,
     1564     &                       iGCM,jGCM,lonGCM,latGCM,
     1565     &                       jjmp1,nlevSTD,clevSTD,
     1566     &                       nbteta, ctetaSTD, dtime,ok_veget,
     1567     &                       type_ocean,iflag_pbl,ok_mensuel,ok_journe,
     1568     &                       ok_hf,ok_instan,ok_LES,ok_ade,ok_aie,
     1569     &                       read_climoz, phys_out_filestations,
     1570     &                       new_aod, aerosol_couple
     1571     &                        )
    15081572c$OMP END MASTER
    15091573c$OMP BARRIER
    … …  
    15251589#endif
    15261590
    1527 cIM 250308bad guide        ecrit_hf2mth = 30*1/ecrit_hf
    15281591         ecrit_hf2mth = ecrit_mth/ecrit_hf
    15291592
    … …  
    15381601         ecrit_reg = ecrit_reg * un_jour
    15391602         ecrit_tra = ecrit_tra * un_jour
    1540          ecrit_ISCCP = ecrit_ISCCP * un_jour
    15411603         ecrit_LES = ecrit_LES * un_jour
    15421604c
    … …  
    15441606     .   ecrit_hf,ecrit_day,ecrit_mth,ecrit_reg,ecrit_tra,ecrit_ISCCP,
    15451607     .   ecrit_hf2mth
    1546 cIM 030306 END
    1547 
    15481608
    15491609cXXXPB Positionner date0 pour initialisation de ORCHIDEE
    … …  
    16031663      END IF
    16041664C$omp end single
     1665c
     1666cIM betaCRF
     1667      pfree=70000. !Pa
     1668      beta_pbl=1.
     1669      beta_free=1.
     1670      lon1_beta=-180.
     1671      lon2_beta=+180.
     1672      lat1_beta=90.
     1673      lat2_beta=-90.
     1674      mskocean_beta=.FALSE.
     1675
     1676      OPEN(99,file='beta_crf.data',status='old',
     1677     $          form='formatted',err=9999)
     1678      READ(99,*,end=9998) pfree
     1679      READ(99,*,end=9998) beta_pbl
     1680      READ(99,*,end=9998) beta_free
     1681      READ(99,*,end=9998) lon1_beta
     1682      READ(99,*,end=9998) lon2_beta
     1683      READ(99,*,end=9998) lat1_beta
     1684      READ(99,*,end=9998) lat2_beta
     1685      READ(99,*,end=9998) mskocean_beta
     16869998  Continue
     1687      CLOSE(99)
     16889999  Continue
     1689      WRITE(*,*)'pfree=',pfree
     1690      WRITE(*,*)'beta_pbl=',beta_pbl
     1691      WRITE(*,*)'beta_free=',beta_free
     1692      WRITE(*,*)'lon1_beta=',lon1_beta
     1693      WRITE(*,*)'lon2_beta=',lon2_beta
     1694      WRITE(*,*)'lat1_beta=',lat1_beta
     1695      WRITE(*,*)'lat2_beta=',lat2_beta
     1696      WRITE(*,*)'mskocean_beta=',mskocean_beta
    16051697      ENDIF
    16061698!
    … …  
    19061998     s     zxtsol,    zxfluxlat, zt2m,    qsat2m,
    19071999     s     d_t_vdf,   d_q_vdf,   d_u_vdf, d_v_vdf,
    1908      s     coefh,     slab_wfbils,               
     2000     s     coefh,     coefm,     slab_wfbils,               
    19092001     d     qsol,      zq2m,      s_pblh,  s_lcl,
    19102002     d     s_capCL,   s_oliqCL,  s_cteiCL,s_pblT,
    … …  
    29333025         mass_solu_aero_pi(:,:) = ccm(:,:,2)
    29343026      END IF
    2935 
     3027c
    29363028      if (ok_newmicro) then
    29373029      CALL newmicro (paprs, pplay,ok_newmicro,
    … …  
    29543046      endif
    29553047c
     3048cIM betaCRF
     3049c
     3050      cldtaurad = cldtau
     3051      cldemirad = cldemi
     3052c
     3053      if(lon1_beta.EQ.-180..AND.lon2_beta.EQ.180..AND.
     3054     $lat1_beta.EQ.90..AND.lat2_beta.EQ.-90.) THEN
     3055c
     3056c global
     3057c
     3058       DO k=1, klev
     3059       DO i=1, klon
     3060        if (pplay(i,k).GE.pfree) THEN
     3061         beta(i,k) = beta_pbl
     3062        else
     3063         beta(i,k) = beta_free
     3064        endif
     3065        if (mskocean_beta) THEN
     3066         beta(i,k) = beta(i,k) * pctsrf(i,is_oce)
     3067        endif
     3068        cldtaurad(i,k) = cldtau(i,k) * beta(i,k)
     3069        cldemirad(i,k) = cldemi(i,k) * beta(i,k)
     3070       ENDDO
     3071       ENDDO
     3072c
     3073      else
     3074c
     3075c regional
     3076c
     3077       DO k=1, klev
     3078       DO i=1,klon
     3079c
     3080        if (rlon(i).ge.lon1_beta.AND.rlon(i).le.lon2_beta.AND.
     3081     $      rlat(i).le.lat1_beta.AND.rlat(i).ge.lat2_beta) THEN
     3082         if (pplay(i,k).GE.pfree) THEN
     3083          beta(i,k) = beta_pbl
     3084         else
     3085          beta(i,k) = beta_free
     3086         endif
     3087         if (mskocean_beta) THEN
     3088          beta(i,k) = beta(i,k) * pctsrf(i,is_oce)
     3089         endif
     3090        cldtaurad(i,k) = cldtau(i,k) * beta(i,k)
     3091        cldemirad(i,k) = cldemi(i,k) * beta(i,k)
     3092        endif
     3093c
     3094       ENDDO
     3095       ENDDO
     3096c
     3097      endif
     3098c
    29563099c Appeler le rayonnement mais calculer tout d'abord l'albedo du sol.
    29573100c
    … …  
    29823125     e        paprs, pplay,zxtsol,albsol1, albsol2, t_seri,q_seri,
    29833126     e        wo(:, :, 1),
    2984      e        cldfra, cldemi, cldtau,
     3127     e        cldfra, cldemirad, cldtaurad,
    29853128     s        heat,heat0,cool,cool0,radsol,albpla,
    29863129     s        topsw,toplw,solsw,sollw,
    … …  
    30003143#endif
    30013144      ELSE
    3002 
     3145c
     3146cIM calcul radiatif pour le cas actuel
     3147c
     3148       RCO2 = RCO2_act
     3149       RCH4 = RCH4_act
     3150       RN2O = RN2O_act
     3151       RCFC11 = RCFC11_act
     3152       RCFC12 = RCFC12_act
     3153c
    30033154         CALL radlwsw
    30043155     e        (dist, rmu0, fract,
    30053156     e        paprs, pplay,zxtsol,albsol1, albsol2,
    30063157     e        t_seri,q_seri,wo,
    3007      e        cldfra, cldemi, cldtau,
     3158     e        cldfra, cldemirad, cldtaurad,
    30083159     e        ok_ade, ok_aie,
    30093160     e        tau_aero, piz_aero, cg_aero,
    … …  
    30233174     o        topswcf_aero, solswcf_aero)
    30243175         
    3025 
     3176c
     3177cIM 2eme calcul radiatif pour le cas perturbe ou au moins un
     3178cIM des taux doit etre different du taux actuel
     3179cIM Par defaut on a les taux perturbes egaux aux taux actuels
     3180c
     3181       if (RCO2_per.NE.RCO2_act.OR.RCH4_per.NE.RCH4_act.OR.
     3182     $RN2O_per.NE.RN2O_act.OR.RCFC11_per.NE.RCFC11_act.OR.
     3183     $RCFC12_per.NE.RCFC12_act) THEN
     3184c
     3185       RCO2 = RCO2_per
     3186       RCH4 = RCH4_per
     3187       RN2O = RN2O_per
     3188       RCFC11 = RCFC11_per
     3189       RCFC12 = RCFC12_per
     3190c
     3191         CALL radlwsw
     3192     e        (dist, rmu0, fract,
     3193     e        paprs, pplay,zxtsol,albsol1, albsol2,
     3194     e        t_seri,q_seri,wo,
     3195     e        cldfra, cldemi, cldtau,
     3196     e        ok_ade, ok_aie,
     3197     e        tau_aero, piz_aero, cg_aero,
     3198     e        cldtaupi,new_aod,
     3199     e        zqsat, flwc, fiwc,
     3200     s        heatp,heat0p,coolp,cool0p,radsolp,albplap,
     3201     s        topswp,toplwp,solswp,sollwp,
     3202     s        sollwdownp,
     3203     s        topsw0p,toplw0p,solsw0p,sollw0p,
     3204     s        lwdn0p, lwdnp, lwup0p, lwupp,
     3205     s        swdn0p, swdnp, swup0p, swupp,
     3206     s        topswad_aerop, solswad_aerop,
     3207     s        topswai_aerop, solswai_aerop,
     3208     o        topswad0_aerop, solswad0_aerop,
     3209     o        topsw_aerop, topsw0_aerop,
     3210     o        solsw_aerop, solsw0_aerop,
     3211     o        topswcf_aerop, solswcf_aerop)
     3212       endif
     3213c
    30263214      ENDIF ! aerosol_couple
    30273215      itaprad = 0
    … …  
    31843372c
    31853373c  ajout des tendances
    3186         CALL add_phys_tend(d_u_hin,d_v_hin,d_t_hin,dq0,dql0,'lif')
     3374        CALL add_phys_tend(d_u_hin,d_v_hin,d_t_hin,dq0,dql0,'hin')
    31873375
    31883376      ENDIF
    … …  
    32603448     $                   prfl(:,1:klev),psfl(:,1:klev),
    32613449     $                   pmflxr(:,1:klev),pmflxs(:,1:klev),
    3262      $                   mr_ozone,cldtau, cldemi)
     3450     $                   mr_ozone,cldtaurad, cldemirad)
    32633451
    32643452!     L          calipso2D,calipso3D,cfadlidar,parasolrefl,atb,betamol,
    … …  
    34163604      vwriteSTD(:,:,5)=vlevSTD(:,:)
    34173605      wwriteSTD(:,:,5)=wlevSTD(:,:)
     3606c
     3607cIM initialisation 6eme fichier de sortie
     3608      twriteSTD(:,:,6)=tlevSTD(:,:)
     3609      qwriteSTD(:,:,6)=qlevSTD(:,:)
     3610      rhwriteSTD(:,:,6)=rhlevSTD(:,:)
     3611      phiwriteSTD(:,:,6)=philevSTD(:,:)
     3612      uwriteSTD(:,:,6)=ulevSTD(:,:)
     3613      vwriteSTD(:,:,6)=vlevSTD(:,:)
     3614      wwriteSTD(:,:,6)=wlevSTD(:,:)
    34183615cIM for NMC files
    34193616      DO n=1, nlevSTD3
Note: See TracChangeset for help on using the changeset viewer.