From 46feb413d36d8fec4b95078d66317477260903a6 Mon Sep 17 00:00:00 2001
From: Jean-Pierre CHABOUREAU <jean-pierre.chaboureau@aero.obs-mip.fr>
Date: Fri, 21 Oct 2022 17:36:56 +0200
Subject: [PATCH] Jean-Pierre 21/10/2022: PREP_REAL_CASE fix for GFS/ERA5
analysis on pressure levels
---
src/MNH/read_all_data_grib_case.f90 | 17 +++++++++--------
src/MNH/ver_prep_gribex_case.f90 | 29 ++++++++++++++---------------
2 files changed, 23 insertions(+), 23 deletions(-)
diff --git a/src/MNH/read_all_data_grib_case.f90 b/src/MNH/read_all_data_grib_case.f90
index 04a8ff5e8..ec83d1da9 100644
--- a/src/MNH/read_all_data_grib_case.f90
+++ b/src/MNH/read_all_data_grib_case.f90
@@ -136,7 +136,8 @@ END MODULE MODI_READ_ALL_DATA_GRIB_CASE
! Q. Rodier 21/04/2020: correction GFS u and v wind component written in the right vertical order
! Q. Rodier 02/09/2020: Read and interpol geopotential height for interpolation on isobaric surface Grid of NCEP
! P. Wautelet 09/03/2021: move some chemistry initializations to ini_nsv
-!JP Chaboureau 02/08/2021: add ERA5 reanlysis in pressure levels
+!JP Chaboureau 02/08/2021: add ERA5 reanalysis in pressure levels
+!JP Chaboureau 18/10/2022: correction on vertical level for GFS and ERA5 reanalyses in pressure levels
!-------------------------------------------------------------------------------
!
!* 0. DECLARATIONS
@@ -925,6 +926,9 @@ ALLOCATE (ZEXNM_G(INI,INLEVEL))
ZEXNM_G(:,1:INLEVEL-1) = (ZEXNF_G(:,1:INLEVEL-1)-ZEXNF_G(:,2:INLEVEL)) / &
(LOG(ZEXNF_G(:,1:INLEVEL-1))-LOG(ZEXNF_G(:,2:INLEVEL)))
ZEXNM_G(:,INLEVEL) = (ZPF_G(:,INLEVEL)/2./XP00)**(XRD/XCPD)
+!
+IF (IMODEL==10.OR.IMODEL==11) ZEXNM_G(:,:)=ZEXNF_G(:,:) ! for GFS and ERA5 on pressure levels
+!
DEALLOCATE (ZEXNF_G)
DEALLOCATE (ZPF_G)
!
@@ -1010,8 +1014,6 @@ ALLOCATE (ZRV_G(INI))
ALLOCATE (ZOUT(INO))
IF (IMODEL/=10) THEN ! others than NCEP
DO JLOOP1=1, INLEVEL
- !WRITE (ILUOUT0,*) 'JLOOP1=',JLOOP1,MINVAL(ZPM_G(:,JLOOP1)),MINVAL(ZT_G(:,JLOOP1)),MINVAL(ZQ_G(:,JLOOP1))
- !WRITE (ILUOUT0,*) ' ',MAXVAL(ZPM_G(:,JLOOP1)),MAXVAL(ZT_G(:,JLOOP1)),MAXVAL(ZQ_G(:,JLOOP1))
!
! Compute Theta V and relative humidity on grib grid
!
@@ -1039,17 +1041,13 @@ IF (IMODEL/=10) THEN ! others than NCEP
CALL ARRAY_1D_TO_2D (INO,ZOUT,IIU,IJU,ZTHV_LS(:,:,JLOOP1))
!
END DO
-ELSE !NCEP
+ELSE !GFS and ERA5 on pressure levels
DO JLOOP1=1, INLEVEL
- !WRITE (ILUOUT0,*) 'JLOOP1=',JLOOP1,MINVAL(ZPM_G(:,JLOOP1)),MINVAL(ZT_G(:,JLOOP1)),MINVAL(ZQ_G(:,JLOOP1))
- !WRITE (ILUOUT0,*) ' ',MAXVAL(ZPM_G(:,JLOOP1)),MAXVAL(ZT_G(:,JLOOP1)),MAXVAL(ZQ_G(:,JLOOP1))
ZH_G(:) =ZQ_G(:,JLOOP1)
ZRV_G(:) = (XRD/XRV)*SM_FOES(ZT_G(:,JLOOP1))*0.01*ZH_G(:) &
/(ZPM_G(:,JLOOP1) -SM_FOES(ZT_G(:,JLOOP1))*0.01*ZH_G(:))
- !WRITE (ILUOUT0,*) ' ',MINVAL(ZRV_G(:)),MAXVAL(ZRV_G(:))
ZTHV_G(:)=ZT_G(:,JLOOP1) * ((XP00/ZPM_G(:,JLOOP1))**(XRD/XCPD)) * &
((1. + XRV*ZRV_G(:)/XRD) / (1. + ZRV_G(:)))
- !WRITE (ILUOUT0,*) ' ',MINVAL(ZTHV_G(:)),MAXVAL(ZTHV_G(:))
!
! Interpolation : H
CALL HORIBL(ZPARAM(3),ZPARAM(4),ZPARAM(5),ZPARAM(6),INT(ZPARAM(2)),IINLO,INI, &
@@ -1140,6 +1138,9 @@ ALLOCATE (ZEXNM_LS(IIU,IJU,INLEVEL))
ZEXNM_LS(:,:,1:INLEVEL-1) = (ZEXNF_LS(:,:,1:INLEVEL-1)-ZEXNF_LS(:,:,2:INLEVEL)) / &
(LOG(ZEXNF_LS(:,:,1:INLEVEL-1))-LOG(ZEXNF_LS(:,:,2:INLEVEL)))
ZEXNM_LS(:,:,INLEVEL) = (ZPF_LS(:,:,INLEVEL)/2./XP00)**(XRD/XCPD)
+!
+IF (IMODEL==10.OR.IMODEL==11) ZEXNM_LS(:,:,:)=ZEXNF_LS(:,:,:) ! for GFS and ERA5 on pressure levels
+!
DEALLOCATE (ZEXNF_LS)
DEALLOCATE (ZPF_LS)
!
diff --git a/src/MNH/ver_prep_gribex_case.f90 b/src/MNH/ver_prep_gribex_case.f90
index 8a6cc54bb..0763ea885 100644
--- a/src/MNH/ver_prep_gribex_case.f90
+++ b/src/MNH/ver_prep_gribex_case.f90
@@ -87,6 +87,8 @@ END MODULE MODI_VER_PREP_GRIBEX_CASE
!! Philippe Wautelet: 05/2016-04/2018: new data structures and calls for I/O
!! Sep, 02, 2020 (Q. Rodier) use of geopotential height instead of
!! height above orography for isobaric surface interpolation
+!! Oct, 18 2022 (J.-P. Chaboureau) correction on vertical level
+!! for input fields on isobaric surface (GFS and ERA5)
!-------------------------------------------------------------------------------
!
!* 0. DECLARATIONS
@@ -142,6 +144,7 @@ REAL,DIMENSION(:,:,:), ALLOCATABLE :: ZTKE_LS
INTEGER :: JRR ! loop counter
INTEGER :: JSV ! loop counter
INTEGER :: JK ! loop counter
+CHARACTER(LEN=4) :: YWLOC_LS ! localisation of vertical wind speed
!-------------------------------------------------------------------------------
!
ILUOUT0 = TLUOUT0%NLU
@@ -195,6 +198,9 @@ END IF
ALLOCATE(ZZMASS_LS(IIU,IJU,ILU+2*JPVEXT))
ALLOCATE(ZSV_LS(IIU,IJU,ILU+2*JPVEXT,SIZE(XSV_LS,4)))
IF (HFILE(1:3)=='ATM') THEN
+ IF (SIZE(XB_LS)==0) THEN ! pressure level
+ XZMASS_LS(:,:,:)=XGH_LS(:,:,:)
+ END IF
ALLOCATE(ZZFLUX_LS(IIU,IJU,ILU+2*JPVEXT))
ALLOCATE(ZPMHP_LS(IIU,IJU,ILU+2*JPVEXT))
ALLOCATE(ZTHV_LS(IIU,IJU,ILU+2*JPVEXT))
@@ -262,27 +268,20 @@ IF (HFILE(1:3)=='ATM') THEN
END IF
!
IF (HFILE(1:3)=='ATM') THEN
-
- IF (SIZE(XB_LS)/=0) THEN ! hybrid level (w at flux points)
- CALL VER_INTERP_TO_MIXED_GRID('ATM ',.TRUE.,XZS_LS,XZSMT_LS, &
+ IF (SIZE(XB_LS)/=0) THEN ! hybrid level (w at flux points)
+ YWLOC_LS='FLUX'
+ ELSE ! isobaric surfaces (w at mass points)
+ YWLOC_LS='MASS'
+ END IF
+ !Warning, for pressure level (GFS and ERA5 only for now), ZZFLUX_LS is not correct (but not used)
+ CALL VER_INTERP_TO_MIXED_GRID('ATM ',.TRUE.,XZS_LS,XZSMT_LS, &
ZZMASS_LS,ZSV_LS, &
ZZFLUX_LS,XPS_LS,ZPMHP_LS, &
ZTHV_LS,ZR_LS, &
ZHU_LS, &
ZTKE_LS, &
ZU_LS,ZV_LS, &
- ZW_LS,'FLUX' )
- ELSE ! isobaric surfaces (w at mass points)
- !Warning, in that case (NCEP only for now) ZZFLUX_LS is not correct (but not used)
- CALL VER_INTERP_TO_MIXED_GRID('ATM ',.TRUE.,XZS_LS,XZSMT_LS, &
- XGH_LS,ZSV_LS, &
- ZZFLUX_LS,XPS_LS,ZPMHP_LS, &
- ZTHV_LS,ZR_LS, &
- ZHU_LS, &
- ZTKE_LS, &
- ZU_LS,ZV_LS, &
- ZW_LS,'MASS' )
- END IF
+ ZW_LS,YWLOC_LS )
ELSE IF (HFILE=='CHEM') THEN
CALL VER_INTERP_TO_MIXED_GRID(HFILE,.TRUE.,XZS_SV_LS,XZS_SV_LS,&
ZZMASS_LS,ZSV_LS )
--
GitLab