From 452b361e9e4d6b4b059c4e43d7953abbf01a551b Mon Sep 17 00:00:00 2001
From: Quentin Rodier <quentin.rodier@meteo.fr>
Date: Fri, 8 Dec 2023 16:21:17 +0100
Subject: [PATCH] Valery M. 08/12/2023: add new option NRAD_AGG to compute
radiation only on an average column of size (NRAD_AGG^2 points). It reduces
by NRAD_AGG^2 the cost of the radiation calls but is not yet fully optimized
when the number of inner domain by cores are lower than 20x20 columns
---
src/MNH/default_desfmn.f90 | 1 +
src/MNH/ini_modeln.f90 | 5 +-
src/MNH/ini_radiations.f90 | 28 +-
src/MNH/ini_radiations_agg.f90 | 202 ++++++++++
src/MNH/modd_param_radn.f90 | 3 +
src/MNH/modd_radiationsn.f90 | 16 +
src/MNH/modn_param_radn.f90 | 8 +-
src/MNH/phys_paramn.f90 | 4 +-
src/MNH/radiations.f90 | 63 ++-
src/MNH/radiations_agg.f90 | 680 +++++++++++++++++++++++++++++++++
src/MNH/read_exsegn.f90 | 19 +
11 files changed, 983 insertions(+), 46 deletions(-)
create mode 100644 src/MNH/ini_radiations_agg.f90
create mode 100644 src/MNH/radiations_agg.f90
diff --git a/src/MNH/default_desfmn.f90 b/src/MNH/default_desfmn.f90
index 07359fc5a..9ccf63878 100644
--- a/src/MNH/default_desfmn.f90
+++ b/src/MNH/default_desfmn.f90
@@ -727,6 +727,7 @@ COPILW = 'EBCU'
XFUDG = 1.
LAERO_FT=.FALSE.
LFIX_DAT=.FALSE.
+NRAD_AGG = 1
!
#ifdef MNH_ECRAD
!* 13bis. SET DEFAULT VALUES FOR MODD_PARAM_ECRAD_n :
diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90
index 8f00cebc5..ce91fc05c 100644
--- a/src/MNH/ini_modeln.f90
+++ b/src/MNH/ini_modeln.f90
@@ -389,6 +389,7 @@ USE MODD_PASPOL
USE MODD_PASPOL_n
USE MODD_PAST_FIELD_n
use modd_precision, only: LFIINT
+USE MODD_PARAM_RAD_n
USE MODD_RADIATIONS_n
USE MODD_RECYCL_PARAM_n
USE MODD_REF
@@ -2519,7 +2520,9 @@ IF (CRAD /= 'NONE') THEN
TDTRAD_FULL,TDTRAD_CLONLY, &
TZINITHALO2D_ll, &
XRADEFF,XSWU,XSWD,XLWU, &
- XLWD,XDTHRADSW,XDTHRADLW )
+ XLWD,XDTHRADSW,XDTHRADLW, &
+ NRAD_AGG,NI_RAD_AGG,NJ_RAD_AGG, &
+ NIOR_RAD_AGG,NJOR_RAD_AGG )
!
IF (GINIRAD) CALL SUNPOS_n(XZENITH,PAZIMSOL=XAZIM)
CALL SURF_SOLAR_GEOM (XZS, XZS_XY)
diff --git a/src/MNH/ini_radiations.f90 b/src/MNH/ini_radiations.f90
index 883179c8b..2b11a5890 100644
--- a/src/MNH/ini_radiations.f90
+++ b/src/MNH/ini_radiations.f90
@@ -15,7 +15,8 @@ INTERFACE
PDTHRAD,PDIRFLASWD,PSCAFLASWD, &
PFLALWD,PDIRSRFSWD,KCLEARCOL_TM1, &
PZENITH, PAZIM, TPDTRAD_FULL,TPDTRAD_CLONLY,TPINITHALO2D_ll, &
- PRADEFF,PSWU,PSWD,PLWU,PLWD,PDTHRADSW,PDTHRADLW )
+ PRADEFF,PSWU,PSWD,PLWU,PLWD,PDTHRADSW,PDTHRADLW, &
+ KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG )
!
USE MODD_ARGSLIST_ll, ONLY : LIST_ll
USE MODD_IO, ONLY : TFILEDATA
@@ -59,6 +60,13 @@ REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDTHRADSW ! dthrad sw
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDTHRADLW ! dthrad lw
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PRADEFF ! effective radius
!
+!
+INTEGER, INTENT(IN) :: KRAD_AGG ! number of aggregated points
+INTEGER, INTENT(OUT) :: KI_RAD_AGG ! reformatted X array size
+INTEGER, INTENT(OUT) :: KJ_RAD_AGG ! reformatted Y array size
+INTEGER, INTENT(OUT) :: KIOR_RAD_AGG ! index of first point of packed array according to current domain
+INTEGER, INTENT(OUT) :: KJOR_RAD_AGG ! index of first point of packed array according to current domain
+
END SUBROUTINE INI_RADIATIONS
!
END INTERFACE
@@ -73,7 +81,8 @@ END MODULE MODI_INI_RADIATIONS
PDTHRAD,PDIRFLASWD,PSCAFLASWD, &
PFLALWD,PDIRSRFSWD,KCLEARCOL_TM1, &
PZENITH, PAZIM, TPDTRAD_FULL,TPDTRAD_CLONLY,TPINITHALO2D_ll, &
- PRADEFF,PSWU,PSWD,PLWU,PLWD,PDTHRADSW,PDTHRADLW )
+ PRADEFF,PSWU,PSWD,PLWU,PLWD,PDTHRADSW,PDTHRADLW, &
+ KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG )
! ####################################################################
!
!!**** *INI_RADIATION_TIME * - initialisation for radiation scheme in the MesoNH framework
@@ -130,6 +139,7 @@ USE MODE_IO_FIELD_READ, only: IO_Field_read
USE MODE_ll
!
USE MODI_SHUMAN
+USE MODI_INI_RADIATIONS_AGG
!
IMPLICIT NONE
!
@@ -172,6 +182,12 @@ REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDTHRADSW ! dthrad sw
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDTHRADLW ! dthrad lw
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PRADEFF ! effective radius
!
+INTEGER, INTENT(IN) :: KRAD_AGG ! number of aggregated points
+INTEGER, INTENT(OUT) :: KI_RAD_AGG ! reformatted X array size
+INTEGER, INTENT(OUT) :: KJ_RAD_AGG ! reformatted Y array size
+INTEGER, INTENT(OUT) :: KIOR_RAD_AGG ! index of first point of packed array according to current domain
+INTEGER, INTENT(OUT) :: KJOR_RAD_AGG ! index of first point of packed array according to current domain
+!
!* 0.2 declarations of local variables
!
INTEGER, DIMENSION(0:11) :: IBIS, INOBIS ! Cumulative number of days per month
@@ -321,6 +337,14 @@ ELSE
CALL IO_Field_read(TPINIFILE,'ZENITH', PZENITH)
CALL IO_Field_read(TPINIFILE,'AZIM', PAZIM)
END IF
+!
+!-------------------------------------------------------------------------------
+!
+!* 10. INITIALIZE COLUMN AGGREGATION FOR RADIATION CALL
+! -------------------------------------------------
+
+CALL INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG)
+!
!-------------------------------------------------------------------------------
!
END SUBROUTINE INI_RADIATIONS
diff --git a/src/MNH/ini_radiations_agg.f90 b/src/MNH/ini_radiations_agg.f90
new file mode 100644
index 000000000..0ec198677
--- /dev/null
+++ b/src/MNH/ini_radiations_agg.f90
@@ -0,0 +1,202 @@
+!MNH_LIC Copyright 1995-2022 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+!-----------------------------------------------------------------
+! ########################
+ MODULE MODI_INI_RADIATIONS_AGG
+! ########################
+!
+INTERFACE
+
+ SUBROUTINE INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG)
+INTEGER, INTENT(IN) :: KRAD_AGG ! number of aggregated points
+INTEGER, INTENT(OUT) :: KI_RAD_AGG ! reformatted X array size
+INTEGER, INTENT(OUT) :: KJ_RAD_AGG ! reformatted Y array size
+INTEGER, INTENT(OUT) :: KIOR_RAD_AGG ! index of first point of packed array according to current domain
+INTEGER, INTENT(OUT) :: KJOR_RAD_AGG ! index of first point of packed array according to current domain
+END SUBROUTINE INI_RADIATIONS_AGG
+
+END INTERFACE
+
+END MODULE MODI_INI_RADIATIONS_AGG
+!
+! ############################################################################
+ SUBROUTINE INI_RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG)
+! ############################################################################
+!
+!!**** *INI_RADIATIONS_AGG * - routine to call the SW and LW radiation calculations
+!!
+!! PURPOSE
+!! -------
+!! The purpose of this routine is to aggregate columns of the temperature, water vapor
+!! liquid water, cloud fraction, ozone profiles for the ECMWF radiation
+!! calculations. There is a great number of available radiative fluxes in
+!! the output, but only the potential temperature radiative tendency and the
+!! SW and LW surface fluxes are provided in the output of the routine.
+!!
+!!** METHOD
+!! ------
+!!
+!! All columns are aggregated according to NRAD_AGG * NRAD_AGG points.
+!! if NRAD_AGG = 1 : no aggregation (default)
+!! if NRAD_AGG = 2 ; aggregation on 4 points, bottom/left is 1st (modulo 2) physical point
+!! of the whole domain (not processor domain)
+!! if NRAD_AGG = 3 ; aggregation on 9 points, bottom/left is 1st (modulo 3) physical point
+!! of the whole domain (not processor domain) ; NHALO must be at least equal to 2
+!! if NRAD_AGG = 4 ; aggregation on 16 points, bottom/left is 1st (modulo 4) physical point
+!! of the whole domain (not processor domain) ; NHALO must be at least equal to 3
+!! if NRAD_AGG = 5 ; aggregation on 25 points, bottom/left is 1st (modulo 5) physical point
+!! of the whole domain (not processor domain) ; NHALO must be at least equal to 4
+!!
+!!
+!! EXTERNAL
+!! --------
+!! Subroutine INI_RADIATIONS_AGG
+!!
+!! IMPLICIT ARGUMENTS
+!! ------------------
+!!
+!! REFERENCE
+!! ---------
+!!
+!! AUTHOR
+!! ------
+!! V. Masson * CNRM *
+!!
+!! MODIFICATIONS
+!! -------------
+!! Original 27/10/23
+!-------------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+USE PARKIND1, ONLY: JPRB
+!
+USE MODD_CST
+USE MODD_LUNIT_n, ONLY : TLUOUT
+USE MODD_CONF, ONLY : NHALO
+USE MODD_LBC_n, ONLY : CLBCX, CLBCY
+USE MODD_PARAMETERS, ONLY : JPHEXT
+!
+USE MODE_ll
+USE MODE_MSG
+USE MODE_MODELN_HANDLER
+USE MODD_VAR_ll, ONLY : IP
+!
+!
+IMPLICIT NONE
+!
+!* 0.1 DECLARATIONS OF DUMMY ARGUMENTS :
+!
+INTEGER, INTENT(IN) :: KRAD_AGG ! number of aggregated points
+INTEGER, INTENT(OUT) :: KI_RAD_AGG ! reformatted X array size
+INTEGER, INTENT(OUT) :: KJ_RAD_AGG ! reformatted Y array size
+INTEGER, INTENT(OUT) :: KIOR_RAD_AGG ! index of first point of packed array according to current domain
+INTEGER, INTENT(OUT) :: KJOR_RAD_AGG ! index of first point of packed array according to current domain
+!
+!
+!* 0.2 DECLARATIONS OF LOCAL VARIABLES
+!
+!
+INTEGER :: IIB ! I index value of the first inner mass point (current domain)
+INTEGER :: IJB ! J index value of the first inner mass point
+INTEGER :: IIE ! I index value of the last inner mass point
+INTEGER :: IJE ! J index value of the last inner mass point
+INTEGER :: IIB_ll ! I index value of the first inner mass point (whole domain)
+INTEGER :: IJB_ll ! J index value of the first inner mass point
+INTEGER :: IIE_ll ! I index value of the last inner mass point
+INTEGER :: IJE_ll ! J index value of the last inner mass point
+INTEGER :: IIMAX_ll ! I index size of the whole physical domain
+INTEGER :: IJMAX_ll ! J index size of the whole physical domain
+!
+INTEGER :: IIOR_RAD_AGG_ll ! index of first point of packed array according to whole domain
+INTEGER :: IJOR_RAD_AGG_ll ! index of first point of packed array according to whole domain
+!
+INTEGER :: IIOR_ll ! index of first point in the processor relative to the whole domain
+INTEGER :: IJOR_ll ! index of first point in the processor relative to the whole domain
+!
+INTEGER :: ILUOUT ! Logical unit
+INTEGER :: IMI
+!-------------------------------------------------------------------------
+!-------------------------------------------------------------------------
+!-------------------------------------------------------------------------
+!
+!
+!* 1. CHECK OF COHERENCE
+! ------------------
+!
+IMI = GET_CURRENT_MODEL_INDEX()
+
+ILUOUT = TLUOUT%NLU
+!
+IF (KRAD_AGG > NHALO+1) THEN
+ WRITE(ILUOUT,*) ' +------------------------------------------------------+'
+ WRITE(ILUOUT,*) ' [ Error in Radiation columns aggregation |'
+ WRITE(ILUOUT,*) ' [ NRAD_AGG = ',KRAD_AGG,' |'
+ WRITE(ILUOUT,*) ' [ NHALO = ',NHALO,' |'
+ WRITE(ILUOUT,*) ' [ NRAD_AGG must be smaller than or equal to NHALO+1 |'
+ WRITE(ILUOUT,*) ' +------------------------------------------------------+'
+ !
+ CALL PRINT_MSG( NVERB_FATAL, 'GEN', 'INI_RADIATIONS_AGG','Incoherence between NRAD_AGG and NHALO' )
+END IF
+!
+!* 2. COMPUTE DIMENSIONS OF ARRAYS AND OTHER INDICES
+! ----------------------------------------------
+!
+! full arrays
+CALL GET_INDICE_ll (IIB,IJB,IIE,IJE) ! this definition must be coherent with
+ ! the one used in ini_radiations routine
+!
+! index of first physical point of the whole domain, relative to the whole domain
+!
+! Physical global domain bounds
+!
+CALL GET_GLOBALDIMS_ll ( IIMAX_ll,IJMAX_ll)
+IF (CLBCX(1)=='CYCL') THEN
+ IIB_ll = 1 + NHALO
+ IIE_ll = IIMAX_ll + NHALO
+ELSE
+ IIB_ll = 1 + JPHEXT
+ IIE_ll = IIMAX_ll + JPHEXT
+END IF
+
+IF (CLBCX(2)=='CYCL') THEN
+ IJB_ll = 1 + NHALO
+ IJE_ll = IJMAX_ll + NHALO
+ELSE
+ IJB_ll = 1 + JPHEXT
+ IJE_ll = IJMAX_ll + JPHEXT
+END IF
+!
+! index of first point in the modelessor relative to the whole domain (includes non-physical moints)
+!
+CALL GET_OR_ll('B',IIOR_ll,IJOR_ll)
+!
+! index of the first PACKED radiative column domain in the current processor relative to the whole domain
+! There must be an entire number of packed domain between first global physical point and the packed domain in this processor
+IIOR_RAD_AGG_ll = KRAD_AGG * ((IIOR_ll-1) / KRAD_AGG) + IIB_ll
+IJOR_RAD_AGG_ll = KRAD_AGG * ((IJOR_ll-1) / KRAD_AGG) + IJB_ll
+!
+!
+! index of the first PACKED radiative column domain in the current processor relative to the current processor domain
+! It can be lower than first physical point index (i.e. inside the HALO points)
+KIOR_RAD_AGG = IIOR_RAD_AGG_ll - IIOR_ll + 1
+KJOR_RAD_AGG = IJOR_RAD_AGG_ll - IJOR_ll + 1
+!
+!
+! removes packed columns that are entirely included in the HALO (no need to compute them on this processor)
+!
+KIOR_RAD_AGG = KIOR_RAD_AGG + ((IIB-KIOR_RAD_AGG)/KRAD_AGG) * KRAD_AGG
+KJOR_RAD_AGG = KJOR_RAD_AGG + ((IJB-KJOR_RAD_AGG)/KRAD_AGG) * KRAD_AGG
+!
+!
+! Number of PACKED radiative subdomains inside current processor domain
+KI_RAD_AGG = (IIE - KIOR_RAD_AGG) / KRAD_AGG + 1
+KJ_RAD_AGG = (IJE - KJOR_RAD_AGG) / KRAD_AGG + 1
+!
+!-------------------------------------------------------------------------------
+!
+END SUBROUTINE INI_RADIATIONS_AGG
+
diff --git a/src/MNH/modd_param_radn.f90 b/src/MNH/modd_param_radn.f90
index 51f4260ec..f429c6938 100644
--- a/src/MNH/modd_param_radn.f90
+++ b/src/MNH/modd_param_radn.f90
@@ -84,6 +84,7 @@ TYPE PARAM_RAD_t
! aerosol and ozone distribution
LOGICAL :: LFIX_DAT ! logical switch to fix the date to a constant
! perpetual day ( diurnal cycle is conserved)
+ INTEGER:: NRAD_AGG ! number of aggregated points in each direction in radiative aggregated columns
!-------------------------------------------------------------------------------
!
@@ -115,6 +116,7 @@ CHARACTER (LEN=4), POINTER :: COPWSW=>NULL()
CHARACTER (LEN=4), POINTER :: COPISW=>NULL()
LOGICAL, POINTER :: LAERO_FT=>NULL()
LOGICAL, POINTER :: LFIX_DAT=>NULL()
+INTEGER, POINTER :: NRAD_AGG=>NULL()
CONTAINS
@@ -142,6 +144,7 @@ COPWSW=>PARAM_RAD_MODEL(KTO)%COPWSW
COPISW=>PARAM_RAD_MODEL(KTO)%COPISW
LAERO_FT=>PARAM_RAD_MODEL(KTO)%LAERO_FT
LFIX_DAT=>PARAM_RAD_MODEL(KTO)%LFIX_DAT
+NRAD_AGG=>PARAM_RAD_MODEL(KTO)%NRAD_AGG
END SUBROUTINE PARAM_RAD_GOTO_MODEL
diff --git a/src/MNH/modd_radiationsn.f90 b/src/MNH/modd_radiationsn.f90
index 743c37074..64be9ce63 100644
--- a/src/MNH/modd_radiationsn.f90
+++ b/src/MNH/modd_radiationsn.f90
@@ -127,6 +127,12 @@ TYPE RADIATIONS_t
REAL, DIMENSION(:,:,:), POINTER :: XDTHRADLW => NULL() ! DTHRAD LW
REAL, DIMENSION(:,:,:), POINTER :: XRADEFF => NULL() ! effective radius
!
+! For aggregation of columns into packed radiative columns
+ INTEGER :: NI_RAD_AGG ! reformatted X array size
+ INTEGER :: NJ_RAD_AGG ! reformatted Y array size
+ INTEGER :: NIOR_RAD_AGG ! index of first point of packed array according to current domain
+ INTEGER :: NJOR_RAD_AGG ! index of first point of packed array according to current domain
+!
END TYPE RADIATIONS_t
TYPE(RADIATIONS_t), DIMENSION(JPMODELMAX), TARGET, SAVE :: RADIATIONS_MODEL
@@ -186,6 +192,11 @@ REAL, DIMENSION(:,:,:), POINTER :: XDTHRADSW=>NULL()
REAL, DIMENSION(:,:,:), POINTER :: XDTHRADLW=>NULL()
REAL, DIMENSION(:,:,:), POINTER :: XRADEFF=>NULL()
+INTEGER, POINTER :: NI_RAD_AGG=>NULL() ! reformatted X array size
+INTEGER, POINTER :: NJ_RAD_AGG=>NULL() ! reformatted Y array size
+INTEGER, POINTER :: NIOR_RAD_AGG=>NULL() ! index of first point of packed array according to current domain
+INTEGER, POINTER :: NJOR_RAD_AGG=>NULL() ! index of first point of packed array according to current domain
+
CONTAINS
SUBROUTINE RADIATIONS_GOTO_MODEL(KFROM, KTO)
@@ -282,6 +293,11 @@ XDTHRADSW=>RADIATIONS_MODEL(KTO)%XDTHRADSW
XDTHRADLW=>RADIATIONS_MODEL(KTO)%XDTHRADLW
XRADEFF=>RADIATIONS_MODEL(KTO)%XRADEFF
+NI_RAD_AGG=>RADIATIONS_MODEL(KTO)%NI_RAD_AGG
+NJ_RAD_AGG=>RADIATIONS_MODEL(KTO)%NJ_RAD_AGG
+NIOR_RAD_AGG=>RADIATIONS_MODEL(KTO)%NIOR_RAD_AGG
+NJOR_RAD_AGG=>RADIATIONS_MODEL(KTO)%NJOR_RAD_AGG
+
END SUBROUTINE RADIATIONS_GOTO_MODEL
END MODULE MODD_RADIATIONS_n
diff --git a/src/MNH/modn_param_radn.f90 b/src/MNH/modn_param_radn.f90
index 87f962d75..477e14df9 100644
--- a/src/MNH/modn_param_radn.f90
+++ b/src/MNH/modn_param_radn.f90
@@ -34,7 +34,8 @@ USE MODD_PARAM_RAD_n, ONLY: &
COPISW_n => COPISW, &
XFUDG_n => XFUDG, &
LAERO_FT_n => LAERO_FT, &
- LFIX_DAT_n => LFIX_DAT
+ LFIX_DAT_n => LFIX_DAT, &
+ NRAD_AGG_n => NRAD_AGG
!
IMPLICIT NONE
@@ -56,10 +57,11 @@ CHARACTER (LEN=4), SAVE :: COPISW
REAL, SAVE :: XFUDG
LOGICAL,SAVE :: LAERO_FT
LOGICAL,SAVE :: LFIX_DAT
+INTEGER,SAVE :: NRAD_AGG
!
NAMELIST/NAM_PARAM_RADn/XDTRAD,XDTRAD_CLONLY,LCLEAR_SKY,NRAD_COLNBR,&
NRAD_DIAG,CLW,CAER,CAOP,CEFRADL,CEFRADI,COPWLW,&
- COPILW,COPWSW,COPISW,XFUDG,LAERO_FT,LFIX_DAT
+ COPILW,COPWSW,COPISW,XFUDG,LAERO_FT,LFIX_DAT, NRAD_AGG
!
CONTAINS
!
@@ -81,6 +83,7 @@ SUBROUTINE INIT_NAM_PARAM_RADn
XFUDG = XFUDG_n
LAERO_FT = LAERO_FT_n
LFIX_DAT = LFIX_DAT_n
+ NRAD_AGG = NRAD_AGG_n
END SUBROUTINE INIT_NAM_PARAM_RADn
SUBROUTINE UPDATE_NAM_PARAM_RADn
@@ -101,6 +104,7 @@ SUBROUTINE UPDATE_NAM_PARAM_RADn
XFUDG_n = XFUDG
LAERO_FT_n = LAERO_FT
LFIX_DAT_n = LFIX_DAT
+ NRAD_AGG_n = NRAD_AGG
END SUBROUTINE UPDATE_NAM_PARAM_RADn
END MODULE MODN_PARAM_RAD_n
diff --git a/src/MNH/phys_paramn.f90 b/src/MNH/phys_paramn.f90
index cc55fe411..8c14197c8 100644
--- a/src/MNH/phys_paramn.f90
+++ b/src/MNH/phys_paramn.f90
@@ -354,7 +354,7 @@ USE MODI_GROUND_PARAM_n
USE MODI_GRADIENT_M
USE MODI_GRADIENT_W
USE MODI_PASPOL
-USE MODI_RADIATIONS
+USE MODI_RADIATIONS_AGG
USE MODI_SALT_FILTER
USE MODI_SEDIM_DUST
USE MODI_SEDIM_SALT
@@ -777,7 +777,7 @@ CALL SUNPOS_n ( XZENITH, ZCOSZEN, ZSINZEN, ZAZIMSOL )
XLWD(:,:,:)=0.0
XDTHRADSW(:,:,:)=0.0
XDTHRADLW(:,:,:)=0.0
- CALL RADIATIONS( TPFILE, &
+ CALL RADIATIONS_AGG(NRAD_AGG,NI_RAD_AGG,NJ_RAD_AGG,NIOR_RAD_AGG,NJOR_RAD_AGG, TPFILE, &
LCLEAR_SKY, OCLOUD_ONLY, NCLEARCOL_TM1, CEFRADL, CEFRADI, COPWSW, COPISW, &
COPWLW, COPILW, XFUDG, &
NDLON, NFLEV, NRAD_DIAG, NFLUX, NRAD, NAER, NSWB_OLD, NSWB_MNH, NLWB_MNH, &
diff --git a/src/MNH/radiations.f90 b/src/MNH/radiations.f90
index f4db08bfc..7917ae111 100644
--- a/src/MNH/radiations.f90
+++ b/src/MNH/radiations.f90
@@ -10,10 +10,11 @@
CONTAINS
!
! ############################################################################
- SUBROUTINE RADIATIONS (TPFILE,OCLEAR_SKY,OCLOUD_ONLY, &
+ SUBROUTINE RADIATIONS (TPFILE,KIB,KIE,KJB,KJE,OCLEAR_SKY,OCLOUD_ONLY, &
KCLEARCOL_TM1,HEFRADL,HEFRADI,HOPWSW,HOPISW,HOPWLW,HOPILW, &
PFUDG, KDLON, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD, &
KSWB_MNH,KLWB_MNH, KSTATM,KRAD_COLNBR,PCOSZEN,PSEA, PCORSOL, &
+ PLAT, PLON, &
PDIR_ALB, PSCA_ALB,PEMIS, PCLDFR, PCCO2, PTSRAD, PSTATM, &
PTHT, PRT, PPABST, POZON, PAER, PDST_WL, PAER_CLIM, PSVT, &
PDTHRAD, PSRFLWD, PSRFSWD_DIR,PSRFSWD_DIF, PRHODREF, PZZ, &
@@ -178,6 +179,10 @@ IMPLICIT NONE
!* 0.1 DECLARATIONS OF DUMMY ARGUMENTS :
!
TYPE(TFILEDATA), INTENT(IN) :: TPFILE ! Output file
+INTEGER, INTENT(IN) :: KIB ! first X physical point in array
+INTEGER, INTENT(IN) :: KIE ! last X physical point in array
+INTEGER, INTENT(IN) :: KJB ! first Y physical point in array
+INTEGER, INTENT(IN) :: KJE ! last Y physical point in array
LOGICAL, INTENT(IN) :: OCLOUD_ONLY! flag for the cloud column
! computations only
LOGICAL, INTENT(IN) :: OCLEAR_SKY !
@@ -214,6 +219,8 @@ CHARACTER (LEN=*), INTENT (IN) :: HOPILW !ice water LW optical properties
REAL, INTENT(IN) :: PFUDG ! subgrid cloud inhomogenity factor
REAL, DIMENSION(:,:), INTENT(IN) :: PCOSZEN ! COS(zenithal solar angle)
REAL, INTENT(IN) :: PCORSOL ! SOLar constant CORrection
+REAL, DIMENSION(:,:), INTENT(IN) :: PLAT ! Latitude
+REAL, DIMENSION(:,:), INTENT(IN) :: PLON ! Longitude
REAL, DIMENSION(:,:), INTENT(IN) :: PSEA ! Land-sea mask
!
REAL, DIMENSION(:,:,:), INTENT(IN) :: PDIR_ALB! Surface direct ALBedo
@@ -229,10 +236,10 @@ REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT ! moist variables at t (humidity
REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! pressure at t
REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT ! scalar variable ( C2R2 and C1R3 particle)
!
-REAL, DIMENSION(:,:,:), POINTER :: POZON ! OZONE field from clim.
-REAL, DIMENSION(:,:,:,:), POINTER :: PAER ! AERosols optical thickness from clim.
-REAL, DIMENSION(:,:,:,:), POINTER :: PDST_WL ! AERosols Extinction by wavelength .
-REAL, DIMENSION(:,:,:,:), POINTER :: PAER_CLIM ! AERosols optical thickness from clim.
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: POZON ! OZONE field from clim.
+REAL, DIMENSION(:,:,:,:), INTENT(INOUT) :: PAER ! AERosols optical thickness from clim.
+REAL, DIMENSION(:,:,:,:), INTENT(INOUT) :: PDST_WL ! AERosols Extinction by wavelength .
+REAL, DIMENSION(:,:,:,:), INTENT(INOUT) :: PAER_CLIM ! AERosols optical thickness from clim.
! note : the vertical dimension of
! these fields include the "radiation levels"
! above domain top
@@ -558,8 +565,11 @@ REAL :: ZCLEAR_COL_ll , ZDLON_ll
!* 1. COMPUTE DIMENSIONS OF ARRAYS AND OTHER INDICES
! ----------------------------------------------
!
-CALL GET_INDICE_ll (IIB,IJB,IIE,IJE) ! this definition must be coherent with
- ! the one used in ini_radiations routine
+IIB = KIB
+IIE = KIE
+IJB = KJB
+IJE = KJE
+
IKU = SIZE(PTHT,3)
IKB = 1 + JPVEXT
IKE = IKU - JPVEXT
@@ -569,38 +579,17 @@ IKUP = IKE-JPVEXT+1
!
ISWB = SIZE(PSRFSWD_DIR,3)
!
-!-------------------------------------------------------------------------------
-!* 1.1 CHECK PRESSURE DECREASING
-! -------------------------
-ZDZPABST(:,:,1:IKU-1) = PPABST(:,:,1:IKU-1) - PPABST(:,:,2:IKU)
-ZDZPABST(:,:,IKU) = ZDZPABST(:,:,IKU-1)
-!
-ZMINVAL=MIN_ll(ZDZPABST,IINFO_ll)
-!
-IF ( ZMINVAL <= 0.0 ) THEN
- ILUOUT = TLUOUT%NLU
- IMINLOC=GMINLOC_ll( ZDZPABST )
- WRITE(ILUOUT,*) ' radiation.f90 STOP :: SOMETHING WRONG WITH PRESSURE , ZDZPABST <= 0.0 '
- WRITE(ILUOUT,*) ' radiation :: ZDZPABST ', ZMINVAL,' located at ', IMINLOC
- FLUSH(unit=ILUOUT)
- call Print_msg( NVERB_FATAL, 'GEN', 'RADIATIONS', 'something wrong with pressure: ZDZPABST <= 0.0' )
-
-ENDIF
!------------------------------------------------------------------------------
ALLOCATE(ZLAT(KDLON))
ALLOCATE(ZLON(KDLON))
-IF(LCARTESIAN) THEN
- ZLAT(:) = XLAT0*(XPI/180.)
- ZLON(:) = XLON0*(XPI/180.)
-ELSE
- DO JJ=IJB,IJE
- DO JI=IIB,IIE
- IIJ = 1 + (JI-IIB) + (IIE-IIB+1)*(JJ-IJB)
- ZLAT(IIJ) = XLAT(JI,JJ)*(XPI/180.)
- ZLON(IIJ) = XLON(JI,JJ)*(XPI/180.)
- END DO
+DO JJ=IJB,IJE
+ DO JI=IIB,IIE
+ IIJ = 1 + (JI-IIB) + (IIE-IIB+1)*(JJ-IJB)
+ ZLAT(IIJ) = PLAT(JI,JJ)
+ ZLON(IIJ) = PLON(JI,JJ)
END DO
-END IF
+END DO
+!
!-------------------------------------------------------------------------------
!
!* 2. INITIALIZES THE MEAN-LAYER VARIABLES
@@ -1215,10 +1204,6 @@ IF(OCLOUD_ONLY .OR. OCLEAR_SKY) THEN
CALL REDUCESUM_ll(ZCLEAR_COL_ll,IINFO_ll)
!ZDLON_ll = KDLON
!CALL REDUCESUM_ll(ZDLON_ll,IINFO_ll)
-
- !IF (IP == 1 )
- !print*,",RADIATIOn COULD_ONLY=OCLOUD_ONLY,OCLEAR_SKY,ZCLEAR_COL_ll,ICLEAR_COL,ICLOUD_COL,KDON,ZDLON_ll,GNOCL=", &
- ! OCLOUD_ONLY,OCLEAR_SKY,ZCLEAR_COL_ll,ICLEAR_COL,ICLOUD_COL,KDLON,ZDLON_ll,GNOCL
!
!!$ IF( ICLEAR_COL /=0 ) THEN ! at least one clear-sky column exists -> average profiles on clear columns
IF( ZCLEAR_COL_ll /= 0.0 ) THEN ! at least one clear-sky column exists -> average profiles on clear columns
diff --git a/src/MNH/radiations_agg.f90 b/src/MNH/radiations_agg.f90
new file mode 100644
index 000000000..bbc5d383f
--- /dev/null
+++ b/src/MNH/radiations_agg.f90
@@ -0,0 +1,680 @@
+!MNH_LIC Copyright 1995-2022 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
+!MNH_LIC for details. version 1.
+!-----------------------------------------------------------------
+! ########################
+ MODULE MODI_RADIATIONS_AGG
+! ########################
+!
+INTERFACE
+
+ SUBROUTINE RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG, &
+ TPFILE,OCLEAR_SKY,OCLOUD_ONLY, &
+ KCLEARCOL_TM1,HEFRADL,HEFRADI,HOPWSW,HOPISW,HOPWLW,HOPILW, &
+ PFUDG, KDLON, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD, &
+ KSWB_MNH,KLWB_MNH, KSTATM,KRAD_COLNBR,PCOSZEN,PSEA, PCORSOL, &
+ PDIR_ALB, PSCA_ALB,PEMIS, PCLDFR, PCCO2, PTSRAD, PSTATM, &
+ PTHT, PRT, PPABST, POZON, PAER, PDST_WL, PAER_CLIM, PSVT, &
+ PDTHRAD, PSRFLWD, PSRFSWD_DIR,PSRFSWD_DIF, PRHODREF, PZZ, &
+ PRADEFF, PSWU, PSWD, PLWU,PLWD, PDTHRADSW, PDTHRADLW )
+!
+USE MODD_IO, ONLY: TFILEDATA
+
+INTEGER, INTENT(IN) :: KRAD_AGG ! number of aggregated points
+INTEGER, INTENT(IN) :: KI_RAD_AGG ! reformatted X array size
+INTEGER, INTENT(IN) :: KJ_RAD_AGG ! reformatted Y array size
+INTEGER, INTENT(IN) :: KIOR_RAD_AGG ! index of first point of packed array according to current domain
+INTEGER, INTENT(IN) :: KJOR_RAD_AGG ! index of first point of packed array according to current domain
+
+TYPE(TFILEDATA), INTENT(IN) :: TPFILE ! Output file
+LOGICAL, INTENT(IN) :: OCLOUD_ONLY! flag for the cloud column
+ ! computations only
+LOGICAL, INTENT(IN) :: OCLEAR_SKY !
+INTEGER, INTENT(IN) :: KDLON ! number of columns where the
+ ! radiation calculations are
+ ! performed
+INTEGER, INTENT(IN) :: KFLEV ! number of vertical levels
+ ! where the radiation
+ ! calculations are performed
+INTEGER, INTENT(IN) :: KRAD_DIAG ! index for the number of
+ ! fields in the output
+INTEGER, INTENT(IN) :: KFLUX ! number of top and ground
+ ! fluxes for the ZFLUX array
+INTEGER, INTENT(IN) :: KRAD ! number of satellite radiances
+ ! for the ZRAD and ZRADCS arrays
+INTEGER, INTENT(IN) :: KAER ! number of AERosol classes
+
+INTEGER, INTENT(IN) :: KSWB_OLD ! number of SW band ECMWF
+INTEGER, INTENT(IN) :: KSWB_MNH ! number of SW band ECRAD
+INTEGER, INTENT(IN) :: KLWB_MNH ! number of LW band ECRAD
+INTEGER, INTENT(IN) :: KSTATM ! index of the standard
+ ! atmosphere level just above
+ ! the model top
+INTEGER, INTENT(IN) :: KRAD_COLNBR ! factor by which the memory
+ ! is split
+ !
+ !Choice of :
+CHARACTER (LEN=*), INTENT (IN) :: HEFRADL !
+CHARACTER (LEN=*), INTENT (IN) :: HEFRADI !
+CHARACTER (LEN=*), INTENT (IN) :: HOPWSW !cloud water SW optical properties
+CHARACTER (LEN=*), INTENT (IN) :: HOPISW !ice water SW optical properties
+CHARACTER (LEN=*), INTENT (IN) :: HOPWLW !cloud water LW optical properties
+CHARACTER (LEN=*), INTENT (IN) :: HOPILW !ice water LW optical properties
+REAL, INTENT(IN) :: PFUDG ! subgrid cloud inhomogenity factor
+REAL, DIMENSION(:,:), INTENT(IN) :: PCOSZEN ! COS(zenithal solar angle)
+REAL, INTENT(IN) :: PCORSOL ! SOLar constant CORrection
+REAL, DIMENSION(:,:), INTENT(IN) :: PSEA ! Land-sea mask
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PDIR_ALB! Surface direct ALBedo
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PSCA_ALB! Surface diffuse ALBedo
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PEMIS ! Surface IR EMISsivity
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PCLDFR ! CLouD FRaction
+REAL, INTENT(IN) :: PCCO2 ! CO2 content
+REAL, DIMENSION(:,:), INTENT(IN) :: PTSRAD ! RADiative Surface Temperature
+REAL, DIMENSION(:,:), INTENT(IN) :: PSTATM ! selected standard atmosphere
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PTHT ! THeta at t
+REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT ! moist variables at t (humidity, cloud water, rain water, ice water)
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! pressure at t
+REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT ! scalar variable ( C2R2 and C1R3 particle)
+!
+REAL, DIMENSION(:,:,:), POINTER :: POZON ! OZONE field from clim.
+REAL, DIMENSION(:,:,:,:), POINTER :: PAER ! AERosols optical thickness from clim.
+REAL, DIMENSION(:,:,:,:), POINTER :: PDST_WL ! AERosols Extinction by wavelength .
+REAL, DIMENSION(:,:,:,:), POINTER :: PAER_CLIM ! AERosols optical thickness from clim.
+ ! note : the vertical dimension of
+ ! these fields include the "radiation levels"
+ ! above domain top
+ !
+
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF ![kg/m3] air density
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ![m] height of layers
+
+INTEGER, DIMENSION(:,:), INTENT(INOUT) :: KCLEARCOL_TM1 ! trace of cloud/clear col
+ ! at the previous radiation step
+!
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PDTHRAD ! THeta RADiative Tendancy
+REAL, DIMENSION(:,:), INTENT(INOUT) :: PSRFLWD ! Downward SuRFace LW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSRFSWD_DIR ! Downward SuRFace SW Flux DIRect
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSRFSWD_DIF ! Downward SuRFace SW Flux DIFfuse
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSWU ! upward SW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSWD ! downward SW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLWU ! upward LW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLWD ! downward LW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PDTHRADSW ! dthrad sw
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PDTHRADLW ! dthradsw
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRADEFF ! effective radius
+
+END SUBROUTINE RADIATIONS_AGG
+!
+END INTERFACE
+!
+END MODULE MODI_RADIATIONS_AGG
+!
+! ############################################################################
+ SUBROUTINE RADIATIONS_AGG (KRAD_AGG,KI_RAD_AGG,KJ_RAD_AGG,KIOR_RAD_AGG,KJOR_RAD_AGG, &
+ TPFILE,OCLEAR_SKY,OCLOUD_ONLY, &
+ KCLEARCOL_TM1,HEFRADL,HEFRADI,HOPWSW,HOPISW,HOPWLW,HOPILW, &
+ PFUDG, KDLON, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD, &
+ KSWB_MNH,KLWB_MNH, KSTATM,KRAD_COLNBR,PCOSZEN,PSEA, PCORSOL, &
+ PDIR_ALB, PSCA_ALB,PEMIS, PCLDFR, PCCO2, PTSRAD, PSTATM, &
+ PTHT, PRT, PPABST, POZON, PAER, PDST_WL, PAER_CLIM, PSVT, &
+ PDTHRAD, PSRFLWD, PSRFSWD_DIR,PSRFSWD_DIF, PRHODREF, PZZ, &
+ PRADEFF, PSWU, PSWD, PLWU,PLWD, PDTHRADSW, PDTHRADLW )
+! ############################################################################
+!
+!!**** *RADIATIONS * - routine to call the SW and LW radiation calculations
+!!
+!! PURPOSE
+!! -------
+!! The purpose of this routine is to aggregate columns of the temperature, water vapor
+!! liquid water, cloud fraction, ozone profiles for the ECMWF radiation
+!! calculations. There is a great number of available radiative fluxes in
+!! the output, but only the potential temperature radiative tendency and the
+!! SW and LW surface fluxes are provided in the output of the routine.
+!!
+!!** METHOD
+!! ------
+!!
+!! All columns are aggregated according to NRAD_AGG * NRAD_AGG points.
+!!
+!!
+!! EXTERNAL
+!! --------
+!! Subroutine RADIATIONS
+!!
+!! IMPLICIT ARGUMENTS
+!! ------------------
+!!
+!! REFERENCE
+!! ---------
+!! Book2 of documentation ( routine RADIATIONS )
+!!
+!! AUTHOR
+!! ------
+!! V. Masson * CNRM *
+!!
+!! MODIFICATIONS
+!! -------------
+!! Original 27/10/23
+!-------------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+USE PARKIND1, ONLY: JPRB
+!
+USE MODD_CONF, ONLY: LCARTESIAN
+USE MODD_CST
+USE MODD_GRID , ONLY: XLAT0, XLON0
+USE MODD_GRID_n , ONLY: XLAT, XLON, XXHAT, XYHAT, XXHAT_ll, XYHAT_ll
+USE MODD_IO, ONLY: TFILEDATA
+USE MODD_LUNIT_n, ONLY: TLUOUT
+USE MODD_LBC_n, ONLY: CLBCX, CLBCY
+USE MODD_PARAMETERS
+!
+USE MODI_RADIATIONS
+USE MODE_ll
+use mode_msg
+USE MODE_SUM_ll, ONLY: MIN_ll
+USE MODE_SUM2_ll, ONLY: GMINLOC_ll
+USE MODD_VAR_ll, ONLY: IP
+USE MODD_ARGSLIST_ll, ONLY : LIST_ll
+!
+!
+IMPLICIT NONE
+!
+!* 0.1 DECLARATIONS OF DUMMY ARGUMENTS :
+!
+INTEGER, INTENT(IN) :: KRAD_AGG ! number of aggregated points
+INTEGER, INTENT(IN) :: KI_RAD_AGG ! reformatted X array size
+INTEGER, INTENT(IN) :: KJ_RAD_AGG ! reformatted Y array size
+INTEGER, INTENT(IN) :: KIOR_RAD_AGG ! index of first point of packed array according to current domain
+INTEGER, INTENT(IN) :: KJOR_RAD_AGG ! index of first point of packed array according to current domain
+
+TYPE(TFILEDATA), INTENT(IN) :: TPFILE ! Output file
+LOGICAL, INTENT(IN) :: OCLOUD_ONLY! flag for the cloud column
+ ! computations only
+LOGICAL, INTENT(IN) :: OCLEAR_SKY !
+INTEGER, INTENT(IN) :: KDLON ! number of columns where the
+ ! radiation calculations are
+ ! performed
+INTEGER, INTENT(IN) :: KFLEV ! number of vertical levels
+ ! where the radiation
+ ! calculations are performed
+INTEGER, INTENT(IN) :: KRAD_DIAG ! index for the number of
+ ! fields in the output
+INTEGER, INTENT(IN) :: KFLUX ! number of top and ground
+ ! fluxes for the ZFLUX array
+INTEGER, INTENT(IN) :: KRAD ! number of satellite radiances
+ ! for the ZRAD and ZRADCS arrays
+INTEGER, INTENT(IN) :: KAER ! number of AERosol classes
+
+INTEGER, INTENT(IN) :: KSWB_OLD ! number of SW band ECMWF
+INTEGER, INTENT(IN) :: KSWB_MNH ! number of SW band ECRAD
+INTEGER, INTENT(IN) :: KLWB_MNH ! number of LW band ECRAD
+INTEGER, INTENT(IN) :: KSTATM ! index of the standard
+ ! atmosphere level just above
+ ! the model top
+INTEGER, INTENT(IN) :: KRAD_COLNBR ! factor by which the memory
+ ! is split
+ !
+ !Choice of :
+CHARACTER (LEN=*), INTENT (IN) :: HEFRADL !
+CHARACTER (LEN=*), INTENT (IN) :: HEFRADI !
+CHARACTER (LEN=*), INTENT (IN) :: HOPWSW !cloud water SW optical properties
+CHARACTER (LEN=*), INTENT (IN) :: HOPISW !ice water SW optical properties
+CHARACTER (LEN=*), INTENT (IN) :: HOPWLW !cloud water LW optical properties
+CHARACTER (LEN=*), INTENT (IN) :: HOPILW !ice water LW optical properties
+REAL, INTENT(IN) :: PFUDG ! subgrid cloud inhomogenity factor
+REAL, DIMENSION(:,:), INTENT(IN) :: PCOSZEN ! COS(zenithal solar angle)
+REAL, INTENT(IN) :: PCORSOL ! SOLar constant CORrection
+REAL, DIMENSION(:,:), INTENT(IN) :: PSEA ! Land-sea mask
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PDIR_ALB! Surface direct ALBedo
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PSCA_ALB! Surface diffuse ALBedo
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PEMIS ! Surface IR EMISsivity
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PCLDFR ! CLouD FRaction
+REAL, INTENT(IN) :: PCCO2 ! CO2 content
+REAL, DIMENSION(:,:), INTENT(IN) :: PTSRAD ! RADiative Surface Temperature
+REAL, DIMENSION(:,:), INTENT(IN) :: PSTATM ! selected standard atmosphere
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PTHT ! THeta at t
+REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PRT ! moist variables at t (humidity, cloud water, rain water, ice water)
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PPABST ! pressure at t
+REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PSVT ! scalar variable ( C2R2 and C1R3 particle)
+!
+REAL, DIMENSION(:,:,:), POINTER :: POZON ! OZONE field from clim.
+REAL, DIMENSION(:,:,:,:), POINTER :: PAER ! AERosols optical thickness from clim.
+REAL, DIMENSION(:,:,:,:), POINTER :: PDST_WL ! AERosols Extinction by wavelength .
+REAL, DIMENSION(:,:,:,:), POINTER :: PAER_CLIM ! AERosols optical thickness from clim.
+ ! note : the vertical dimension of
+ ! these fields include the "radiation levels"
+ ! above domain top
+ !
+
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRHODREF ![kg/m3] air density
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PZZ ![m] height of layers
+
+INTEGER, DIMENSION(:,:), INTENT(INOUT) :: KCLEARCOL_TM1 ! trace of cloud/clear col
+ ! at the previous radiation step
+!
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PDTHRAD ! THeta RADiative Tendancy
+REAL, DIMENSION(:,:), INTENT(INOUT) :: PSRFLWD ! Downward SuRFace LW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSRFSWD_DIR ! Downward SuRFace SW Flux DIRect
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSRFSWD_DIF ! Downward SuRFace SW Flux DIFfuse
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSWU ! upward SW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PSWD ! downward SW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLWU ! upward LW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLWD ! downward LW Flux
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PDTHRADSW ! dthrad sw
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PDTHRADLW ! dthradsw
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRADEFF ! effective radius
+!
+!
+!* 0.2 DECLARATIONS OF LOCAL VARIABLES
+!
+INTEGER :: IIB ! I index value of the first inner mass point
+INTEGER :: IJB ! J index value of the first inner mass point
+INTEGER :: IKB ! K index value of the first inner mass point
+INTEGER :: IIE ! I index value of the last inner mass point
+INTEGER :: IJE ! J index value of the last inner mass point
+INTEGER :: IKE ! K index value of the last inner mass point
+INTEGER :: IIU ! array size for the first index
+INTEGER :: IJU ! array size for the second index
+INTEGER :: IKU ! array size for the third index
+INTEGER :: IIMAX ! last X possible point to consider in aggregation
+INTEGER :: IJMAX ! last Y possible point to consider in aggregation
+!
+INTEGER :: JIP ! X packed array index
+INTEGER :: JJP ! Y packed array index
+INTEGER :: JI ! X full array index in current processor
+INTEGER :: JJ ! Y full array index in current processor
+INTEGER :: IXORP ! Index of left X point of packed domain being treated
+INTEGER :: IYORP ! Index of bottom Y point of packed domain being treated
+!
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ZLAT ! Latitude
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ZLON ! Longitude
+!
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ZCOSZEN ! COS(zenithal solar angle)
+REAL :: ZCORSOL ! SOLar constant CORrection
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ZSEA ! Land-sea mask
+!
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PDIR_ALB,3)) :: ZDIR_ALB! Surface direct ALBedo
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PSCA_ALB,3)) :: ZSCA_ALB! Surface diffuse ALBedo
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PEMIS,3)) :: ZEMIS ! Surface IR EMISsivity
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PCLDFR,3)) :: ZCLDFR ! CLouD FRaction
+REAL :: ZCCO2 ! CO2 content
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ZTSRAD ! RADiative Surface Temperature
+!
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PTHT,3)) :: ZTHT ! THeta at t
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PRT,3),SIZE(PRT,4)) :: ZRT ! moist variables at t (humidity, cloud water, rain water, ice water)
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PPABST,3)) :: ZPABST ! pressure at t
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PSVT,3),SIZE(PSVT,4)) :: ZSVT ! scalar variable ( C2R2 and C1R3 particle)
+!
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(POZON,3)) :: ZOZON ! OZONE field from clim.
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PAER,3),SIZE(PAER,4)) :: ZAER ! AERosols optical thickness from clim.
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PDST_WL,3),SIZE(PDST_WL,4)) :: ZDST_WL ! AERosols Extinction by wavelength .
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PAER_CLIM,3),SIZE(PAER_CLIM,4)) :: ZAER_CLIM ! AERosols optical thickness from clim.
+ ! note : the vertical dimension of
+ ! these fields include the "radiation levels"
+ ! above domain top
+ !
+
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PRHODREF,3)) :: ZRHODREF ![kg/m3] air density
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PZZ,3)) :: ZZZ ![m] height of layers
+
+INTEGER, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ICLEARCOL_TM1 ! trace of cloud/clear col
+ ! at the previous radiation step
+!
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PDTHRAD,3)) :: ZDTHRAD ! THeta RADiative Tendancy
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG) :: ZSRFLWD ! Downward SuRFace LW Flux
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PSRFSWD_DIR,3)) :: ZSRFSWD_DIR ! Downward SuRFace SW Flux DIRect
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PSRFSWD_DIF,3)) :: ZSRFSWD_DIF ! Downward SuRFace SW Flux DIFfuse
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PSWU,3)) :: ZSWU ! upward SW Flux
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PSWD,3)) :: ZSWD ! downward SW Flux
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PLWU,3)) :: ZLWU ! upward LW Flux
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PLWD,3)) :: ZLWD ! downward LW Flux
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PDTHRADSW,3)) :: ZDTHRADSW ! dthrad sw
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PDTHRADLW,3)) :: ZDTHRADLW ! dthradsw
+REAL, DIMENSION(KI_RAD_AGG,KJ_RAD_AGG,SIZE(PRADEFF,3)) :: ZRADEFF ! effective radius
+!
+
+REAL, DIMENSION(SIZE(PCOSZEN,1),SIZE(PCOSZEN,2)) :: ZZLAT ! Latitude
+REAL, DIMENSION(SIZE(PCOSZEN,1),SIZE(PCOSZEN,2)) :: ZZLON ! Longitude
+!
+!
+REAL, DIMENSION(SIZE(PTHT,1),SIZE(PTHT,2),SIZE(PTHT,3)) :: ZDZPABST
+INTEGER :: ILUOUT ! Logical unit number for output-listing
+REAL :: ZMINVAL
+INTEGER, DIMENSION(3) :: IMINLOC
+INTEGER :: IINFO_ll
+!-------------------------------------------------------------------------
+!-------------------------------------------------------------------------
+!-------------------------------------------------------------------------
+!
+!* 1. COMPUTE DIMENSIONS OF ARRAYS AND OTHER INDICES
+! ----------------------------------------------
+!
+! full arrays
+CALL GET_DIM_EXT_ll ('B',IIU,IJU)
+CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
+!
+
+! Limit of points to integrate into the aggregation on X and Y limits
+IIMAX = IIU ; IF (LEAST_ll() .AND. CLBCX(2)/='CYCL') IIMAX = IIE
+IJMAX = IJU ; IF (LNORTH_ll() .AND. CLBCY(2)/='CYCL') IJMAX = IJE
+!
+IKU = SIZE(PTHT,3)
+IKB = 1 + JPVEXT
+IKE = IKU - JPVEXT
+!
+!
+!* 1.1 CHECK PRESSURE DECREASING
+! -------------------------
+ZDZPABST(:,:,1:IKU-1) = PPABST(:,:,1:IKU-1) - PPABST(:,:,2:IKU)
+ZDZPABST(:,:,IKU) = ZDZPABST(:,:,IKU-1)
+!
+ZMINVAL=MIN_ll(ZDZPABST,IINFO_ll)
+!
+IF ( ZMINVAL <= 0.0 ) THEN
+ ILUOUT = TLUOUT%NLU
+ IMINLOC=GMINLOC_ll( ZDZPABST )
+ WRITE(ILUOUT,*) ' radiation.f90 STOP :: SOMETHING WRONG WITH PRESSURE , ZDZPABST <= 0.0 '
+ WRITE(ILUOUT,*) ' radiation :: ZDZPABST ', ZMINVAL,' located at ', IMINLOC
+ FLUSH(unit=ILUOUT)
+ call Print_msg( NVERB_FATAL, 'GEN', 'RADIATIONS', 'something wrong with pressure: ZDZPABST <= 0.0' )
+END IF
+!
+!-------------------------------------------------------------------------------
+!
+! special case of no packing: avoids transfer of data from one array to another (but duplicates the code to call the radaitions routine)
+! -------------------------
+!
+!--------------------
+IF (KRAD_AGG==1) THEN
+!--------------------
+
+ IF(LCARTESIAN) THEN
+ ZZLAT(:,:) = XLAT0 *(XPI/180.)
+ ZZLON(:,:) = XLON0 *(XPI/180.)
+ ELSE
+ ZZLAT = XLAT *(XPI/180.)
+ ZZLON = XLON *(XPI/180.)
+ END IF
+!
+ CALL RADIATIONS( TPFILE, IIB, IIE,IJB,IJE, &
+ OCLEAR_SKY, OCLOUD_ONLY, KCLEARCOL_TM1, HEFRADL, HEFRADI, HOPWSW, HOPISW, &
+ HOPWLW, HOPILW, PFUDG, &
+ KDLON, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD, KSWB_MNH, KLWB_MNH, &
+ KSTATM, KRAD_COLNBR, PCOSZEN, PSEA, PCORSOL, ZZLAT, ZZLON, &
+ PDIR_ALB, PSCA_ALB, PEMIS, PCLDFR, PCCO2, PTSRAD, PSTATM, PTHT, PRT, &
+ PPABST, POZON, PAER,PDST_WL, PAER_CLIM, PSVT, &
+ PDTHRAD, PSRFLWD, PSRFSWD_DIR, PSRFSWD_DIF, PRHODREF, PZZ , &
+ PRADEFF, PSWU, PSWD, PLWU, PLWD, PDTHRADSW, PDTHRADLW )
+ RETURN
+!--------------------
+END IF
+!--------------------
+!
+!-------------------------------------------------------------------------------
+!
+! 2. Packs arrays
+! ------------
+!
+! latitude and longitude definition & packing
+!
+! Only the middle point of each packed subdomain is kept, there is nopt averaging
+! (to avoid issues at poles and longitude sign change).
+!
+IF(LCARTESIAN) THEN
+ ZLAT(:,:) = XLAT0
+ ZLON(:,:) = XLON0
+ELSE
+ CALL PACK_RAD_AGG_MID(XLAT,ZLAT,.FALSE.)
+ CALL PACK_RAD_AGG_MID(XLON,ZLON,.FALSE.)
+END IF
+!
+ZLAT = ZLAT *(XPI/180.)
+ZLON = ZLON *(XPI/180.)
+!
+!
+! packing of input fileds
+!
+CALL PACK_RAD_AGG_2D(PCOSZEN,ZCOSZEN,.TRUE.)
+CALL PACK_RAD_AGG_2D(PSEA,ZSEA,.TRUE.)
+CALL PACK_RAD_AGG_3D(PDIR_ALB,ZDIR_ALB,.TRUE.)
+CALL PACK_RAD_AGG_3D(PSCA_ALB,ZSCA_ALB,.TRUE.)
+CALL PACK_RAD_AGG_3D(PEMIS,ZEMIS,.TRUE.)
+CALL PACK_RAD_AGG_3D(PCLDFR,ZCLDFR,.FALSE.)
+CALL PACK_RAD_AGG_2D(PTSRAD,ZTSRAD,.TRUE.)
+CALL PACK_RAD_AGG_3D(PTHT,ZTHT,.FALSE.)
+CALL PACK_RAD_AGG_4D(PRT,ZRT,.FALSE.)
+CALL PACK_RAD_AGG_3D(PPABST,ZPABST,.FALSE.)
+CALL PACK_RAD_AGG_4D(PSVT,ZSVT,.FALSE.)
+CALL PACK_RAD_AGG_3D(POZON,ZOZON,.TRUE.)
+CALL PACK_RAD_AGG_4D(PAER,ZAER,.TRUE.)
+CALL PACK_RAD_AGG_4D(PDST_WL,ZDST_WL,.TRUE.)
+CALL PACK_RAD_AGG_4D(PAER_CLIM,ZAER_CLIM,.TRUE.)
+CALL PACK_RAD_AGG_3D(PRHODREF,ZRHODREF,.FALSE.)
+CALL PACK_RAD_AGG_3D(PZZ,ZZZ,.FALSE.)
+CALL PACK_RAD_AGG_I2(KCLEARCOL_TM1,ICLEARCOL_TM1,.TRUE.)
+CALL PACK_RAD_AGG_3D(PDTHRAD,ZDTHRAD,.TRUE.)
+CALL PACK_RAD_AGG_2D(PSRFLWD,ZSRFLWD,.TRUE.)
+CALL PACK_RAD_AGG_3D(PSRFSWD_DIR,ZSRFSWD_DIR,.TRUE.)
+CALL PACK_RAD_AGG_3D(PSRFSWD_DIF,ZSRFSWD_DIF,.TRUE.)
+CALL PACK_RAD_AGG_3D(PSWU,ZSWU,.TRUE.)
+CALL PACK_RAD_AGG_3D(PSWD,ZSWD,.TRUE.)
+CALL PACK_RAD_AGG_3D(PLWU,ZLWU,.TRUE.)
+CALL PACK_RAD_AGG_3D(PLWD,ZLWD,.TRUE.)
+CALL PACK_RAD_AGG_3D(PDTHRADSW,ZDTHRADSW,.TRUE.)
+CALL PACK_RAD_AGG_3D(PDTHRADLW,ZDTHRADLW,.TRUE.)
+CALL PACK_RAD_AGG_3D(PRADEFF,ZRADEFF,.TRUE.)
+!-------------------------------------------------------------------------------
+!
+!-------------------------------------------------------------------------------
+!
+! 3. Call radiation on packed columns
+! --------------------------------
+!
+ CALL RADIATIONS( TPFILE, 1, KI_RAD_AGG,1,KJ_RAD_AGG, &
+ OCLEAR_SKY, OCLOUD_ONLY, ICLEARCOL_TM1, HEFRADL, HEFRADI, HOPWSW, HOPISW, &
+ HOPWLW, HOPILW, PFUDG, &
+ KI_RAD_AGG*KJ_RAD_AGG, KFLEV, KRAD_DIAG, KFLUX, KRAD, KAER, KSWB_OLD, KSWB_MNH, KLWB_MNH, &
+ KSTATM, KRAD_COLNBR, ZCOSZEN, ZSEA, PCORSOL, ZLAT, ZLON, &
+ ZDIR_ALB, ZSCA_ALB, ZEMIS, ZCLDFR, PCCO2, ZTSRAD, PSTATM, ZTHT, ZRT, &
+ ZPABST, ZOZON, ZAER,ZDST_WL, ZAER_CLIM, ZSVT, &
+ ZDTHRAD, ZSRFLWD, ZSRFSWD_DIR, ZSRFSWD_DIF, ZRHODREF, ZZZ , &
+ ZRADEFF, ZSWU, ZSWD, ZLWU, ZLWD, ZDTHRADSW, ZDTHRADLW )
+!
+!-------------------------------------------------------------------------------
+!
+! 4. Unpacks arrays
+! --------------
+!
+CALL UNPACK_RAD_AGG_3D(PDTHRAD,ZDTHRAD)
+CALL UNPACK_RAD_AGG_2D(PSRFLWD,ZSRFLWD)
+CALL UNPACK_RAD_AGG_3D(PSRFSWD_DIR,ZSRFSWD_DIR)
+CALL UNPACK_RAD_AGG_3D(PSRFSWD_DIF,ZSRFSWD_DIF)
+CALL UNPACK_RAD_AGG_3D(PSWU,ZSWU)
+CALL UNPACK_RAD_AGG_3D(PSWD,ZSWD)
+CALL UNPACK_RAD_AGG_3D(PLWU,ZLWU)
+CALL UNPACK_RAD_AGG_3D(PLWD,ZLWD)
+CALL UNPACK_RAD_AGG_3D(PDTHRADSW,ZDTHRADSW)
+CALL UNPACK_RAD_AGG_3D(PDTHRADLW,ZDTHRADLW)
+CALL UNPACK_RAD_AGG_3D(PRADEFF,ZRADEFF)
+!
+!-------------------------------------------------------------------------------
+!
+CONTAINS
+!
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE PACK_RAD_AGG_2D(PFULL,PPACK,OEXCH)
+REAL, DIMENSION(:,:), INTENT(IN) :: PFULL
+REAL, DIMENSION(:,:), INTENT(OUT) :: PPACK
+LOGICAL :: OEXCH
+INTEGER :: ICOUNT
+TYPE(LIST_ll), POINTER :: TZFIELDS_ll ! list of fields to exchange
+REAL, DIMENSION(SIZE(PFULL,1),SIZE(PFULL,2)) :: ZFULL
+
+ NULLIFY(TZFIELDS_ll)
+ ZFULL = PFULL
+ IF (OEXCH) THEN
+ CALL ADD2DFIELD_ll( TZFIELDS_ll, ZFULL, 'RADIATION_AGG: PACK_RAD_AGG_2D' )
+ CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
+ CALL CLEANLIST_ll(TZFIELDS_ll)
+ END IF
+
+ PPACK = 0.
+
+ DO JJP=1,KJ_RAD_AGG
+ DO JIP=1,KI_RAD_AGG
+ ICOUNT = 0
+ IXORP = KIOR_RAD_AGG + (JIP-1) * KRAD_AGG
+ IYORP = KJOR_RAD_AGG + (JJP-1) * KRAD_AGG
+ DO JJ=IYORP,MIN(IYORP+KRAD_AGG-1,IJMAX)
+ DO JI=IXORP,MIN(IXORP+KRAD_AGG-1,IIMAX)
+ PPACK(JIP,JJP) = PPACK(JIP,JJP) + ZFULL(JI,JJ)
+ ICOUNT = ICOUNT + 1
+ END DO
+ END DO
+ PPACK(JIP,JJP) = PPACK(JIP,JJP) / ICOUNT
+ END DO
+ END DO
+END SUBROUTINE PACK_RAD_AGG_2D
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE PACK_RAD_AGG_3D(PFULL,PPACK,OEXCH)
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PFULL
+REAL, DIMENSION(:,:,:), INTENT(OUT) :: PPACK
+LOGICAL :: OEXCH
+INTEGER :: JK
+TYPE(LIST_ll), POINTER :: TZFIELDS_ll ! list of fields to exchange
+REAL, DIMENSION(SIZE(PFULL,1),SIZE(PFULL,2),SIZE(PFULL,3)) :: ZFULL
+
+NULLIFY(TZFIELDS_ll)
+ZFULL = PFULL
+IF (OEXCH) THEN
+ DO JK=1,SIZE(PFULL,3)
+ CALL ADD2DFIELD_ll( TZFIELDS_ll, ZFULL(:,:,JK), 'RADIATION_AGG: PACK_RAD_AGG_3D' )
+ END DO
+ CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
+ CALL CLEANLIST_ll(TZFIELDS_ll)
+END IF
+
+DO JK=1,SIZE(PFULL,3)
+ CALL PACK_RAD_AGG_2D(ZFULL(:,:,JK),PPACK(:,:,JK),.FALSE.)
+END DO
+
+END SUBROUTINE PACK_RAD_AGG_3D
+
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE PACK_RAD_AGG_4D(PFULL,PPACK,OEXCH)
+REAL, DIMENSION(:,:,:,:), INTENT(IN) :: PFULL
+REAL, DIMENSION(:,:,:,:), INTENT(OUT) :: PPACK
+LOGICAL :: OEXCH
+INTEGER :: JK, JL
+TYPE(LIST_ll), POINTER :: TZFIELDS_ll ! list of fields to exchange
+REAL, DIMENSION(SIZE(PFULL,1),SIZE(PFULL,2),SIZE(PFULL,3),SIZE(PFULL,4)) :: ZFULL
+
+NULLIFY(TZFIELDS_ll)
+ZFULL = PFULL
+IF (OEXCH) THEN
+ DO JL=1,SIZE(PFULL,4)
+ DO JK=1,SIZE(PFULL,3)
+ CALL ADD2DFIELD_ll( TZFIELDS_ll, ZFULL(:,:,JK,JL), 'RADIATION_AGG: PACK_RAD_AGG_4D' )
+ END DO
+ END DO
+ CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
+ CALL CLEANLIST_ll(TZFIELDS_ll)
+END IF
+
+DO JL=1,SIZE(PFULL,4)
+ DO JK=1,SIZE(PFULL,3)
+ CALL PACK_RAD_AGG_2D(ZFULL(:,:,JK,JL),PPACK(:,:,JK,JL),.FALSE.)
+ END DO
+END DO
+
+END SUBROUTINE PACK_RAD_AGG_4D
+
+
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE PACK_RAD_AGG_I2(KFULL,KPACK,OEXCH)
+INTEGER, DIMENSION(:,:), INTENT(IN) :: KFULL
+INTEGER, DIMENSION(:,:), INTENT(OUT) :: KPACK
+LOGICAL :: OEXCH
+REAL, DIMENSION(SIZE(KFULL,1),SIZE(KFULL,2)) :: ZFULL
+REAL, DIMENSION(SIZE(KPACK,1),SIZE(KPACK,2)) :: ZPACK
+
+ZFULL = FLOAT(KFULL)
+CALL PACK_RAD_AGG_2D(ZFULL(:,:),ZPACK(:,:),OEXCH)
+!
+KPACK=0
+WHERE ( ZPACK>=1.-0.5/FLOAT(KRAD_AGG*KRAD_AGG) ) KPACK = 1
+! This ensures that the output is coherent with the clear_sky columns definitions in radiations:
+! If the averaged value is 1 (allowing calculation error), this means that all columns are clear, so we keep this information as clear
+! If any column was not clear, then the averaged value is below one
+! (in fact, necessarily below 1-1/FLOAT(KRAD_AGG**2), because entry values are 0 or 1 ),
+! and then the packed column is partly cloudy, so we keep the value zero.
+!
+END SUBROUTINE PACK_RAD_AGG_I2
+
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE PACK_RAD_AGG_MID(PFULL,PPACK,OEXCH)
+REAL, DIMENSION(:,:), INTENT(IN) :: PFULL
+REAL, DIMENSION(:,:), INTENT(OUT) :: PPACK
+LOGICAL :: OEXCH
+
+! DO JJP=1,MIN(KJ_RAD_AGG,IJMAX)
+! DO JIP=1,MIN(KI_RAD_AGG,IIMAX)
+ DO JJP=1,KJ_RAD_AGG
+ DO JIP=1,KI_RAD_AGG
+ IXORP = KIOR_RAD_AGG + (JIP-1) * KRAD_AGG
+ IYORP = KJOR_RAD_AGG + (JJP-1) * KRAD_AGG
+ PPACK(JIP,JJP) = PFULL(MIN(IXORP + KRAD_AGG/2,IImax),MIN(IYORP + KRAD_AGG/2,Ijmax) )
+ END DO
+ END DO
+END SUBROUTINE PACK_RAD_AGG_MID
+!
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE UNPACK_RAD_AGG_2D(PFULL,PPACK)
+REAL, DIMENSION(:,:), INTENT(OUT) :: PFULL
+REAL, DIMENSION(:,:), INTENT(IN) :: PPACK
+
+ DO JJP=1,KJ_RAD_AGG
+ DO JIP=1,KI_RAD_AGG
+ IXORP = KIOR_RAD_AGG + (JIP-1) * KRAD_AGG
+ IYORP = KJOR_RAD_AGG + (JJP-1) * KRAD_AGG
+ DO JJ=IYORP,MIN(IYORP+KRAD_AGG-1,IJMAX)
+ DO JI=IXORP,MIN(IXORP+KRAD_AGG-1,IIMAX)
+ PFULL(JI,JJ) = PPACK(JIP,JJP)
+ END DO
+ END DO
+ END DO
+ END DO
+
+END SUBROUTINE UNPACK_RAD_AGG_2D
+!-------------------------------------------------------------------------------
+!
+SUBROUTINE UNPACK_RAD_AGG_3D(PFULL,PPACK)
+REAL, DIMENSION(:,:,:), INTENT(OUT) :: PFULL
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PPACK
+INTEGER :: JK
+
+DO JK=1,SIZE(PFULL,3)
+ CALL UNPACK_RAD_AGG_2D(PFULL(:,:,JK),PPACK(:,:,JK))
+END DO
+
+END SUBROUTINE UNPACK_RAD_AGG_3D
+
+
+!-------------------------------------------------------------------------------
+!
+END SUBROUTINE RADIATIONS_AGG
+
diff --git a/src/MNH/read_exsegn.f90 b/src/MNH/read_exsegn.f90
index cab38bf98..8fed63820 100644
--- a/src/MNH/read_exsegn.f90
+++ b/src/MNH/read_exsegn.f90
@@ -2251,6 +2251,25 @@ IF ( CRAD /= 'NONE' .AND. CPROGRAM=='MESONH' ) THEN
WRITE(UNIT=ILUOUT,FMT=*) 'YOU MUST USE LAERO_FT=T WITH CAER=TEGE IF CCONF=RESTA IN ALL SEGMENTS'
WRITE(UNIT=ILUOUT,FMT=*) 'TO UPDATE THE OZONE AND AEROSOLS CLIMATOLOGY USED BY THE RADIATION CODE;'
END IF
+!
+! impossible to have aggregation of radiative columns and LCEAR_SKy or CLOUD_ONLY options
+! because of parallelization issues when aggregating all clear sky columns between processors
+!
+! Please note that XDTRAD_CLONLY and XDTRAD are both initialized per default to 60. sec.
+! This is because the default time step is 60. sec, which is almost never the case.
+! This should be corrected (put a XUNDEF value for default for XDTRAD_CLONLY would be nice)
+!
+ IF (NRAD_AGG>1 .AND. (LCLEAR_SKY .OR. (XDTRAD<XDTRAD_CLONLY) ) ) THEN
+ WRITE(UNIT=ILUOUT,FMT=9003) KMI
+ WRITE(UNIT=ILUOUT,FMT=*) 'AGGREGATION OF RADIATIVE COLUMNS CANNOT BE DONE IF LCLEAR_SKY OPTION,'
+ WRITE(UNIT=ILUOUT,FMT=*) 'OR CLOUD_ONLY OPTION ARE ACTIVATED'
+ WRITE(UNIT=ILUOUT,FMT=*) 'AGGREGATION OF RADIATIVE COLUMNS : NRAD_AGG =', NRAD_AGG
+ WRITE(UNIT=ILUOUT,FMT=*) 'LCLEAR_SKY OPTION : LCLEAR_SKY =', LCLEAR_SKY
+ WRITE(UNIT=ILUOUT,FMT=*) 'CLOUD_ONLY OPTION, XDTRAD /= XDTRAD_CLONLY : '
+ WRITE(UNIT=ILUOUT,FMT=*) 'CLOUD_ONLY OPTION : XDTRAD =', XDTRAD
+ WRITE(UNIT=ILUOUT,FMT=*) 'CLOUD_ONLY OPTION : XDTRAD_CLONLY =', XDTRAD_CLONLY
+ CALL PRINT_MSG(NVERB_FATAL,'GEN','READ_EXSEG_n','NRAD_AGG OPTION NOT IMPLEMENTED FOR LCLEAR_SKY AND CLOUD_ONLY OPTIONS')
+ END IF
END IF
!
! 3.6 check the initialization of the deep convection scheme
--
GitLab