Ignore:
Timestamp:
Nov 19, 2021, 4:58:59 PM (3 years ago)
Author:
lguez
Message:

Sync latest trunk changes to Ocean_skin

Location:
LMDZ6/branches/Ocean_skin
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/branches/Ocean_skin

  • LMDZ6/branches/Ocean_skin/libf/phylmd/radlwsw_m.F90

    r3798 r4013  
    1616   t,q,wo,&
    1717   cldfra, cldemi, cldtaupd,&
    18    ok_ade, ok_aie, ok_volcan, flag_aerosol,&
     18   ok_ade, ok_aie, ok_volcan, flag_volc_surfstrat, flag_aerosol,&
    1919   flag_aerosol_strat, flag_aer_feedback, &
    2020   tau_aero, piz_aero, cg_aero,&
    21    tau_aero_sw_rrtm, piz_aero_sw_rrtm, cg_aero_sw_rrtm,& ! rajoute par OB pour RRTM
    22    tau_aero_lw_rrtm, &                                   ! rajoute par C. Kleinschmitt pour RRTM
     21   tau_aero_sw_rrtm, piz_aero_sw_rrtm, cg_aero_sw_rrtm,& ! rajoute par OB RRTM
     22   tau_aero_lw_rrtm, &              ! rajoute par C.Kleinschmitt pour RRTM
    2323   cldtaupi, &
    2424   qsat, flwc, fiwc, &
    … …  
    4545   ZSWFT0_i, ZFSDN0, ZFSUP0)
    4646
    47 
    48 
     47! Modules necessaires
    4948  USE DIMPHY
    5049  USE assert_m, ONLY : assert
    5150  USE infotrac_phy, ONLY : type_trac
    5251  USE write_field_phy
     52
    5353#ifdef REPROBUS
    5454  USE CHEM_REP, ONLY : solaireTIME, ok_SUNTIME, ndimozon
    5555#endif
     56
    5657#ifdef CPP_RRTM
    5758!    modules necessaires au rayonnement
    5859!    -----------------------------------------
    59 !     USE YOMCST   , ONLY : RG       ,RD       ,RTT      ,RPI
    60 !     USE YOERAD   , ONLY : NSW      ,LRRTM    ,LINHOM   , LCCNL,LCCNO,
    61 !     USE YOERAD   , ONLY : NSW      ,LRRTM    ,LCCNL    ,LCCNO ,&
    62 ! NSW mis dans .def MPL 20140211
    63 ! NLW ajoute par OB
    6460      USE YOERAD   , ONLY : NLW, LRRTM    ,LCCNL    ,LCCNO ,&
    6561          NRADIP   , NRADLP , NICEOPT, NLIQOPT ,RCCNLND  , RCCNSEA
    … …  
    7369          RFLDD1   ,RFLDD2   ,RFLDD3   ,RFUETA   ,RASWCA,&
    7470          RASWCB   ,RASWCC   ,RASWCD   ,RASWCE   ,RASWCF
    75 !    &    RASWCB   ,RASWCC   ,RASWCD   ,RASWCE   ,RASWCF, RLINLI
    7671      USE YOERDU   , ONLY : NUAER  ,NTRAER ,REPLOG ,REPSC  ,REPSCW ,DIFF
    77 !      USE YOETHF   , ONLY : RTICE
    7872      USE YOERRTWN , ONLY : DELWAVE   ,TOTPLNK     
    7973      USE YOMPHY3  , ONLY : RII0
    … …  
    8175      USE aero_mod
    8276
     77! AI 02.2021
     78! Besoin pour ECRAD de pctsrf, zmasq, longitude, altitude
     79#ifdef CPP_ECRAD
     80      USE geometry_mod, ONLY: latitude, longitude
     81      USE phys_state_var_mod, ONLY: pctsrf
     82      USE indice_sol_mod
     83      USE time_phylmdz_mod, only: current_time
     84      USE phys_cal_mod, only: day_cur
     85#endif
     86
    8387  !======================================================================
    8488  ! Auteur(s): Z.X. Li (LMD/CNRS) date: 19960719
    8589  ! Objet: interface entre le modele et les rayonnements
    8690  ! Arguments:
    87   ! dist-----input-R- distance astronomique terre-soleil
    88   ! rmu0-----input-R- cosinus de l'angle zenithal
    89   ! fract----input-R- duree d'ensoleillement normalisee
    90   ! co2_ppm--input-R- concentration du gaz carbonique (en ppm)
    91   ! paprs----input-R- pression a inter-couche (Pa)
    92   ! pplay----input-R- pression au milieu de couche (Pa)
    93   ! tsol-----input-R- temperature du sol (en K)
    94   ! alb1-----input-R- albedo du sol(entre 0 et 1) dans l'interval visible
    95   ! alb2-----input-R- albedo du sol(entre 0 et 1) dans l'interval proche infra-rouge   
    96   ! t--------input-R- temperature (K)
    97   ! q--------input-R- vapeur d'eau (en kg/kg)
    98   ! cldfra---input-R- fraction nuageuse (entre 0 et 1)
    99   ! cldtaupd---input-R- epaisseur optique des nuages dans le visible (present-day value)
    100   ! cldemi---input-R- emissivite des nuages dans l'IR (entre 0 et 1)
    101   ! ok_ade---input-L- apply the Aerosol Direct Effect or not?
    102   ! ok_aie---input-L- apply the Aerosol Indirect Effect or not?
    103   ! ok_volcan-input-L- activate volcanic diags (SW heat & LW cool rate, SW & LW flux)
    104   ! flag_aerosol-input-I- aerosol flag from 0 to 6
    105   ! flag_aerosol_strat-input-I- use stratospheric aerosols flag (0, 1, 2)
    106   ! flag_aer_feedback-input-I- activate aerosol radiative feedback (T, F)
    107   ! tau_ae, piz_ae, cg_ae-input-R- aerosol optical properties (calculated in aeropt.F)
    108   ! cldtaupi-input-R- epaisseur optique des nuages dans le visible
     91  !                  INPUTS
     92  ! dist----- input-R- distance astronomique terre-soleil
     93  ! rmu0----- input-R- cosinus de l'angle zenithal
     94  ! fract---- input-R- duree d'ensoleillement normalisee
     95  ! co2_ppm-- input-R- concentration du gaz carbonique (en ppm)
     96  ! paprs---- input-R- pression a inter-couche (Pa)
     97  ! pplay---- input-R- pression au milieu de couche (Pa)
     98  ! tsol----- input-R- temperature du sol (en K)
     99  ! alb1----- input-R- albedo du sol(entre 0 et 1) dans l'interval visible
     100  ! alb2----- input-R- albedo du sol(entre 0 et 1) dans l'interval proche infra-rouge   
     101  ! t-------- input-R- temperature (K)
     102  ! q-------- input-R- vapeur d'eau (en kg/kg)
     103  ! cldfra--- input-R- fraction nuageuse (entre 0 et 1)
     104  ! cldtaupd- input-R- epaisseur optique des nuages dans le visible (present-day value)
     105  ! cldemi--- input-R- emissivite des nuages dans l'IR (entre 0 et 1)
     106  ! ok_ade--- input-L- apply the Aerosol Direct Effect or not?
     107  ! ok_aie--- input-L- apply the Aerosol Indirect Effect or not?
     108  ! ok_volcan input-L- activate volcanic diags (SW heat & LW cool rate, SW & LW flux)
     109  ! flag_volc_surfstrat input-I- activate volcanic surf cooling or strato heating (or nothing)
     110  ! flag_aerosol input-I- aerosol flag from 0 to 6
     111  ! flag_aerosol_strat input-I- use stratospheric aerosols flag (0, 1, 2)
     112  ! flag_aer_feedback  input-I- activate aerosol radiative feedback (T, F)
     113  ! tau_ae, piz_ae, cg_ae input-R- aerosol optical properties (calculated in aeropt.F)
     114  ! cldtaupi  input-R- epaisseur optique des nuages dans le visible
    109115  !                   calculated for pre-industrial (pi) aerosol concentrations, i.e. with smaller
    110116  !                   droplet concentration, thus larger droplets, thus generally cdltaupi cldtaupd
    111117  !                   it is needed for the diagnostics of the aerosol indirect radiative forcing     
    112118  !
     119  !                  OUTPUTS
    113120  ! heat-----output-R- echauffement atmospherique (visible) (K/jour)
    114121  ! cool-----output-R- refroidissement dans l'IR (K/jour)
    … …  
    177184  !
    178185  ! ====================================================================
     186
     187! ==============
     188! DECLARATIONS
     189! ==============
    179190  include "YOETHF.h"
    180191  include "YOMCST.h"
    … …  
    200211  LOGICAL, INTENT(in)  :: ok_ade, ok_aie                                 ! switches whether to use aerosol direct (indirect) effects or not
    201212  LOGICAL, INTENT(in)  :: ok_volcan                                      ! produce volcanic diags (SW/LW heat flux and rate)
    202   LOGICAL              :: lldebug
     213  INTEGER, INTENT(in)  :: flag_volc_surfstrat                            ! allow to impose volcanic cooling rate at surf or heating in strato
     214  LOGICAL              :: lldebug=.false.
    203215  INTEGER, INTENT(in)  :: flag_aerosol                                   ! takes value 0 (no aerosol) or 1 to 6 (aerosols)
    204216  INTEGER, INTENT(in)  :: flag_aerosol_strat                             ! use stratospheric aerosols
    … …  
    286298  REAL(KIND=8) PTAVE(kdlon,kflev)
    287299  REAL(KIND=8) PWV(kdlon,kflev), PQS(kdlon,kflev)
     300
     301!!!!!!! Declarations specifiques pour ECRAD !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     302! AI 02.2021
     303#ifdef CPP_ECRAD
     304! ATTENTION les dimensions klon, kdlon ???
     305! INPUTS
     306  REAL, DIMENSION(kdlon,kflev+1) :: ZSWFT0_ii, ZLWFT0_ii
     307  REAL(KIND=8) ZEMISW(klon), &              ! LW emissivity inside the window region
     308               ZEMIS(klon)                  ! LW emissivity outside the window region
     309  REAL(KIND=8) ZGELAM(klon), &              ! longitudes en rad
     310               ZGEMU(klon)                  ! sin(latitude)
     311  REAL(KIND=8) ZCO2(klon,klev), &           ! CO2 mass mixing ratios on full levels
     312               ZCH4(klon,klev), &           ! CH4 mass mixing ratios on full levels
     313               ZN2O(klon,klev), &           ! N2O mass mixing ratios on full levels
     314               ZNO2(klon,klev), &           ! NO2 mass mixing ratios on full levels
     315               ZCFC11(klon,klev), &         ! CFC11
     316               ZCFC12(klon,klev), &         ! CFC12
     317               ZHCFC22(klon,klev), &        ! HCFC22
     318               ZCCL4(klon,klev)           ! CCL4
     319!               ZO3_DP(klon,klev), ZO3_DP_i(klon,klev)            ! Ozone
     320  REAL(KIND=8) ZQ_RAIN(klon,klev), &        ! Rain cloud mass mixing ratio (kg/kg) ?
     321               ZQ_SNOW(klon,klev)           ! Snow cloud mass mixing ratio (kg/kg) ?
     322  REAL(KIND=8) ZAEROSOL_OLD(KLON,6,KLEV), &  !
     323               ZAEROSOL(KLON,KLEV,naero_tot) !
     324! OUTPUTS
     325  REAL(KIND=8) ZFLUX_DIR(klon), &           ! Direct compt of surf flux into horizontal plane
     326               ZFLUX_DIR_CLEAR(klon), &     ! CS Direct
     327               ZFLUX_DIR_INTO_SUN(klon), &  !
     328               ZFLUX_UV(klon), &            ! UV flux
     329               ZFLUX_PAR(klon), &           ! photosynthetically active radiation similarly
     330               ZFLUX_PAR_CLEAR(klon), &     ! CS photosynthetically
     331               ZFLUX_SW_DN_TOA(klon), &     ! DN SW flux at TOA
     332               ZEMIS_OUT(klon)              ! effective broadband emissivity
     333  REAL(KIND=8) ZLWDERIVATIVE(klon,klev+1)   ! LW derivatives
     334  REAL(KIND=8) ZSWDIFFUSEBAND(klon,NSW), &  ! SW DN flux in diffuse albedo band
     335               ZSWDIRECTBAND(klon,NSW)      ! SW DN flux in direct albedo band
     336#endif
     337!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    288338
    289339  REAL(kind=8) POZON(kdlon, kflev, size(wo, 3)) ! mass fraction of ozone
    … …  
    317367  REAL(KIND=8) ztopswaiaero(kdlon), zsolswaiaero(kdlon)     ! dito, indirect
    318368!--NL
    319   REAL(KIND=8) zswadaero(kdlon,kflev+1)                       ! SW Aerosol direct forcing
    320   REAL(KIND=8) zlwadaero(kdlon,kflev+1)                       ! LW Aerosol direct forcing
     369  REAL(KIND=8) zswadaero(kdlon,kflev+1)                     ! SW Aerosol direct forcing
     370  REAL(KIND=8) zlwadaero(kdlon,kflev+1)                     ! LW Aerosol direct forcing
     371  REAL(KIND=8) volmip_solsw(kdlon)                          ! SW clear sky in the case of VOLMIP
    321372!-LW by CK
    322373  REAL(KIND=8) ztoplwadaero(kdlon), zsollwadaero(kdlon)     ! LW Aerosol direct forcing at TOAand surface
    … …  
    401452  REAL zdir, zdif
    402453
     454! =========  INITIALISATIONS ==============================================
     455 IF (lldebug) THEN
     456  print*,'Entree dans radlwsw '
     457  print*,'************* INITIALISATIONS *****************************'
     458  print*,'klon, kdlon, klev, kflev =',klon, kdlon, klev, kflev
     459 ENDIF
     460
    403461  CALL assert(size(wo, 1) == klon, size(wo, 2) == klev, "radlwsw wo")
    404   ! initialisation
     462 
    405463  ist=1
    406464  iend=klon
    407465  ktdia=1
    408466  kmode=ist
     467! Aeros
    409468  tauaero(:,:,:,:)=0.
    410469  pizaero(:,:,:,:)=0.
    411470  cgaero(:,:,:,:)=0.
    412   lldebug=.FALSE.
     471!  lldebug=.FALSE.
    413472
    414473  ztopsw_aero(:,:)  = 0. !ym missing init : warning : not initialized in SW_AEROAR4
    … …  
    462521  ENDIF
    463522
     523 IF (lldebug) THEN
     524  print*,'************** Debut boucle de 1 a ', nb_gr
     525 ENDIF
     526
    464527  DO j = 1, nb_gr
    465528    iof = kdlon*(j-1)
    466529    DO i = 1, kdlon
    467530      zfract(i) = fract(iof+i)
    468 !     zfract(i) = 1.     !!!!!!  essai MPL 19052010
    469531      zrmu0(i) = rmu0(iof+i)
    470532
    471533
    472 !albedo SB >>>
    473 !
    474534      IF (iflag_rrtm==0) THEN
    475 !
     535!     Albedo
    476536        PALBD(i,1)=alb_dif(iof+i,1)
    477537        PALBD(i,2)=alb_dif(iof+i,2)
    478538        PALBP(i,1)=alb_dir(iof+i,1)
    479539        PALBP(i,2)=alb_dir(iof+i,2)
    480 !
    481       ELSEIF (iflag_rrtm==1) THEn
    482 !
     540! AI 02.2021 cas iflag_rrtm=1 et 2
     541       ELSEIF (iflag_rrtm==1.OR.iflag_rrtm==2) THEN
    483542        DO kk=1,NSW
    484543          PALBD_NEW(i,kk)=alb_dif(iof+i,kk)
    … …  
    488547      ENDIF
    489548!albedo SB <<<
    490 
    491549
    492550      PEMIS(i) = 1.0    !!!!! A REVOIR (MPL)
    … …  
    569627      ENDDO
    570628    ENDDO
     629!
     630! AI 02.2021
     631#ifdef CPP_ECRAD
     632  ZEMIS = 1.0
     633  ZEMISW = 1.0
     634  ZGELAM = longitude
     635  ZGEMU = sin(latitude)
     636  ZCO2 = RCO2
     637  ZCH4 = RCH4
     638  ZN2O = RN2O
     639  ZNO2 = 0.0
     640  ZCFC11 = RCFC11
     641  ZCFC12 = RCFC12
     642  ZHCFC22 = 0.0
     643  ZCCL4 = 0.0
     644  ZQ_RAIN = 0.0
     645  ZQ_SNOW = 0.0
     646  ZAEROSOL_OLD = 0.0
     647  ZAEROSOL = 0.0
     648#endif
    571649!
    572650!===== iflag_rrtm ================================================
    … …  
    693771       ENDDO 
    694772!
    695     ELSE
     773    ELSE IF (iflag_rrtm == 1) then
    696774#ifdef CPP_RRTM
    697775!      if (prt_level.gt.10)write(lunout,*)'CPP_RRTM=.T.'
    … …  
    804882            ENDDO
    805883         ENDDO
     884
    806885!       print *,'RADLWSW: avant RECMWFL, RI0,rmu0=',solaire,rmu0
    807886
    … …  
    819898! RII0 = RIP0M15 ! =rip0m if Morcrette non-each time step call.
    820899         RII0=solaire/zdist/zdist
    821 !print*,'+++ radlwsw: solaire ,RII0',solaire,RII0
    822900!  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    823901! Ancien appel a RECMWF (celui du cy25)
    … …  
    852930         PALBD_NEW,PALBP_NEW, paprs_i , pplay_i , RCO2   , cldfra_i,&
    853931         POZON_i  , PAER_i  , PDP_i   , PEMIS   , rmu0   ,&
    854           q_i     , qsat_i  , fiwc_i  , flwc_i  , zmasq  , t_i  ,tsol,&
     932         q_i     , qsat_i  , fiwc_i  , flwc_i  , zmasq  , t_i  ,tsol,&
    855933         ref_liq_i, ref_ice_i, &
    856934         ref_liq_pi_i, ref_ice_pi_i, &   ! rajoute par OB pour diagnostiquer effet indirect
    … …  
    873951         ZTOPLWAIAERO,ZSOLLWAIAERO, &
    874952         ZLWADAERO, & !--NL
     953         volmip_solsw, flag_volc_surfstrat, & !--VOLMIP
    875954         ok_ade, ok_aie, ok_volcan, flag_aerosol,flag_aerosol_strat, flag_aer_feedback) ! flags aerosols
     955
     956!--OB diagnostics
     957! & PTOPSWAIAERO,PSOLSWAIAERO,&
     958! & PTOPSWCFAERO,PSOLSWCFAERO,&
     959! & PSWADAERO,& !--NL
     960!!--LW diagnostics CK
     961! & PTOPLWADAERO,PSOLLWADAERO,&
     962! & PTOPLWAD0AERO,PSOLLWAD0AERO,&
     963! & PTOPLWAIAERO,PSOLLWAIAERO,&
     964! & PLWADAERO,& !--NL
     965!!..end
     966! & ok_ade, ok_aie, ok_volcan, flag_aerosol,flag_aerosol_strat,&
     967! & flag_aer_feedback)
     968
    876969           
    877970!        print *,'RADLWSW: apres RECMWF'
    … …  
    902995        CALL writefield_phy('zfcup_i',ZFCUP_i,klev+1)
    903996      ENDIF
    904 ! --------- output RECMWFL
    905 !  ZEMTD        (KPROMA,KLEV+1)  ; TOTAL DOWNWARD LONGWAVE EMISSIVITY
    906 !  ZEMTU        (KPROMA,KLEV+1)  ; TOTAL UPWARD   LONGWAVE EMISSIVITY
    907 !  ZTRSO        (KPROMA,KLEV+1)  ; TOTAL SHORTWAVE TRANSMISSIVITY
    908 !  ZTH          (KPROMA,KLEV+1)  ; HALF LEVEL TEMPERATURE
    909 !  ZCTRSO       (KPROMA,2)       ; CLEAR-SKY SHORTWAVE TRANSMISSIVITY
    910 !  ZCEMTR       (KPROMA,2)       ; CLEAR-SKY NET LONGWAVE EMISSIVITY
    911 !  ZTRSOD       (KPROMA)         ; TOTAL-SKY SURFACE SW TRANSMISSITY
    912 !  ZLWFC        (KPROMA,2)       ; CLEAR-SKY LONGWAVE FLUXES
    913 !  ZLWFT        (KPROMA,KLEV+1)  ; TOTAL-SKY LONGWAVE FLUXES
    914 !  ZSWFC        (KPROMA,2)       ; CLEAR-SKY SHORTWAVE FLUXES
    915 !  ZSWFT        (KPROMA,KLEV+1)  ; TOTAL-SKY SHORTWAVE FLUXES
    916 !  PPIZA_TOT    (KPROMA,KLEV,NSW); Single scattering albedo of total aerosols
    917 !  PCGA_TOT     (KPROMA,KLEV,NSW); Assymetry factor for total aerosols
    918 !  PTAU_TOT     (KPROMA,KLEV,NSW); Optical depth of total aerosols
    919 !  PPIZA_NAT    (KPROMA,KLEV,NSW); Single scattering albedo of natural aerosols
    920 !  PCGA_NAT     (KPROMA,KLEV,NSW); Assymetry factor for natural aerosols
    921 !  PTAU_NAT     (KPROMA,KLEV,NSW); Optical depth of natiral aerosols
    922 !  PTAU_LW_TOT  (KPROMA,KLEV,NLW); LW Optical depth of total aerosols 
    923 !  PTAU_LW_NAT  (KPROMA,KLEV,NLW); LW Optical depth of natural aerosols 
    924 !  PSFSWDIR     (KPROMA,NSW)     ;
    925 !  PSFSWDIF     (KPROMA,NSW)     ;
    926 !  PFSDNN       (KPROMA)         ;
    927 !  PFSDNV       (KPROMA)         ;
     997
    928998! ---------
    929999! ---------
    … …  
    9831053      ZSOLSWCF_AERO(:,3)=ZSOLSWCF_AERO(:,3)*fract(:)
    9841054
    985 !     print*,'SW_RRTM ZFSDN0 1 , klev:',ZFSDN0(1:klon,1),ZFSDN0(1:klon,klev)
    986 !     print*,'SW_RRTM ZFSUP0 1 , klev:',ZFSUP0(1:klon,1),ZFSUP0(1:klon,klev)
    987 !     print*,'SW_RRTM ZFSDN  1 , klev:',ZFSDN(1:klon,1),ZFSDN(1:klon,klev)
    988 !     print*,'SW_RRTM ZFSUP  1 , klev:',ZFSUP(1:klon,1),ZFSUP(1:klon,klev)     
    989 !     print*,'OK1'
    9901055! ---------
    9911056! ---------
    … …  
    10341099!     print*,'OK2'
    10351100
     1101!--add VOLMIP (surf cool or strat heat activate)
     1102      IF (flag_volc_surfstrat > 0) THEN
     1103         DO i = 1, kdlon
     1104            zsolsw(i)    = volmip_solsw(i)*fract(i)
     1105         ENDDO
     1106      ENDIF
     1107
    10361108! extrait de SW_AR4
    10371109!     DO k = 1, KFLEV
    … …  
    10611133    call abort_physic(modname, abort_message, 1)
    10621134#endif
    1063     ENDIF ! iflag_rrtm
     1135!======================================================================
     1136! AI fev 2021
     1137    ELSE IF(iflag_rrtm == 2) THEN
     1138    print*,'Traitement cas iflag_rrtm = ',iflag_rrtm
     1139!    print*,'Mise a zero des flux '
     1140#ifdef CPP_ECRAD
     1141      DO k = 1, kflev+1
     1142      DO i = 1, kdlon
     1143        ZEMTD_i(i,k)=0.
     1144        ZEMTU_i(i,k)=0.
     1145        ZTRSO_i(i,k)=0.
     1146        ZTH_i(i,k)=0.
     1147        ZLWFT_i(i,k)=0.
     1148        ZSWFT_i(i,k)=0.
     1149        ZFLUX_i(i,1,k)=0.
     1150        ZFLUX_i(i,2,k)=0.
     1151        ZFLUC_i(i,1,k)=0.
     1152        ZFLUC_i(i,2,k)=0.
     1153        ZFSDWN_i(i,k)=0.
     1154        ZFCDWN_i(i,k)=0.
     1155        ZFCCDWN_i(i,k)=0.
     1156        ZFSUP_i(i,k)=0.
     1157        ZFCUP_i(i,k)=0.
     1158        ZFCCUP_i(i,k)=0.
     1159        ZFLCCDWN_i(i,k)=0.
     1160        ZFLCCUP_i(i,k)=0.
     1161      ENDDO
     1162      ENDDO
     1163!
     1164! AI ATTENTION Aerosols A REVOIR
     1165!      DO i = 1, kdlon
     1166!      DO k = 1, kflev
     1167!      DO kk=1, NSW
     1168!
     1169!      PTAU_TOT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,2,kk)
     1170!      PPIZA_TOT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,2,kk)
     1171!      PCGA_TOT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,2,kk)
     1172!
     1173!      PTAU_NAT(i,kflev+1-k,kk)=tau_aero_sw_rrtm(i,k,1,kk)
     1174!      PPIZA_NAT(i,kflev+1-k,kk)=piz_aero_sw_rrtm(i,k,1,kk)
     1175!      PCGA_NAT(i,kflev+1-k,kk)=cg_aero_sw_rrtm(i,k,1,kk)
     1176!
     1177!      ENDDO
     1178!      ENDDO
     1179!      ENDDO
     1180!-end OB
     1181!
     1182!      DO i = 1, kdlon
     1183!      DO k = 1, kflev
     1184!      DO kk=1, NLW
     1185!
     1186!      PTAU_LW_TOT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,2,kk)
     1187!      PTAU_LW_NAT(i,kflev+1-k,kk)=tau_aero_lw_rrtm(i,k,1,kk)
     1188!
     1189!      ENDDO
     1190!      ENDDO
     1191!      ENDDO
     1192!-end C. Kleinschmitt
     1193!     
     1194      DO i = 1, kdlon
     1195      ZCTRSO(i,1)=0.
     1196      ZCTRSO(i,2)=0.
     1197      ZCEMTR(i,1)=0.
     1198      ZCEMTR(i,2)=0.
     1199      ZTRSOD(i)=0.
     1200      ZLWFC(i,1)=0.
     1201      ZLWFC(i,2)=0.
     1202      ZSWFC(i,1)=0.
     1203      ZSWFC(i,2)=0.
     1204      PFSDNN(i)=0.
     1205      PFSDNV(i)=0.
     1206      DO kk = 1, NSW
     1207        PSFSWDIR(i,kk)=0.
     1208        PSFSWDIF(i,kk)=0.
     1209      ENDDO
     1210      ENDDO
     1211!----- Fin des mises a zero des tableaux output -------------------             
     1212
     1213! On met les donnees dans l'ordre des niveaux ecrad
     1214!         print*,'On inverse sur la verticale '
     1215         paprs_i(:,1)=paprs(:,klev+1)
     1216         DO k=1,klev
     1217            paprs_i(1:klon,k+1) =paprs(1:klon,klev+1-k)
     1218            pplay_i(1:klon,k)   =pplay(1:klon,klev+1-k)
     1219            cldfra_i(1:klon,k)  =cldfra(1:klon,klev+1-k)
     1220            PDP_i(1:klon,k)     =PDP(1:klon,klev+1-k)
     1221            t_i(1:klon,k)       =t(1:klon,klev+1-k)
     1222            q_i(1:klon,k)       =q(1:klon,klev+1-k)
     1223            qsat_i(1:klon,k)    =qsat(1:klon,klev+1-k)
     1224            flwc_i(1:klon,k)    =flwc(1:klon,klev+1-k)
     1225            fiwc_i(1:klon,k)    =fiwc(1:klon,klev+1-k)
     1226            ref_liq_i(1:klon,k) =ref_liq(1:klon,klev+1-k)
     1227            ref_ice_i(1:klon,k) =ref_ice(1:klon,klev+1-k)
     1228!-OB
     1229            ref_liq_pi_i(1:klon,k) =ref_liq_pi(1:klon,klev+1-k)
     1230            ref_ice_pi_i(1:klon,k) =ref_ice_pi(1:klon,klev+1-k)
     1231         ENDDO
     1232         DO k=1,kflev
     1233            POZON_i(1:klon,k,:)=POZON(1:klon,kflev+1-k,:)
     1234!            ZO3_DP_i(1:klon,k)=ZO3_DP(1:klon,kflev+1-k)
     1235!            DO i=1,6
     1236            PAER_i(1:klon,k,:)=PAER(1:klon,kflev+1-k,:)
     1237!            ENDDO
     1238         ENDDO
     1239! AI 02.2021
     1240! Calcul de ZTH_i (temp aux interfaces 1:klev+1)
     1241         DO K=2,KLEV
     1242            ZTH_i(:,K)=&
     1243              & (t_i(:,K-1)*pplay_i(:,K-1)*(pplay_i(:,K)-paprs_i(:,K))&
     1244              & +t_i(:,K)*pplay_i(:,K)*(paprs_i(:,K)-pplay_i(:,K-1)))&
     1245              & *(1.0/(paprs_i(:,K)*(pplay_i(:,K)-pplay_i(:,K-1))))
     1246         ENDDO
     1247            ZTH_i(:,KLEV+1)=tsol(:)
     1248            ZTH_i(:,1)=t_i(:,1)-pplay_i(:,1)*(t_i(:,1)-ZTH_i(:,2))&
     1249                      & /(pplay_i(:,1)-paprs_i(:,2))
     1250
     1251      print *,'RADLWSW: avant RADIATION_SCHEME '
     1252      IF (lldebug) THEN
     1253        CALL writefield_phy('rmu0',rmu0,1)
     1254        CALL writefield_phy('tsol',tsol,1)
     1255        CALL writefield_phy('emissiv_out',ZEMIS,1)
     1256        CALL writefield_phy('emissiv_in',ZEMISW,1)
     1257        CALL writefield_phy('pctsrf_ter',pctsrf(:,is_ter),1)
     1258        CALL writefield_phy('pctsrf_oce',pctsrf(:,is_oce),1)
     1259        CALL writefield_phy('ZGELAM',ZGELAM,1)
     1260        CALL writefield_phy('ZGEMU',ZGEMU,1)
     1261        CALL writefield_phy('zmasq',zmasq,1)
     1262        CALL writefield_phy('paprs_i',paprs_i,klev+1)
     1263        CALL writefield_phy('pplay_i',pplay_i,klev)
     1264        CALL writefield_phy('t_i',t_i,klev)
     1265        CALL writefield_phy('ZTH_i',ZTH_i,klev+1)
     1266        CALL writefield_phy('cldfra_i',cldfra_i,klev)
     1267        CALL writefield_phy('paer_i',PAER_i,klev)
     1268        CALL writefield_phy('q_i',q_i,klev)
     1269        CALL writefield_phy('fiwc_i',fiwc_i,klev)
     1270        CALL writefield_phy('flwc_i',flwc_i,klev)
     1271        CALL writefield_phy('palbd_new',PALBD_NEW,NSW)
     1272        CALL writefield_phy('palbp_new',PALBP_NEW,NSW)
     1273!        CALL writefield_phy('ZO3_DP',ZO3_DP,klev)
     1274      ENDIF
     1275 
     1276      CALL RADIATION_SCHEME &
     1277      & (ist, iend, klon, klev, naero_tot, NSW, &
     1278! ??? naero_tot
     1279      & day_cur, current_time, &
     1280!      & solaire, &
     1281      & PSCT, &
     1282      & rmu0, tsol, PALBD_NEW,PALBP_NEW, &   
     1283!       PEMIS_WINDOW (???), &
     1284      & ZEMIS, ZEMISW, &
     1285!       PCCN_LAND, PCCN_SEA, & ???
     1286      & pctsrf(:,is_ter), pctsrf(:,is_oce), &
     1287!       longitude(rad), sin(latitude), PMASQ_ ???
     1288      & ZGELAM, ZGEMU, zmasq, &
     1289!       pression et temp aux milieux
     1290      & pplay_i, t_i, &
     1291!       PTEMPERATURE_H ?,
     1292      & paprs_i, ZTH_i, q_i, qsat_i, &
     1293!       Gas
     1294       & ZCO2, ZCH4, ZN2O, ZNO2, ZCFC11, ZCFC12, ZHCFC22, ZCCL4, POZON_i(:,:,1), &
     1295!       nuages :
     1296      & cldfra_i, flwc_i, fiwc_i, ZQ_RAIN, ZQ_SNOW, & 
     1297      & ref_liq_i, ref_ice_i, &
     1298!       aerosols
     1299      & ZAEROSOL_OLD, ZAEROSOL, &
     1300! Outputs
     1301!       Net flux :
     1302      & ZSWFT_i, ZLWFT_i, ZSWFT0_ii, ZLWFT0_ii, &
     1303!       DWN flux :
     1304      & ZFSDWN_i, ZFLUX_i(:,2,:), ZFCDWN_i, ZFLUC_i(:,2,:), &
     1305!       UP flux :
     1306      & ZFSUP_i, ZFLUX_i(:,1,:), ZFCUP_i, ZFLUC_i(:,1,:), &
     1307!       Surf Direct flux : ATTENTION
     1308      & ZFLUX_DIR, ZFLUX_DIR_CLEAR, ZFLUX_DIR_INTO_SUN, &
     1309!       UV and para flux
     1310      & ZFLUX_UV, ZFLUX_PAR, ZFLUX_PAR_CLEAR, &
     1311!      & ZFLUX_SW_DN_TOA,
     1312      & ZEMIS_OUT, ZLWDERIVATIVE, &
     1313      & PSFSWDIF, PSFSWDIR)
     1314
     1315      print *,'========= RADLWSW: apres RADIATION_SCHEME ==================== '
     1316
     1317      IF (lldebug) THEN
     1318        CALL writefield_phy('zlwft_i',ZLWFT_i,klev+1)
     1319        CALL writefield_phy('zlwft0_ii',ZLWFT0_ii,klev+1)
     1320        CALL writefield_phy('zswft_i',ZSWFT_i,klev+1)
     1321        CALL writefield_phy('zswft0_i',ZSWFT0_ii,klev+1)
     1322        CALL writefield_phy('zfsdwn_i',ZFSDWN_i,klev+1)
     1323        CALL writefield_phy('zflux2_i',ZFLUX_i(:,2,:),klev+1)
     1324        CALL writefield_phy('zfcdwn_i',ZFCDWN_i,klev+1)
     1325        CALL writefield_phy('zfluc2_i',ZFLUC_i(:,2,:),klev+1)
     1326        CALL writefield_phy('psfswdir',PSFSWDIR,6)
     1327        CALL writefield_phy('psfswdif',PSFSWDIF,6)
     1328        CALL writefield_phy('zflux1_i',ZFLUX_i(:,1,:),klev+1)
     1329        CALL writefield_phy('zfluc1_i',ZFLUC_i(:,1,:),klev+1)
     1330        CALL writefield_phy('zfsup_i',ZFSUP_i,klev+1)
     1331        CALL writefield_phy('zfcup_i',ZFCUP_i,klev+1)
     1332      ENDIF
     1333! ---------
     1334! On retablit l'ordre des niveaux lmd pour les tableaux de sortie
     1335! D autre part, on multiplie les resultats SW par fract pour etre coherent
     1336! avec l ancien rayonnement AR4. Si nuit, fract=0 donc pas de
     1337! rayonnement SW. (MPL 260609)
     1338      print*,'On retablit l ordre des niveaux verticaux pour LMDZ'
     1339      print*,'On multiplie les flux SW par fract et LW dwn par -1'
     1340      DO k=0,klev
     1341         DO i=1,klon
     1342         ZEMTD(i,k+1)  = ZEMTD_i(i,klev+1-k)
     1343         ZEMTU(i,k+1)  = ZEMTU_i(i,klev+1-k)
     1344         ZTRSO(i,k+1)  = ZTRSO_i(i,klev+1-k)
     1345!         ZTH(i,k+1)    = ZTH_i(i,klev+1-k)
     1346! AI ATTENTION
     1347          ZLWFT(i,k+1)  = ZLWFT_i(i,klev+1-k)
     1348          ZSWFT(i,k+1)  = ZSWFT_i(i,klev+1-k)*fract(i)
     1349          ZSWFT0_i(i,k+1) = ZSWFT0_ii(i,klev+1-k)*fract(i)
     1350          ZLWFT0_i(i,k+1) = ZLWFT0_ii(i,klev+1-k)
     1351!
     1352         ZFLUP(i,k+1)  = ZFLUX_i(i,1,klev+1-k)
     1353         ZFLDN(i,k+1)  = -1.*ZFLUX_i(i,2,klev+1-k)
     1354         ZFLUP0(i,k+1) = ZFLUC_i(i,1,klev+1-k)
     1355         ZFLDN0(i,k+1) = -1.*ZFLUC_i(i,2,klev+1-k)
     1356         ZFSDN(i,k+1)  = ZFSDWN_i(i,klev+1-k)*fract(i)
     1357         ZFSDN0(i,k+1) = ZFCDWN_i(i,klev+1-k)*fract(i)
     1358         ZFSDNC0(i,k+1)= ZFCCDWN_i(i,klev+1-k)*fract(i)
     1359         ZFSUP (i,k+1) = ZFSUP_i(i,klev+1-k)*fract(i)
     1360         ZFSUP0(i,k+1) = ZFCUP_i(i,klev+1-k)*fract(i)
     1361         ZFSUPC0(i,k+1)= ZFCCUP_i(i,klev+1-k)*fract(i)
     1362         ZFLDNC0(i,k+1)= -1.*ZFLCCDWN_i(i,klev+1-k)
     1363         ZFLUPC0(i,k+1)= ZFLCCUP_i(i,klev+1-k)
     1364         IF (ok_volcan) THEN
     1365            ZSWADAERO(i,k+1)=ZSWADAERO(i,klev+1-k)*fract(i) !--NL
     1366         ENDIF
     1367         
     1368!   Nouveau calcul car visiblement ZSWFT et ZSWFC sont nuls dans RRTM cy32
     1369!   en sortie de radlsw.F90 - MPL 7.01.09
     1370! AI ATTENTION
     1371!         ZSWFT(i,k+1)  = (ZFSDWN_i(i,k+1)-ZFSUP_i(i,k+1))*fract(i)
     1372!         ZSWFT0_i(i,k+1) = (ZFCDWN_i(i,k+1)-ZFCUP_i(i,k+1))*fract(i)
     1373!         ZLWFT(i,k+1) =-ZFLUX_i(i,2,k+1)-ZFLUX_i(i,1,k+1)
     1374!         ZLWFT0_i(i,k+1)=-ZFLUC_i(i,2,k+1)-ZFLUC_i(i,1,k+1)
     1375         ENDDO
     1376      ENDDO
     1377
     1378!--ajout OB
     1379      ZTOPSWADAERO(:) =ZTOPSWADAERO(:) *fract(:)
     1380      ZSOLSWADAERO(:) =ZSOLSWADAERO(:) *fract(:)
     1381      ZTOPSWAD0AERO(:)=ZTOPSWAD0AERO(:)*fract(:)
     1382      ZSOLSWAD0AERO(:)=ZSOLSWAD0AERO(:)*fract(:)
     1383      ZTOPSWAIAERO(:) =ZTOPSWAIAERO(:) *fract(:)
     1384      ZSOLSWAIAERO(:) =ZSOLSWAIAERO(:) *fract(:)
     1385      ZTOPSWCF_AERO(:,1)=ZTOPSWCF_AERO(:,1)*fract(:)
     1386      ZTOPSWCF_AERO(:,2)=ZTOPSWCF_AERO(:,2)*fract(:)
     1387      ZTOPSWCF_AERO(:,3)=ZTOPSWCF_AERO(:,3)*fract(:)
     1388      ZSOLSWCF_AERO(:,1)=ZSOLSWCF_AERO(:,1)*fract(:)
     1389      ZSOLSWCF_AERO(:,2)=ZSOLSWCF_AERO(:,2)*fract(:)
     1390      ZSOLSWCF_AERO(:,3)=ZSOLSWCF_AERO(:,3)*fract(:)
     1391
     1392! ---------
     1393! On renseigne les champs LMDz, pour avoir la meme chose qu'en sortie de
     1394! LW_LMDAR4 et SW_LMDAR4
     1395
     1396      !--fraction of diffuse radiation in surface SW downward radiation
     1397      DO i = 1, kdlon
     1398       IF (fract(i).GT.0.0) THEN
     1399         zdir=SUM(PSFSWDIR(i,:))
     1400         zdif=SUM(PSFSWDIF(i,:))
     1401         zsolswfdiff(i) = zdif/(zdir+zdif)
     1402       ELSE  !--night
     1403         zsolswfdiff(i) = 1.0
     1404       ENDIF
     1405      ENDDO
     1406!
     1407      DO i = 1, kdlon
     1408         zsolsw(i)    = ZSWFT(i,1)
     1409         zsolsw0(i)   = ZSWFT0_i(i,1)
     1410         ztopsw(i)    = ZSWFT(i,klev+1)
     1411         ztopsw0(i)   = ZSWFT0_i(i,klev+1)
     1412         zsollw(i)    = ZLWFT(i,1)
     1413         zsollw0(i)   = ZLWFT0_i(i,1)
     1414         ztoplw(i)    = ZLWFT(i,klev+1)*(-1)
     1415         ztoplw0(i)   = ZLWFT0_i(i,klev+1)*(-1)
     1416!         
     1417         zsollwdown(i)= -1.*ZFLDN(i,1)
     1418      ENDDO
     1419
     1420      DO k=1,kflev
     1421         DO i=1,kdlon
     1422           zheat(i,k)=(ZSWFT(i,k+1)-ZSWFT(i,k))*RDAY*RG/RCPD/PDP(i,k)
     1423           zheat0(i,k)=(ZSWFT0_i(i,k+1)-ZSWFT0_i(i,k))*RDAY*RG/RCPD/PDP(i,k)
     1424           zcool(i,k)=(ZLWFT(i,k)-ZLWFT(i,k+1))*RDAY*RG/RCPD/PDP(i,k)
     1425           zcool0(i,k)=(ZLWFT0_i(i,k)-ZLWFT0_i(i,k+1))*RDAY*RG/RCPD/PDP(i,k)
     1426           IF (ok_volcan) THEN
     1427              zheat_volc(i,k)=(ZSWADAERO(i,k+1)-ZSWADAERO(i,k))*RG/RCPD/PDP(i,k) !NL
     1428              zcool_volc(i,k)=(ZLWADAERO(i,k)-ZLWADAERO(i,k+1))*RG/RCPD/PDP(i,k) !NL
     1429           ENDIF
     1430         ENDDO
     1431      ENDDO
     1432#endif 
     1433  print*,'Fin traitement ECRAD'
     1434! Fin ECRAD
     1435  ENDIF        ! iflag_rrtm
     1436! ecrad
    10641437!======================================================================
    10651438
    … …  
    11021475          solswad_aero(iof+i) = zsolswadaero(i)
    11031476          solswad0_aero(iof+i) = zsolswad0aero(i)
    1104 ! MS the following lines seem to be wrong, why is iof on right hand side???
    1105 !          topsw_aero(iof+i,:) = ztopsw_aero(iof+i,:)
    1106 !          topsw0_aero(iof+i,:) = ztopsw0_aero(iof+i,:)
    1107 !          solsw_aero(iof+i,:) = zsolsw_aero(iof+i,:)
    1108 !          solsw0_aero(iof+i,:) = zsolsw0_aero(iof+i,:)
    11091477          topsw_aero(iof+i,:) = ztopsw_aero(i,:)
    11101478          topsw0_aero(iof+i,:) = ztopsw0_aero(i,:)
    … …  
    11711539 ENDDO ! j = 1, nb_gr
    11721540
     1541IF (lldebug) THEN
     1542 if (0.eq.1) then
     1543! Verifs dans le cas 1D
     1544 print*,'================== Sortie de radlw ================='
     1545 print*,'******** LW LW LW *******************'
     1546 print*,'ZLWFT =',ZLWFT
     1547 print*,'ZLWFT0_i =',ZLWFT0_i
     1548 print*,'ZFLUP0 =',ZFLUP0
     1549 print*,'ZFLDN0 =',ZFLDN0
     1550 print*,'ZFLDNC0 =',ZFLDNC0
     1551 print*,'ZFLUPC0 =',ZFLUPC0
     1552
     1553 print*,'******** SW SW SW *******************'
     1554 print*,'ZSWFT =',ZSWFT
     1555 print*,'ZSWFT0_i =',ZSWFT0_i
     1556 print*,'ZFSDN =',ZFSDN
     1557 print*,'ZFSDN0 =',ZFSDN0
     1558 print*,'ZFSDNC0 =',ZFSDNC0
     1559 print*,'ZFSUP =',ZFSUP
     1560 print*,'ZFSUP0 =',ZFSUP0
     1561 print*,'ZFSUPC0 =',ZFSUPC0
     1562
     1563 print*,'******** LMDZ  *******************'
     1564 print*,'cool = ', cool
     1565 print*,'heat = ', heat
     1566 print*,'topsw = ', topsw
     1567 print*,'toplw = ', toplw
     1568 print*,'sollw = ', sollw
     1569 print*,'solsw = ', solsw
     1570 print*,'lwdn = ', lwdn
     1571 print*,'lwup = ', lwup
     1572 print*,'swdn = ', swdn
     1573 print*,'swup =', swup
     1574 endif
     1575ENDIF
     1576
    11731577END SUBROUTINE radlwsw
    11741578
Note: See TracChangeset for help on using the changeset viewer.