diff --git a/src/LIB/SURCOUCHE/src/mode_field.f90 b/src/LIB/SURCOUCHE/src/mode_field.f90
index a0df6d83f884d7c1e18f0f206f64e06440115d28..71eefd21106cbbf7c0ff13b7881697d371173faa 100644
--- a/src/LIB/SURCOUCHE/src/mode_field.f90
+++ b/src/LIB/SURCOUCHE/src/mode_field.f90
@@ -511,6 +511,30 @@ call Add_field2list( TFIELDDATA( &
   NDIMS      = 1,                &
   LTIMEDEP   = .FALSE.           ) )
 
+call Add_field2list( TFIELDDATA( &
+  CMNHNAME   = 'HAT_BOUND', &
+  CSTDNAME   = '',               &
+  CLONGNAME  = 'HAT_BOUND', &
+  CUNITS     = 'm',              &
+  CDIR       = '--',             &
+  CCOMMENT   = 'Boundaries of domain in the conformal or cartesian plane at u and v points', &
+  NGRID      = 0,                &
+  NTYPE      = TYPEREAL,         &
+  NDIMS      = 1,                &
+  LTIMEDEP   = .FALSE.           ) )
+
+call Add_field2list( TFIELDDATA(  &
+  CMNHNAME   = 'HATM_BOUND', &
+  CSTDNAME   = '',                &
+  CLONGNAME  = 'HATM_BOUND', &
+  CUNITS     = 'm',               &
+  CDIR       = '--',              &
+  CCOMMENT   = 'Boundaries of domain in the conformal or cartesian plane at mass points', &
+  NGRID      = 0,                 &
+  NTYPE      = TYPEREAL,          &
+  NDIMS      = 1,                 &
+  LTIMEDEP   = .FALSE.            ) )
+
 call Add_field2list( TFIELDDATA(      &
   CMNHNAME   = 'ZTOP',                &
   CSTDNAME   = 'altitude_at_top_of_atmosphere_model', &
@@ -3632,6 +3656,8 @@ call Goto_model_1field( 'XHATM', kfrom, kto, XXHATM )
 call Goto_model_1field( 'YHATM', kfrom, kto, XYHATM )
 call Goto_model_1field( 'ZHAT',  kfrom, kto, XZHAT  )
 call Goto_model_1field( 'ZHATM', kfrom, kto, XZHATM )
+call Goto_model_1field( 'HAT_BOUND',  kfrom, kto, XHAT_BOUND  )
+call Goto_model_1field( 'HATM_BOUND', kfrom, kto, XHATM_BOUND )
 call Goto_model_1field( 'ZTOP',  kfrom, kto, XZTOP  )
 call Goto_model_1field( 'DXHAT', kfrom, kto, XDXHAT )
 call Goto_model_1field( 'DYHAT', kfrom, kto, XDYHAT )
diff --git a/src/MNH/ini_modeln.f90 b/src/MNH/ini_modeln.f90
index 9ef32c017305e261ea7a5f88c259e9d920e1d725..a3de7b35a14b3b5c4ac30d3999d57c32a64ae7c2 100644
--- a/src/MNH/ini_modeln.f90
+++ b/src/MNH/ini_modeln.f90
@@ -292,7 +292,7 @@ END MODULE MODI_INI_MODEL_n
 !  S. Riette      04/2020: XHL* fields
 !  F. Auguste     02/2021: add IBM
 !  T.Nigel        02/2021: add turbulence recycling
-! J.L.Redelsperger 06/2011: OCEAN case
+! J.L.Redelsperger 06/2021: OCEAN case
 !---------------------------------------------------------------------------------
 !
 !*       0.    DECLARATIONS
@@ -1013,6 +1013,8 @@ ALLOCATE(XDIRCOSXW(IIU,IJU))
 ALLOCATE(XDIRCOSYW(IIU,IJU))
 ALLOCATE(XCOSSLOPE(IIU,IJU))
 ALLOCATE(XSINSLOPE(IIU,IJU))
+ALLOCATE(XHAT_BOUND (NHAT_BOUND_SIZE))
+ALLOCATE(XHATM_BOUND(NHAT_BOUND_SIZE))
 !
 ALLOCATE(XDXX(IIU,IJU,IKU))
 ALLOCATE(XDYY(IIU,IJU,IKU))
@@ -1858,9 +1860,10 @@ END IF
 CALL SET_GRID( KMI, TPINIFILE, IKU, NIMAX_ll, NJMAX_ll,                        &
                XTSTEP, XSEGLEN,                                                &
                XLONORI, XLATORI, XLON, XLAT,                                   &
-               XXHAT, XYHAT, XDXHAT, XDYHAT, XXHATM, XYHATM, XMAP,             &
-               XZS, XZZ, XZHAT, XZHATM, XZTOP, LSLEVE, XLEN1, XLEN2, XZSMT,    &
-               ZJ,                                                             &
+               XXHAT, XYHAT, XDXHAT, XDYHAT, XXHATM, XYHATM,                   &
+               XHAT_BOUND, XHATM_BOUND,                                        &
+               XMAP, XZS, XZZ, XZHAT, XZHATM, XZTOP, LSLEVE,                   &
+               XLEN1, XLEN2, XZSMT, ZJ,                                        &
                TDTMOD, TDTCUR, NSTOP, NBAK_NUMB, NOUT_NUMB, TBACKUPN, TOUTPUTN )
 !
 CALL METRICS(XMAP,XDXHAT,XDYHAT,XZZ,XDXX,XDYY,XDZX,XDZY,XDZZ)
diff --git a/src/MNH/ini_spectren.f90 b/src/MNH/ini_spectren.f90
index c17715ec618e5a4d91a68047c5c881c1890515a6..a98810782ef1e36442b8fd4bfb6e60150676cbcd 100644
--- a/src/MNH/ini_spectren.f90
+++ b/src/MNH/ini_spectren.f90
@@ -269,6 +269,8 @@ ALLOCATE(XZSMT(IIU,IJU))
 ALLOCATE(XZZ(IIU,IJU,IKU))
 ALLOCATE(XZHAT(IKU))
 ALLOCATE(XZHATM(IKU))
+ALLOCATE(XHAT_BOUND (NHAT_BOUND_SIZE))
+ALLOCATE(XHATM_BOUND(NHAT_BOUND_SIZE))
 !
 ALLOCATE(XDXX(IIU,IJU,IKU))
 ALLOCATE(XDYY(IIU,IJU,IKU))
@@ -688,9 +690,10 @@ CALL INI_BIKHARDT_n (NDXRATIO_ALL(KMI),NDYRATIO_ALL(KMI),KMI)
 CALL SET_GRID( KMI, TPINIFILE, IKU, NIMAX_ll, NJMAX_ll,                        &
                XTSTEP, XSEGLEN,                                                &
                XLONORI, XLATORI, XLON, XLAT,                                   &
-               XXHAT, XYHAT, XDXHAT, XDYHAT, XXHATM, XYHATM, XMAP,             &
-               XZS, XZZ, XZHAT, XZHATM, XZTOP, LSLEVE, XLEN1, XLEN2, XZSMT,    &
-               ZJ,                                                             &
+               XXHAT, XYHAT, XDXHAT, XDYHAT, XXHATM, XYHATM,                   &
+               XHAT_BOUND, XHATM_BOUND,                                        &
+               XMAP, XZS, XZZ, XZHAT, XZHATM, XZTOP, LSLEVE,                   &
+               XLEN1, XLEN2, XZSMT, ZJ,                                        &
                TDTMOD, TDTCUR, NSTOP, NBAK_NUMB, NOUT_NUMB, TBACKUPN, TOUTPUTN )
 !
 CALL METRICS(XMAP,XDXHAT,XDYHAT,XZZ,XDXX,XDYY,XDZX,XDZY,XDZZ)
diff --git a/src/MNH/modd_gridn.f90 b/src/MNH/modd_gridn.f90
index 6437cfc8bcff5cdbfcfc63c283990283b7c8a06b..1f0c5d895d330437f0578dd5e6f598aac62a2ff2 100644
--- a/src/MNH/modd_gridn.f90
+++ b/src/MNH/modd_gridn.f90
@@ -31,20 +31,39 @@
 !!
 !!    MODIFICATIONS
 !!    -------------
-!!      Original    05/05/94                      
-!!      J. Stein    15/11/95  add the slope angle
-!!      V. Ducrocq   13/08/98  // : add XLATOR_ll and XLONOR_ll       
-!!      V. Masson   nov 2004  supress XLATOR,XLONOR,XLATOR_ll,XLONOR_ll
-!!  Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
-!  P. Wautelet 02/09/2022: add XXHATM, XYHATM and XZHATM
+!  Original    05/05/94
+!  J. Stein    15/11/95:  add the slope angle
+!  V. Ducrocq  13/08/98: //: add XLATOR_ll and XLONOR_ll
+!  V. Masson      11/2004: supress XLATOR, XLONOR, XLATOR_ll, XLONOR_ll
+!  P. Wautelet 05/2016-04/2018: new data structures and calls for I/O
+!  P. Wautelet    09/2022: add XXHATM, XYHATM, XZHATM, XHAT_BOUND and XHATM_BOUND
 !-------------------------------------------------------------------------------
 !
 !*       0.   DECLARATIONS
 !             ------------
 !
 USE MODD_PARAMETERS, ONLY: JPMODELMAX
+
 IMPLICIT NONE
 
+SAVE
+
+! Parameters for XHAT_BOUND and XHATM_BOUND
+INTEGER, PARAMETER :: NHAT_BOUND_SIZE = 12
+INTEGER, PARAMETER :: NPHYS_XMIN = 1  ! Position of minimum position in physical domain in x direction
+INTEGER, PARAMETER :: NPHYS_XMAX = 2  ! Position of maximum position in physical domain in x direction
+INTEGER, PARAMETER :: NPHYS_YMIN = 3  ! Position of minimum position in physical domain in y direction
+INTEGER, PARAMETER :: NPHYS_YMAX = 4  ! Position of maximum position in physical domain in y direction
+INTEGER, PARAMETER :: NPHYS_ZMIN = 5  ! Position of minimum position in physical domain in z direction
+INTEGER, PARAMETER :: NPHYS_ZMAX = 6  ! Position of maximum position in physical domain in z direction
+INTEGER, PARAMETER :: NEXTE_XMIN = 7  ! Position of minimum position in extended domain in x direction
+INTEGER, PARAMETER :: NEXTE_XMAX = 8  ! Position of maximum position in extended domain in x direction
+INTEGER, PARAMETER :: NEXTE_YMIN = 9  ! Position of minimum position in extended domain in y direction
+INTEGER, PARAMETER :: NEXTE_YMAX = 10 ! Position of maximum position in extended domain in y direction
+INTEGER, PARAMETER :: NEXTE_ZMIN = 11 ! Position of minimum position in extended domain in z direction
+INTEGER, PARAMETER :: NEXTE_ZMAX = 12 ! Position of maximum position in extended domain in z direction
+
+
 REAL, DIMENSION(:,:),  POINTER :: XLON=>NULL(),XLAT=>NULL() ! Longitude and latitude  
 REAL, DIMENSION(:),    POINTER :: XXHAT=>NULL()             ! Position x in the conformal or cartesian plane
 REAL, DIMENSION(:),    POINTER :: XYHAT=>NULL()             ! Position y in the conformal or cartesian plane
@@ -67,5 +86,7 @@ LOGICAL,               POINTER  :: LSLEVE=>NULL()            ! Logical for SLEVE
 REAL,                  POINTER  :: XLEN1=>NULL()             ! Decay scale for smooth topography
 REAL,                  POINTER  :: XLEN2=>NULL()             ! Decay scale for small-scale topography deviation
 REAL, DIMENSION(:,:),  POINTER  :: XZSMT=>NULL()             ! smooth orography for SLEVE coordinate
+REAL, DIMENSION(:),    POINTER :: XHAT_BOUND  => NULL() ! Boundaries of global domain at u and v points
+REAL, DIMENSION(:),    POINTER :: XHATM_BOUND => NULL() ! Boundaries of global domain at mass points
 
 END MODULE MODD_GRID_n
diff --git a/src/MNH/prep_ideal_case.f90 b/src/MNH/prep_ideal_case.f90
index 6be3241c1ab29860a2f18525367c3374d245c5db..8d4cdc38599eb2b3b185c3f41762975ccbdbc5a2 100644
--- a/src/MNH/prep_ideal_case.f90
+++ b/src/MNH/prep_ideal_case.f90
@@ -379,7 +379,7 @@ USE MODE_ll
 USE MODE_MODELN_HANDLER
 use mode_field,            only: Alloc_field_scalars, Ini_field_list, Ini_field_scalars
 USE MODE_MSG
-USE MODE_SET_GRID,         only: INTERP_HORGRID_TO_MASSPOINTS
+USE MODE_SET_GRID,         only: INTERP_HORGRID_TO_MASSPOINTS, STORE_HORGRID_BOUNDS
 !
 USE MODI_DEFAULT_DESFM_n    ! Interface modules
 USE MODI_DEFAULT_EXPRE
@@ -1242,6 +1242,7 @@ ELSE
 !
   ALLOCATE( XXHAT(NIU),  XYHAT(NJU)  )
   ALLOCATE( XXHATM(NIU), XYHATM(NJU) )
+  ALLOCATE( XHAT_BOUND (NHAT_BOUND_SIZE), XHATM_BOUND(NHAT_BOUND_SIZE) )
 !
 ! define the grid localization at the earth surface by the central point
 ! coordinates
@@ -1254,14 +1255,10 @@ ELSE
 ! conformal coordinates (0,0). This is to allow the centering of the model in
 ! a non-cyclic  configuration regarding to XLATCEN or XLONCEN.
 !
-      ALLOCATE(ZXHAT_ll(NIMAX_ll+2*JPHEXT),ZYHAT_ll(NJMAX_ll+2*JPHEXT))
-      ZXHAT_ll=0.
-      ZYHAT_ll=0.
       CALL SM_LATLON(XLATCEN,XLONCEN,                     &
                        -XDELTAX*(NIMAX_ll/2-0.5+JPHEXT),  &
                        -XDELTAY*(NJMAX_ll/2-0.5+JPHEXT),  &
                        XLATORI,XLONORI)
-        DEALLOCATE(ZXHAT_ll,ZYHAT_ll)
 !
       WRITE(NLUOUT,FMT=*) 'PREP_IDEAL_CASE : XLATORI=' , XLATORI, &
                           ' XLONORI= ', XLONORI
@@ -1288,6 +1285,11 @@ ELSE
 
   ! Interpolations of positions to mass points
   CALL INTERP_HORGRID_TO_MASSPOINTS( XXHAT, XYHAT, XXHATM, XYHATM )
+
+  ! Collect global domain boundaries
+  CALL STORE_HORGRID_BOUNDS( XXHAT,  XYHAT,  XHAT_BOUND  )
+  CALL STORE_HORGRID_BOUNDS( XXHATM, XYHATM, XHATM_BOUND )
+
 END IF
 !
 !*       5.1.2  Orography and Gal-Chen Sommerville transformation :
diff --git a/src/MNH/read_hgridn.f90 b/src/MNH/read_hgridn.f90
index 39586f9795394e7fa9ce49a917b2a306473c789d..7bdf08ec4db7f3a70fb3cbd3a2cbb99ce17c0e5d 100644
--- a/src/MNH/read_hgridn.f90
+++ b/src/MNH/read_hgridn.f90
@@ -93,7 +93,7 @@ USE MODE_IO,            only: IO_Pack_set
 USE MODE_IO_FIELD_READ, only: IO_Field_read
 USE MODE_MSG
 USE MODE_MODELN_HANDLER
-USE MODE_SET_GRID,       only: INTERP_HORGRID_TO_MASSPOINTS
+USE MODE_SET_GRID,       only: INTERP_HORGRID_TO_MASSPOINTS, STORE_HORGRID_BOUNDS
 use MODE_TOOLS_ll,       only: GET_DIM_EXT_ll, GET_DIM_PHYS_ll, GET_INDICE_ll
 !
 IMPLICIT NONE
@@ -254,10 +254,16 @@ CALL IO_Field_read(TPFMFILE,'YHAT',XYHAT)
 
 IF ( .NOT. ASSOCIATED(XXHATM) ) ALLOCATE( XXHATM(SIZE( XXHAT )) )
 IF ( .NOT. ASSOCIATED(XYHATM) ) ALLOCATE( XYHATM(SIZE( XYHAT )) )
+IF ( .NOT. ASSOCIATED(XHAT_BOUND)  ) ALLOCATE( XHAT_BOUND (NHAT_BOUND_SIZE) )
+IF ( .NOT. ASSOCIATED(XHATM_BOUND) ) ALLOCATE( XHATM_BOUND(NHAT_BOUND_SIZE) )
 
 ! Interpolations of positions to mass points
 CALL INTERP_HORGRID_TO_MASSPOINTS( XXHAT, XYHAT, XXHATM, XYHATM )
 
+! Collect global domain boundaries
+CALL STORE_HORGRID_BOUNDS( XXHAT,  XYHAT,  XHAT_BOUND  )
+CALL STORE_HORGRID_BOUNDS( XXHATM, XYHATM, XHATM_BOUND )
+
 !JUAN REALZ
 IF ( CPROGRAM .EQ. "REAL  " ) THEN
 IF (.NOT. (ASSOCIATED(XZS))) ALLOCATE(XZS(IIU,IJU))
diff --git a/src/MNH/read_ver_grid.f90 b/src/MNH/read_ver_grid.f90
index 653bbec0203753a923431ebdea30a16a6a06e093..7fd47bf64b613692b1a085a4757fe8563d1d54d3 100644
--- a/src/MNH/read_ver_grid.f90
+++ b/src/MNH/read_ver_grid.f90
@@ -111,7 +111,7 @@ USE MODD_PARAMETERS
 !
 USE MODE_MSG
 USE MODE_POS
-USE MODE_SET_GRID, ONLY: INTERP_VERGRID_TO_MASSPOINTS
+USE MODE_SET_GRID, ONLY: INTERP_VERGRID_TO_MASSPOINTS, STORE_VERGRID_BOUNDS
 !
 USE MODI_DEFAULT_SLEVE
 !
@@ -331,6 +331,10 @@ XZTOP = XZHAT(IKU)
 ! Interpolations of positions to mass points
 CALL INTERP_VERGRID_TO_MASSPOINTS( XZHAT, XZHATM )
 
+! Collect global domain boundaries
+CALL STORE_VERGRID_BOUNDS( XZHAT,  XHAT_BOUND  )
+CALL STORE_VERGRID_BOUNDS( XZHATM, XHATM_BOUND )
+
 !-------------------------------------------------------------------------------
 !
 !*       5.    TEST ON STRETCHING :
diff --git a/src/MNH/set_grid.f90 b/src/MNH/set_grid.f90
index 1f972fd213b8db63d4b6108dc928da1b78817fd2..25d645e5c6f64e5ef02b131a33ef4db16676cb0e 100644
--- a/src/MNH/set_grid.f90
+++ b/src/MNH/set_grid.f90
@@ -19,6 +19,8 @@ PUBLIC :: SET_GRID
 
 PUBLIC :: INTERP_HORGRID_1DIR_TO_MASSPOINTS, INTERP_HORGRID_TO_MASSPOINTS, INTERP_VERGRID_TO_MASSPOINTS
 
+PUBLIC :: STORE_GRID_BOUNDS, STORE_HORGRID_BOUNDS, STORE_VERGRID_BOUNDS
+
 CONTAINS
 
 !     #####################################################################
@@ -27,6 +29,7 @@ CONTAINS
                            PTSTEP, PSEGLEN,                               &
                            PLONORI, PLATORI, PLON, PLAT,                  &
                            PXHAT, PYHAT, PDXHAT, PDYHAT, PXHATM, PYHATM,  &
+                           PHAT_BOUND, PHATM_BOUND,                       &
                            PMAP, PZS, PZZ, PZHAT, PZHATM, PZTOP, OSLEVE,  &
                            PLEN1, PLEN2, PZSMT, PJ,                       &
                            TPDTMOD, TPDTCUR, KSTOP,                       &
@@ -204,6 +207,8 @@ REAL, DIMENSION(:),     INTENT(OUT) :: PDXHAT    ! horizontal stretching in x
 REAL, DIMENSION(:),     INTENT(OUT) :: PDYHAT    ! horizontal stretching in y
 REAL, DIMENSION(:),     INTENT(OUT) :: PXHATM    ! Position x in the conformal plane or on the cartesian plane at mass points
 REAL, DIMENSION(:),     INTENT(OUT) :: PYHATM    ! Position y in the conformal plane or on the cartesian plane at mass points
+REAL, DIMENSION(:),     INTENT(INOUT) :: PHAT_BOUND  ! Boundaries of global domain in the conformal or cartesian plane
+REAL, DIMENSION(:),     INTENT(INOUT) :: PHATM_BOUND ! Boundaries of global domain in the conformal or cartesian plane at mass pts
 REAL, DIMENSION(:,:),   INTENT(OUT) :: PMAP      ! Map factor
 !
 REAL, DIMENSION(:,:),   INTENT(OUT) :: PZS       ! orography
@@ -346,6 +351,10 @@ CALL IO_Field_read(TPINIFILE,'DTSEG',TDTSEG)
 CALL INTERP_HORGRID_TO_MASSPOINTS( PXHAT, PYHAT, PXHATM, PYHATM )
 CALL INTERP_VERGRID_TO_MASSPOINTS( PZHAT, PZHATM )
 
+! Collect global domain boundaries
+CALL STORE_GRID_BOUNDS( PXHAT,  PYHAT,  PZHAT,  PHAT_BOUND  )
+CALL STORE_GRID_BOUNDS( PXHATM, PYHATM, PZHATM, PHATM_BOUND )
+
 IF (LCARTESIAN) THEN
   CALL SM_GRIDCART(PXHAT,PYHAT,PZHAT,PZS,OSLEVE,PLEN1,PLEN2,PZSMT,PDXHAT,PDYHAT,PZZ,PJ)
 ELSE
@@ -430,75 +439,186 @@ CALL SM_PRINT_TIME(TDTSEG,TLUOUT,YTITLE)
 END SUBROUTINE SET_GRID
 
 
+!-----------------------------------------------------------------
 SUBROUTINE INTERP_HORGRID_1DIR_TO_MASSPOINTS( HDIR, PHAT, PHATM )
-! Interpolate 1 direction of horizontal grid to mass points
+  ! Interpolate 1 direction of horizontal grid to mass points
 
-USE MODD_ARGSLIST_ll, ONLY: LIST1D_ll
+  USE MODD_ARGSLIST_ll, ONLY: LIST1D_ll
 
-USE MODE_ARGSLIST_ll, ONLY: ADD1DFIELD_ll, CLEANLIST1D_ll
-USE MODE_EXCHANGE_ll, ONLY: UPDATE_1DHALO_ll
-USE MODE_MSG
+  USE MODE_ARGSLIST_ll, ONLY: ADD1DFIELD_ll, CLEANLIST1D_ll
+  USE MODE_EXCHANGE_ll, ONLY: UPDATE_1DHALO_ll
+  USE MODE_MSG
 
-IMPLICIT NONE
+  IMPLICIT NONE
 
-CHARACTER(LEN=1),   INTENT(IN)  :: HDIR  ! Direction ('X' or 'Y')
-REAL, DIMENSION(:), INTENT(IN)  :: PHAT  ! Position x or y in the conformal or cartesian plane
-REAL, DIMENSION(:), INTENT(OUT) :: PHATM ! Position x or y in the conformal or cartesian plane at mass points
+  CHARACTER(LEN=1),   INTENT(IN)  :: HDIR  ! Direction ('X' or 'Y')
+  REAL, DIMENSION(:), INTENT(IN)  :: PHAT  ! Position x or y in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(OUT) :: PHATM ! Position x or y in the conformal or cartesian plane at mass points
 
-CHARACTER(LEN=:), ALLOCATABLE :: YNAME
-INTEGER                       :: IINFO_ll ! return code
-TYPE(LIST1D_ll),  POINTER     :: TZLIST   ! pointer for the list of 1D fields to be communicated
+  CHARACTER(LEN=:), ALLOCATABLE :: YNAME
+  INTEGER                       :: IINFO_ll ! return code
+  TYPE(LIST1D_ll),  POINTER     :: TZLIST   ! pointer for the list of 1D fields to be communicated
 
 
-! Interpolate inside subdomain
-PHATM( : UBOUND(PHATM,1)-1 ) = 0.5 * PHAT( : UBOUND(PHAT,1)-1 ) + 0.5 * PHAT( LBOUND(PHAT,1)+1 : UBOUND(PHAT,1) )
-PHATM( UBOUND(PHATM,1)     ) = 1.5 * PHAT( UBOUND(PHAT,1)     ) - 0.5 * PHAT( UBOUND(PHAT,1)-1 )
+  ! Interpolate inside subdomain
+  PHATM( : UBOUND(PHATM,1)-1 ) = 0.5 * PHAT( : UBOUND(PHAT,1)-1 ) + 0.5 * PHAT( LBOUND(PHAT,1)+1 : UBOUND(PHAT,1) )
+  PHATM( UBOUND(PHATM,1)     ) = 1.5 * PHAT( UBOUND(PHAT,1)     ) - 0.5 * PHAT( UBOUND(PHAT,1)-1 )
 
-! Update data between subdomains
-NULLIFY( TZLIST )
+  ! Update data between subdomains
+  NULLIFY( TZLIST )
 
-IF ( HDIR == 'X' ) THEN
-  YNAME = 'XHATM'
-ELSE IF ( HDIR == 'Y' ) THEN
-  YNAME = 'YHATM'
-ELSE
-  CALL PRINT_MSG( NVERB_ERROR, 'GEN', 'INTERP_HORGRID_1DIR_TO_MASSPOINTS', 'invalid direction (valid: X or Y)' )
-END IF
+  IF ( HDIR == 'X' ) THEN
+    YNAME = 'XHATM'
+  ELSE IF ( HDIR == 'Y' ) THEN
+    YNAME = 'YHATM'
+  ELSE
+    CALL PRINT_MSG( NVERB_ERROR, 'GEN', 'INTERP_HORGRID_1DIR_TO_MASSPOINTS', 'invalid direction (valid: X or Y)' )
+  END IF
 
-CALL ADD1DFIELD_ll( HDIR, TZLIST, PHATM, YNAME )
-CALL UPDATE_1DHALO_ll( TZLIST, IINFO_ll )
-CALL CLEANLIST1D_ll( TZLIST )
+  CALL ADD1DFIELD_ll( HDIR, TZLIST, PHATM, YNAME )
+  CALL UPDATE_1DHALO_ll( TZLIST, IINFO_ll )
+  CALL CLEANLIST1D_ll( TZLIST )
 
 END SUBROUTINE INTERP_HORGRID_1DIR_TO_MASSPOINTS
+!-----------------------------------------------------------------
 
 
+!-----------------------------------------------------------------
 SUBROUTINE INTERP_HORGRID_TO_MASSPOINTS( PXHAT, PYHAT, PXHATM, PYHATM )
 
-IMPLICIT NONE
+  IMPLICIT NONE
 
-REAL, DIMENSION(:), INTENT(IN)  :: PXHAT  ! Position x in the conformal or cartesian plane
-REAL, DIMENSION(:), INTENT(IN)  :: PYHAT  ! Position y in the conformal or cartesian plane
-REAL, DIMENSION(:), INTENT(OUT) :: PXHATM ! Position x in the conformal or cartesian plane at mass points
-REAL, DIMENSION(:), INTENT(OUT) :: PYHATM ! Position y in the conformal or cartesian plane at mass points
+  REAL, DIMENSION(:), INTENT(IN)  :: PXHAT  ! Position x in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(IN)  :: PYHAT  ! Position y in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(OUT) :: PXHATM ! Position x in the conformal or cartesian plane at mass points
+  REAL, DIMENSION(:), INTENT(OUT) :: PYHATM ! Position y in the conformal or cartesian plane at mass points
 
-CALL INTERP_HORGRID_1DIR_TO_MASSPOINTS( 'X', PXHAT, PXHATM )
-CALL INTERP_HORGRID_1DIR_TO_MASSPOINTS( 'Y', PYHAT, PYHATM )
+  CALL INTERP_HORGRID_1DIR_TO_MASSPOINTS( 'X', PXHAT, PXHATM )
+  CALL INTERP_HORGRID_1DIR_TO_MASSPOINTS( 'Y', PYHAT, PYHATM )
 
 END SUBROUTINE INTERP_HORGRID_TO_MASSPOINTS
+!-----------------------------------------------------------------
 
 
+!-----------------------------------------------------------------
 PURE SUBROUTINE INTERP_VERGRID_TO_MASSPOINTS( PZHAT, PZHATM )
-! Interpolate vertical grid to mass points
+  ! Interpolate vertical grid to mass points
 
-IMPLICIT NONE
+  IMPLICIT NONE
 
-REAL, DIMENSION(:), INTENT(IN)  :: PZHAT  ! Position z in the conformal or cartesian plane
-REAL, DIMENSION(:), INTENT(OUT) :: PZHATM ! Position z in the conformal or cartesian plane at mass points
+  REAL, DIMENSION(:), INTENT(IN)  :: PZHAT  ! Position z in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(OUT) :: PZHATM ! Position z in the conformal or cartesian plane at mass points
 
-PZHATM( : UBOUND(PZHATM,1)-1 ) = 0.5 * PZHAT( : UBOUND(PZHAT,1)-1 ) + 0.5 * PZHAT( LBOUND(PZHAT,1)+1 : UBOUND(PZHAT,1) )
-PZHATM( UBOUND(PZHATM,1)     ) = 1.5 * PZHAT( UBOUND(PZHAT,1)     ) - 0.5 * PZHAT( UBOUND(PZHAT,1)-1 )
+  PZHATM( : UBOUND(PZHATM,1)-1 ) = 0.5 * PZHAT( : UBOUND(PZHAT,1)-1 ) + 0.5 * PZHAT( LBOUND(PZHAT,1)+1 : UBOUND(PZHAT,1) )
+  PZHATM( UBOUND(PZHATM,1)     ) = 1.5 * PZHAT( UBOUND(PZHAT,1)     ) - 0.5 * PZHAT( UBOUND(PZHAT,1)-1 )
 
 END SUBROUTINE INTERP_VERGRID_TO_MASSPOINTS
+!-----------------------------------------------------------------
+
+
+!-----------------------------------------------------------------
+SUBROUTINE STORE_GRID_1DIR_BOUNDS( HDIR, PHAT, PHAT_BOUND )
+
+  USE MODD_GRID_n
+  USE MODD_PARAMETERS,     ONLY: JPHEXT, JPVEXT
+
+  USE MODE_ALLOCBUFFER_ll, ONLY: ALLOCBUFFER_ll
+  USE MODE_GATHER_ll,      ONLY: GATHERALL_FIELD_ll
+  USE MODE_MSG
+
+  IMPLICIT NONE
+
+  CHARACTER(LEN=1),                 INTENT(IN)    :: HDIR       ! Direction ('X', 'Y' or 'Z')
+  REAL, DIMENSION(:), TARGET,       INTENT(IN)    :: PHAT       ! Position x, y or z in the conformal or cartesian plane
+  REAL, DIMENSION(NHAT_BOUND_SIZE), INTENT(INOUT) :: PHAT_BOUND ! Boundaries of global domain in the conformal or cartesian plane
+
+  INTEGER                     :: IERR
+  LOGICAL                     :: GALLOC
+  REAL, DIMENSION(:), POINTER :: ZHAT_GLOB
+
+  ZHAT_GLOB => NULL()
+
+  SELECT CASE (HDIR)
+    CASE ( 'X' )
+      CALL ALLOCBUFFER_ll( ZHAT_GLOB, PHAT, 'XX', GALLOC )
+      CALL GATHERALL_FIELD_ll( 'XX', PHAT, ZHAT_GLOB, IERR )
+      PHAT_BOUND(NEXTE_XMIN) = ZHAT_GLOB( 1 )
+      PHAT_BOUND(NEXTE_XMAX) = ZHAT_GLOB( UBOUND( ZHAT_GLOB, 1 ) )
+      PHAT_BOUND(NPHYS_XMIN) = ZHAT_GLOB( JPHEXT + 1 )
+      PHAT_BOUND(NPHYS_XMAX) = ZHAT_GLOB( UBOUND( ZHAT_GLOB, 1 ) - JPHEXT )
+
+    CASE ( 'Y' )
+      CALL ALLOCBUFFER_ll( ZHAT_GLOB, PHAT, 'YY', GALLOC )
+      CALL GATHERALL_FIELD_ll( 'YY', PHAT, ZHAT_GLOB, IERR )
+      PHAT_BOUND(NEXTE_YMIN) = ZHAT_GLOB( 1 )
+      PHAT_BOUND(NEXTE_YMAX) = ZHAT_GLOB( UBOUND( ZHAT_GLOB, 1 ) )
+      PHAT_BOUND(NPHYS_YMIN) = ZHAT_GLOB( JPHEXT + 1 )
+      PHAT_BOUND(NPHYS_YMAX) = ZHAT_GLOB( UBOUND( ZHAT_GLOB, 1 ) - JPHEXT )
+
+    CASE ( 'Z' )
+      ZHAT_GLOB => PHAT
+      PHAT_BOUND(NEXTE_ZMIN) = ZHAT_GLOB( 1 )
+      PHAT_BOUND(NEXTE_ZMAX) = ZHAT_GLOB( UBOUND( ZHAT_GLOB, 1 ) )
+      PHAT_BOUND(NPHYS_ZMIN) = ZHAT_GLOB( JPVEXT + 1 )
+      PHAT_BOUND(NPHYS_ZMAX) = ZHAT_GLOB( UBOUND( ZHAT_GLOB, 1 ) - JPVEXT )
+
+    CASE DEFAULT
+      CALL PRINT_MSG( NVERB_ERROR, 'GEN', 'STORE_GRID_1DIR_BOUNDS', 'invalid direction (valid: X, Y or Z)' )
+
+  END SELECT
+
+  IF ( GALLOC ) DEALLOCATE( ZHAT_GLOB )
+
+END SUBROUTINE STORE_GRID_1DIR_BOUNDS
+!-----------------------------------------------------------------
+
+
+!-----------------------------------------------------------------
+SUBROUTINE STORE_GRID_BOUNDS( PXHAT, PYHAT, PZHAT, PHAT_BOUND )
+
+  IMPLICIT NONE
+
+  REAL, DIMENSION(:), INTENT(IN)    :: PXHAT  ! Position x in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(IN)    :: PYHAT  ! Position y in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(IN)    :: PZHAT  ! Position y in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(INOUT) :: PHAT_BOUND ! Boundaries of global domain in the conformal or cartesian plane
+
+  CALL STORE_GRID_1DIR_BOUNDS( 'X', PXHAT, PHAT_BOUND )
+  CALL STORE_GRID_1DIR_BOUNDS( 'Y', PYHAT, PHAT_BOUND )
+  CALL STORE_GRID_1DIR_BOUNDS( 'Z', PZHAT, PHAT_BOUND )
+
+END SUBROUTINE STORE_GRID_BOUNDS
+!-----------------------------------------------------------------
+
+
+!-----------------------------------------------------------------
+SUBROUTINE STORE_HORGRID_BOUNDS( PXHAT, PYHAT, PHAT_BOUND )
+
+  IMPLICIT NONE
+
+  REAL, DIMENSION(:), INTENT(IN)    :: PXHAT  ! Position x in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(IN)    :: PYHAT  ! Position y in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(INOUT) :: PHAT_BOUND ! Boundaries of global domain in the conformal or cartesian plane
+
+  CALL STORE_GRID_1DIR_BOUNDS( 'X', PXHAT, PHAT_BOUND )
+  CALL STORE_GRID_1DIR_BOUNDS( 'Y', PYHAT, PHAT_BOUND )
+
+END SUBROUTINE STORE_HORGRID_BOUNDS
+!-----------------------------------------------------------------
+
+
+!-----------------------------------------------------------------
+SUBROUTINE STORE_VERGRID_BOUNDS( PZHAT, PHAT_BOUND )
+
+  IMPLICIT NONE
+
+  REAL, DIMENSION(:), INTENT(IN)    :: PZHAT  ! Position y in the conformal or cartesian plane
+  REAL, DIMENSION(:), INTENT(INOUT) :: PHAT_BOUND ! Boundaries of global domain in the conformal or cartesian plane
+
+  CALL STORE_GRID_1DIR_BOUNDS( 'Z', PZHAT, PHAT_BOUND )
+
+END SUBROUTINE STORE_VERGRID_BOUNDS
+!-----------------------------------------------------------------
 
 
 END MODULE MODE_SET_GRID
diff --git a/src/MNH/spawn_grid2.f90 b/src/MNH/spawn_grid2.f90
index 6909c693e6822981f98800333cc3df91f07dc70e..564efb71650d842c3e266ec0d9e8b542275bd0f8 100644
--- a/src/MNH/spawn_grid2.f90
+++ b/src/MNH/spawn_grid2.f90
@@ -11,6 +11,7 @@ INTERFACE
 !
      SUBROUTINE SPAWN_GRID2( KXOR, KYOR, KXEND, KYEND, KDXRATIO, KDYRATIO,                &
                              PLONOR, PLATOR, PXHAT, PYHAT, PZHAT, PXHATM, PYHATM, PZHATM, &
+                             PHAT_BOUND, PHATM_BOUND,                                     &
                              PZTOP, OSLEVE, PLEN1, PLEN2,                                 &
                              PZS, PZSMT, PZS_LS, PZSMT_LS,                                &
                              TPDTMOD, TPDTCUR                                             )
@@ -29,6 +30,8 @@ REAL, DIMENSION(:),   INTENT(OUT) :: PXHAT,PYHAT,PZHAT ! positions x,y,z in the
                                      ! conformal plane or on the cartesian plane
 REAL, DIMENSION(:),   INTENT(OUT) :: PXHATM, PYHATM, PZHATM ! positions x,y in the
                                      ! conformal plane or on the cartesian plane at mass points
+REAL, DIMENSION(:),   INTENT(INOUT) :: PHAT_BOUND  ! Boundaries of global domain in the conformal or cartesian plane
+REAL, DIMENSION(:),   INTENT(INOUT) :: PHATM_BOUND ! Boundaries of global domain in the conformal or cartesian plane at mass points
 REAL,                 INTENT(OUT)   :: PZTOP             ! model top (m)
 LOGICAL,              INTENT(OUT)   :: OSLEVE            ! flag for SLEVE coordinate
 REAL,                 INTENT(OUT)   :: PLEN1             ! Decay scale for smooth topography
@@ -51,6 +54,7 @@ END MODULE MODI_SPAWN_GRID2
 !    ######################################################################################
      SUBROUTINE SPAWN_GRID2( KXOR, KYOR, KXEND, KYEND, KDXRATIO, KDYRATIO,                &
                              PLONOR, PLATOR, PXHAT, PYHAT, PZHAT, PXHATM, PYHATM, PZHATM, &
+                             PHAT_BOUND, PHATM_BOUND,                                     &
                              PZTOP, OSLEVE, PLEN1, PLEN2,                                 &
                              PZS, PZSMT, PZS_LS, PZSMT_LS,                                &
                              TPDTMOD, TPDTCUR                                             )
@@ -167,7 +171,7 @@ USE MODD_BIKHARDT_n
 USE MODD_VAR_ll
 use mode_bikhardt
 USE MODE_ll
-USE MODE_SET_GRID,   only: INTERP_HORGRID_TO_MASSPOINTS
+USE MODE_SET_GRID,   only: INTERP_HORGRID_TO_MASSPOINTS, STORE_GRID_BOUNDS
 USE MODE_TIME
 USE MODE_GRIDPROJ
 !
@@ -192,6 +196,8 @@ REAL, DIMENSION(:),   INTENT(OUT) :: PXHAT,PYHAT,PZHAT ! positions x,y,z in the
                                      ! conformal plane or on the cartesian plane
 REAL, DIMENSION(:),   INTENT(OUT) :: PXHATM, PYHATM, PZHATM ! positions x,y in the
                                      ! conformal plane or on the cartesian plane at mass points
+REAL, DIMENSION(:),   INTENT(INOUT) :: PHAT_BOUND  ! Boundaries of global domain in the conformal or cartesian plane
+REAL, DIMENSION(:),   INTENT(INOUT) :: PHATM_BOUND ! Boundaries of global domain in the conformal or cartesian plane at mass points
 REAL,                 INTENT(OUT)   :: PZTOP             ! model top (m)
 LOGICAL,              INTENT(OUT)   :: OSLEVE            ! flag for SLEVE coordinate
 REAL,                 INTENT(OUT)   :: PLEN1             ! Decay scale for smooth topography
@@ -457,6 +463,10 @@ PLEN2    = XLEN21
   ! Interpolations of positions to mass points
   CALL INTERP_HORGRID_TO_MASSPOINTS( PXHAT, PYHAT, PXHATM, PYHATM )
 
+  ! Collect global domain boundaries
+  CALL STORE_GRID_BOUNDS( PXHAT,  PYHAT,  PZHAT,  PHAT_BOUND  )
+  CALL STORE_GRID_BOUNDS( PXHATM, PYHATM, PZHATM, PHATM_BOUND )
+
 !!$=======
 !!$  IXSIZE1=SIZE(XXHAT1)
 !!$  ALLOCATE(ZXHAT_EXTENDED(IXSIZE1+1))
diff --git a/src/MNH/spawn_model2.f90 b/src/MNH/spawn_model2.f90
index 8837733b60b27b1b23673ddffbd5e506afbe8ce6..a824c074b9c2042cada031c2d87887068071d243 100644
--- a/src/MNH/spawn_model2.f90
+++ b/src/MNH/spawn_model2.f90
@@ -736,6 +736,7 @@ END IF
 !
 ALLOCATE(XXHAT(IIU),XYHAT(IJU),XZHAT(IKU))
 ALLOCATE(XXHATM(IIU),XYHATM(IJU),XZHATM(IKU))
+ALLOCATE( XHAT_BOUND(NHAT_BOUND_SIZE), XHATM_BOUND(NHAT_BOUND_SIZE) )
 ALLOCATE(XZTOP)
 ALLOCATE(XMAP(IIU,IJU))
 ALLOCATE(XLAT(IIU,IJU))
@@ -1061,6 +1062,7 @@ ENDIF
 XZS=0.
 CALL SPAWN_GRID2( NXOR, NYOR, NXEND, NYEND, NDXRATIO, NDYRATIO,                  &
                   XLONORI, XLATORI, XXHAT, XYHAT, XZHAT, XXHATM, XYHATM, XZHATM, &
+                  XHAT_BOUND, XHATM_BOUND,                                       &
                   XZTOP, LSLEVE, XLEN1, XLEN2,                                   &
                   XZS, XZSMT, ZZS_LS, ZZSMT_LS, TDTMOD, TDTCUR                   )
 !