From cda9ec857ff4189c751d1606a090ea8519d1d0a1 Mon Sep 17 00:00:00 2001
From: Quentin Rodier <quentin.rodier@meteo.fr>
Date: Mon, 2 Dec 2019 15:44:44 +0100
Subject: [PATCH] Robert S. 02/12/2019: Cartesian coordinates allowed for
stations and passive polutants (paspol). Use of LCARTESIAN instead of LATSTAT
flag for multi-proc
---
src/MNH/ini_surfstationn.f90 | 23 +++++++--
src/MNH/paspol.f90 | 7 ++-
src/MNH/read_grid_time_mesonh_case.f90 | 38 ++++++++-------
src/MNH/stationn.f90 | 64 +++++---------------------
src/MNH/write_stationn.f90 | 25 ++++++++++
5 files changed, 81 insertions(+), 76 deletions(-)
diff --git a/src/MNH/ini_surfstationn.f90 b/src/MNH/ini_surfstationn.f90
index e90f679a3..23f7526ec 100644
--- a/src/MNH/ini_surfstationn.f90
+++ b/src/MNH/ini_surfstationn.f90
@@ -66,7 +66,7 @@ END MODULE MODI_INI_SURFSTATION_n
!! P. Tulet 15/01/2002
!! A. Lemonsu 19/11/2002
!! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
-!!
+!! R. Schoetter : 11/2019 : work for cartesian coordinates + parallel.
!! --------------------------------------------------------------------------
!
!* 0. DECLARATIONS
@@ -84,9 +84,13 @@ USE MODD_TYPE_DATE
USE MODE_GRIDPROJ
USE MODE_IO_ll
USE MODE_ll
+USE MODE_GATHER_ll
USE MODE_MSG
!
USE MODI_INI_STATION_N
+USE MODD_VAR_ll, ONLY: IP
+USE MODD_SHADOWS_n
+USE MODD_DIM_n
!
IMPLICIT NONE
!
@@ -109,6 +113,7 @@ REAL, INTENT(IN) :: PLONOR ! longitude of origine point
!
INTEGER :: ISTORE ! number of storage instants
INTEGER :: ILUOUT ! logical unit
+INTEGER :: IIU_ll,IJU_ll,IRESP
!
!----------------------------------------------------------------------------
ILUOUT = TLUOUT%NLU
@@ -240,7 +245,6 @@ INTEGER :: JII !
INTEGER :: IIU, IJU !
!
IF ( ALL(TSTATION%LAT(:)/=XUNDEF) .AND. ALL(TSTATION%LON(:)/=XUNDEF) ) THEN
- LSTATLAT = .TRUE.
DO JII=1,NUMBSTAT
CALL GET_DIM_EXT_ll ('B',IIU,IJU)
CALL SM_XYHAT(PLATOR,PLONOR, &
@@ -248,11 +252,20 @@ IF ( ALL(TSTATION%LAT(:)/=XUNDEF) .AND. ALL(TSTATION%LON(:)/=XUNDEF) ) THEN
TSTATION%X(JII), TSTATION%Y(JII) )
ENDDO
ELSE
- LSTATLAT = .FALSE.
DO JII=1,NUMBSTAT
- TSTATION%X(JII) = XXHAT(TSTATION%I(JII))
- TSTATION%Y(JII) = XYHAT(TSTATION%I(JII))
CALL GET_DIM_EXT_ll ('B',IIU,IJU)
+ IIU_ll=NIMAX_ll + 2 * JPHEXT
+ IJU_ll=NJMAX_ll + 2 * JPHEXT
+ ALLOCATE(XXHAT_ll (IIU_ll))
+ ALLOCATE(XYHAT_ll (IJU_ll))
+ !
+ CALL GATHERALL_FIELD_ll('XX',XXHAT,XXHAT_ll,IRESP)
+ CALL GATHERALL_FIELD_ll('YY',XYHAT,XYHAT_ll,IRESP)
+ TSTATION%X(JII) = XXHAT_ll(TSTATION%I(JII))
+ TSTATION%Y(JII) = XYHAT_ll(TSTATION%J(JII))
+ IF (LCARTESIAN) THEN
+ XRPK = -1
+ ENDIF
CALL SM_LATLON(PLATOR,PLONOR, &
TSTATION%X(JII), TSTATION%Y(JII), &
TSTATION%LAT(JII), TSTATION%LON(JII) )
diff --git a/src/MNH/paspol.f90 b/src/MNH/paspol.f90
index 99af6a559..11e660d34 100644
--- a/src/MNH/paspol.f90
+++ b/src/MNH/paspol.f90
@@ -216,7 +216,12 @@ IF (GPPFIRSTCALL) THEN
! puis les indices fractionnaires (ZSRCI,ZSRCJ) et entiers
! (IPIGI,IPIGJ) du point de rejet dans le domaine de travail global.
!
- CALL SM_XYHAT(XLATORI,XLONORI,XPPLAT(JSV),XPPLON(JSV),ZSRCX,ZSRCY)
+ IF (LCARTESIAN) THEN !En cartésien écriture dans la namelist des coordonées X,Y et non LAT,LON
+ ZSRCX = XPPLAT(JSV)
+ ZSRCY = XPPLON(JSV)
+ ELSE
+ CALL SM_XYHAT(XLATORI,XLONORI,XPPLAT(JSV),XPPLON(JSV),ZSRCX,ZSRCY)
+ END IF
II=MAX(MIN(COUNT(XXHAT(:)<ZSRCX),IIU-1),1)
IJ=MAX(MIN(COUNT(XYHAT(:)<ZSRCY),IJU-1),1)
ZSRCI=(ZSRCX-XXHAT(II))/(XXHAT(II+1)-XXHAT(II))+FLOAT(II)
diff --git a/src/MNH/read_grid_time_mesonh_case.f90 b/src/MNH/read_grid_time_mesonh_case.f90
index 9c2c6669e..6b27c684e 100644
--- a/src/MNH/read_grid_time_mesonh_case.f90
+++ b/src/MNH/read_grid_time_mesonh_case.f90
@@ -174,23 +174,27 @@ CALL IO_READ_FIELD(TZFMFILE,'RPK', ZRPK_LS)
CALL IO_READ_FIELD(TZFMFILE,'LAT0',ZLAT0_LS)
CALL IO_READ_FIELD(TZFMFILE,'BETA',ZBETA_LS)
!
-IF ( (ABS(ZLAT0_LS-XLAT0)>ZEPS*MAX(1.,ABS(XLAT0))) &
- .OR. (ABS(ZLON0_LS-XLON0)>ZEPS*MAX(1.,ABS(XLON0))) &
- .OR. (ABS(ABS(ZRPK_LS)-ABS(XRPK))>ZEPS*MAX(1.,ABS(XRPK))) &
- .OR. (ABS(ZBETA_LS-XBETA)>ZEPS*MAX(1.,ABS(XBETA))) ) THEN
-!
- WRITE(ILUOUT0,FMT=*) ' '
- WRITE(ILUOUT0,FMT=*) '***************************************************************'
- WRITE(ILUOUT0,FMT=*) 'Projection are different between MESONH input file and PGD file'
- WRITE(ILUOUT0,FMT=*) 'You must recompute a PGD file with PREP_PGD,'
- WRITE(ILUOUT0,FMT=*) 'using the input MESONH file to define its domain.'
- WRITE(ILUOUT0,FMT=*) '***************************************************************'
- WRITE(ILUOUT0,FMT=*) ' '
- WRITE(ILUOUT0,FMT=*) ' input file physiographic data'
- WRITE(ILUOUT0,1) 'LAT0 ',ZLAT0_LS, ' ',XLAT0
- WRITE(ILUOUT0,1) 'LON0 ',ZLON0_LS, ' ',XLON0
- WRITE(ILUOUT0,1) 'RPK ',ZRPK_LS, ' ',XRPK
- WRITE(ILUOUT0,1) 'BETA ',ZBETA_LS, ' ',XBETA
+IF(.NOT.LCARTESIAN) THEN
+ !
+ IF ( (ABS(ZLAT0_LS-XLAT0)>ZEPS*MAX(1.,ABS(XLAT0))) &
+ .OR. (ABS(ZLON0_LS-XLON0)>ZEPS*MAX(1.,ABS(XLON0))) &
+ .OR. (ABS(ABS(ZRPK_LS)-ABS(XRPK))>ZEPS*MAX(1.,ABS(XRPK))) &
+ .OR. (ABS(ZBETA_LS-XBETA)>ZEPS*MAX(1.,ABS(XBETA))) ) THEN
+ !
+ WRITE(ILUOUT0,FMT=*) ' '
+ WRITE(ILUOUT0,FMT=*) '***************************************************************'
+ WRITE(ILUOUT0,FMT=*) 'Projection are different between MESONH input file and PGD file'
+ WRITE(ILUOUT0,FMT=*) 'You must recompute a PGD file with PREP_PGD,'
+ WRITE(ILUOUT0,FMT=*) 'using the input MESONH file to define its domain.'
+ WRITE(ILUOUT0,FMT=*) '***************************************************************'
+ WRITE(ILUOUT0,FMT=*) ' '
+ WRITE(ILUOUT0,FMT=*) ' input file physiographic data'
+ WRITE(ILUOUT0,1) 'LAT0 ',ZLAT0_LS, ' ',XLAT0
+ WRITE(ILUOUT0,1) 'LON0 ',ZLON0_LS, ' ',XLON0
+ WRITE(ILUOUT0,1) 'RPK ',ZRPK_LS, ' ',XRPK
+ WRITE(ILUOUT0,1) 'BETA ',ZBETA_LS, ' ',XBETA
+ END IF
+ !
END IF
!
!* 2.2 Horizontal grid:
diff --git a/src/MNH/stationn.f90 b/src/MNH/stationn.f90
index f1d20df7e..0e1c54e49 100644
--- a/src/MNH/stationn.f90
+++ b/src/MNH/stationn.f90
@@ -86,6 +86,7 @@ END MODULE MODI_STATION_n
!! C.Lac 04/2013 : Add I/J positioning
!! P.Wautelet 28/03/2018 : Replace TEMPORAL_DIST by DATETIME_DISTANCE
!! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
+!! R.Schoetter 11/2019 : use LCARTESIAN instead of LSTATLAT for multiproc in cartesian
!!
!! --------------------------------------------------------------------------
!
@@ -292,7 +293,6 @@ IF (GSTATFIRSTCALL) THEN
!
!* 4.4 Computations only on correct processor
! --------------------------------------
- IF ( LSTATLAT ) THEN
ZXCOEF(I) = 0.
ZYCOEF(I) = 0.
ZUCOEF(I) = 0.
@@ -331,7 +331,6 @@ IF (GSTATFIRSTCALL) THEN
!
END IF
- END IF
ENDDO
END IF
!----------------------------------------------------------------------------
@@ -360,18 +359,16 @@ IF (GSTORE) THEN
ENDIF
END IF
!
- ZGAM = (XRPK * (TSTATION%LON(I) - XLON0) - XBETA)*(XPI/180.)
- IF ( LSTATLAT ) THEN
- ZU_STAT = STATION_INTERP_2D_U(PU(:,:,J))
- ZV_STAT = STATION_INTERP_2D_V(PV(:,:,J))
+ IF (LCARTESIAN) THEN
+ TSTATION%ZON (IN,I) = STATION_INTERP_2D_U(PU(:,:,J))
+ TSTATION%MER (IN,I) = STATION_INTERP_2D_V(PV(:,:,J))
ELSE
- ZU_STAT = PU(TSTATION%I(I),TSTATION%J(I),J)
- ZV_STAT = PV(TSTATION%I(I),TSTATION%J(I),J)
- END IF
- !
- TSTATION%ZON (IN,I) = ZU_STAT * COS(ZGAM) + ZV_STAT * SIN(ZGAM)
- TSTATION%MER (IN,I) = - ZU_STAT * SIN(ZGAM) + ZV_STAT * COS(ZGAM)
- IF ( LSTATLAT ) THEN
+ ZU_STAT = STATION_INTERP_2D_U(PU(:,:,J))
+ ZV_STAT = STATION_INTERP_2D_V(PV(:,:,J))
+ ZGAM = (XRPK * (TSTATION%LON(I) - XLON0) - XBETA)*(XPI/180.)
+ TSTATION%ZON (IN,I) = ZU_STAT * COS(ZGAM) + ZV_STAT * SIN(ZGAM)
+ TSTATION%MER (IN,I) = - ZU_STAT * SIN(ZGAM) + ZV_STAT * COS(ZGAM)
+ ENDIF
TSTATION%W (IN,I) = STATION_INTERP_2D(PW(:,:,J))
TSTATION%TH (IN,I) = STATION_INTERP_2D(PTH(:,:,J))
TSTATION%P (IN,I) = STATION_INTERP_2D(PP(:,:,J))
@@ -410,46 +407,7 @@ IF (GSTORE) THEN
ENDIF
TSTATION%SFCO2 (IN,I) = STATION_INTERP_2D(XCURRENT_SFCO2 )
ENDIF
- ELSE
- TSTATION%W (IN,I) = PW(TSTATION%I(I),TSTATION%J(I),J)
- TSTATION%TH (IN,I) = PTH(TSTATION%I(I),TSTATION%J(I),J)
- TSTATION%P (IN,I) = PP(TSTATION%I(I),TSTATION%J(I),J)
- !
- DO JSV=1,SIZE(PR,4)
- TSTATION%R (IN,I,JSV) = PR(TSTATION%I(I),TSTATION%J(I),J,JSV)
- END DO
- !
- DO JSV=1,SIZE(PSV,4)
- TSTATION%SV (IN,I,JSV) = PSV(TSTATION%I(I),TSTATION%J(I),J,JSV)
- END DO
- !
- IF (SIZE(PTKE)>0) TSTATION%TKE (IN,I) = PTKE(TSTATION%I(I),TSTATION%J(I),J)
- IF (SIZE(PTS) >0) TSTATION%TSRAD(IN,I) = PTS(TSTATION%I(I),TSTATION%J(I))
- TSTATION%ZS(I) = PZ(TSTATION%I(I),TSTATION%J(I),1+JPVEXT)
- !
- IF (LDIAG_IN_RUN) THEN
- TSTATION%ZON10M(IN,I) = XCURRENT_ZON10M(TSTATION%I(I),TSTATION%J(I))
- TSTATION%MER10M(IN,I) = XCURRENT_MER10M(TSTATION%I(I),TSTATION%J(I))
- TSTATION%T2M (IN,I) = XCURRENT_T2M(TSTATION%I(I),TSTATION%J(I))
- TSTATION%Q2M (IN,I) = XCURRENT_Q2M(TSTATION%I(I),TSTATION%J(I))
- TSTATION%HU2M (IN,I) = XCURRENT_HU2M(TSTATION%I(I),TSTATION%J(I))
- TSTATION%RN (IN,I) = XCURRENT_RN(TSTATION%I(I),TSTATION%J(I))
- TSTATION%H (IN,I) = XCURRENT_H(TSTATION%I(I),TSTATION%J(I))
- TSTATION%LE (IN,I) = XCURRENT_LE(TSTATION%I(I),TSTATION%J(I))
- TSTATION%LEI (IN,I) = XCURRENT_LEI(TSTATION%I(I),TSTATION%J(I))
- TSTATION%GFLUX (IN,I) = XCURRENT_GFLUX(TSTATION%I(I),TSTATION%J(I))
- IF (CRAD /= 'NONE') THEN
- TSTATION%SWD (IN,I) = XCURRENT_SWD(TSTATION%I(I),TSTATION%J(I))
- TSTATION%SWU (IN,I) = XCURRENT_SWU(TSTATION%I(I),TSTATION%J(I))
- TSTATION%LWD (IN,I) = XCURRENT_LWD(TSTATION%I(I),TSTATION%J(I))
- TSTATION%LWU (IN,I) = XCURRENT_LWU(TSTATION%I(I),TSTATION%J(I))
- TSTATION%SWDIR (IN,I) = XCURRENT_SWDIR(TSTATION%I(I),TSTATION%J(I))
- TSTATION%SWDIFF(IN,I) = XCURRENT_SWDIFF(TSTATION%I(I),TSTATION%J(I))
- TSTATION%DSTAOD(IN,I) = XCURRENT_DSTAOD(TSTATION%I(I),TSTATION%J(I))
- ENDIF
- TSTATION%SFCO2 (IN,I) = XCURRENT_SFCO2(TSTATION%I(I),TSTATION%J(I))
- ENDIF
- ENDIF
+
!
END IF
!
diff --git a/src/MNH/write_stationn.f90 b/src/MNH/write_stationn.f90
index 48d810835..22d45f61b 100644
--- a/src/MNH/write_stationn.f90
+++ b/src/MNH/write_stationn.f90
@@ -89,6 +89,8 @@ USE MODE_SALT_PSD
USE MODE_AERO_PSD
!
USE MODI_WRITE_DIACHRO
+USE MODD_PASPOL
+USE MODD_CONF
!
IMPLICIT NONE
!
@@ -200,6 +202,20 @@ YUNIT (JPROC) = 'degree'
YCOMMENT (JPROC) = 'Latitude'
ZWORK6 (1,1,1,:,1,JPROC) = TSTATION%LAT(II)
!
+IF (LCARTESIAN) THEN
+ JPROC = JPROC + 1
+ YTITLE (JPROC) = 'X'
+ YUNIT (JPROC) = 'm'
+ YCOMMENT (JPROC) = 'X Pos'
+ ZWORK6 (1,1,1,:,1,JPROC) = TSTATION%X(II)
+ !
+ JPROC = JPROC + 1
+ YTITLE (JPROC) = 'Y'
+ YUNIT (JPROC) = 'm'
+ YCOMMENT (JPROC) = 'Y Pos'
+ ZWORK6 (1,1,1,:,1,JPROC) = TSTATION%Y(II)
+ENDIF
+!
JPROC = JPROC + 1
YTITLE (JPROC) = 'ZON_WIND'
YUNIT (JPROC) = 'm s-1'
@@ -367,6 +383,15 @@ IF (SIZE(TSTATION%TKE,1)>0) THEN
END IF
!
!
+IF (LPASPOL) THEN
+ JSV=1
+ JPROC = JPROC+1
+ WRITE (YTITLE(JPROC),FMT='(A2,I3.3)') 'Sv',JSV
+ YUNIT (JPROC) = 'kg kg-1'
+ YCOMMENT (JPROC) = ' '
+ ZWORK6 (1,1,1,:,1,JPROC) = TSTATION%SV(:,II,JSV)
+ENDIF
+!
IF (SIZE(TSTATION%SV,3)>=1) THEN
! User scalar variables
DO JSV = 1,NSV_USER
--
GitLab