From 6cd5e75e8bacd229978c807dad6f3677a3e41837 Mon Sep 17 00:00:00 2001
From: Philippe WAUTELET <philippe.wautelet@aero.obs-mip.fr>
Date: Tue, 23 Jan 2018 17:17:42 +0100
Subject: [PATCH] Philippe 23/01/2018: IO: * add latitude(_x) and longitude(_x)
fields in NC4 files * set the coordinates attributes for fields * commented
out CSTDNAME for LAT and LON fields to not have 2 times the same
standard_name * corrected NGRID for DX/YHAT fields * modify condition to
detect if program is PREP_SURFEX
---
src/LIB/SURCOUCHE/src/mode_field.f90 | 114 ++++++++++++++++++-
src/LIB/SURCOUCHE/src/mode_netcdf.f90 | 156 ++++++++++++++++++++++----
2 files changed, 245 insertions(+), 25 deletions(-)
diff --git a/src/LIB/SURCOUCHE/src/mode_field.f90 b/src/LIB/SURCOUCHE/src/mode_field.f90
index 3492d7b41..b9a239245 100644
--- a/src/LIB/SURCOUCHE/src/mode_field.f90
+++ b/src/LIB/SURCOUCHE/src/mode_field.f90
@@ -530,7 +530,7 @@ TFIELDLIST(IDX)%CLONGNAME = 'MesoNH: DXHAT'
TFIELDLIST(IDX)%CUNITS = 'm'
TFIELDLIST(IDX)%CDIR = 'XX'
TFIELDLIST(IDX)%CCOMMENT = 'Horizontal stretching in x'
-TFIELDLIST(IDX)%NGRID = 4
+TFIELDLIST(IDX)%NGRID = 2
TFIELDLIST(IDX)%NTYPE = TYPEREAL
TFIELDLIST(IDX)%NDIMS = 1
ALLOCATE(TFIELDLIST(IDX)%TFIELD_X1D(IMODEL))
@@ -544,7 +544,7 @@ TFIELDLIST(IDX)%CLONGNAME = 'MesoNH: DYHAT'
TFIELDLIST(IDX)%CUNITS = 'm'
TFIELDLIST(IDX)%CDIR = 'YY'
TFIELDLIST(IDX)%CCOMMENT = 'Horizontal stretching in y'
-TFIELDLIST(IDX)%NGRID = 4
+TFIELDLIST(IDX)%NGRID = 3
TFIELDLIST(IDX)%NTYPE = TYPEREAL
TFIELDLIST(IDX)%NDIMS = 1
ALLOCATE(TFIELDLIST(IDX)%TFIELD_X1D(IMODEL))
@@ -643,8 +643,113 @@ ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
IDX = IDX+1
!
IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
-TFIELDLIST(IDX)%CMNHNAME = 'LAT'
+TFIELDLIST(IDX)%CMNHNAME = 'latitude'
TFIELDLIST(IDX)%CSTDNAME = 'latitude'
+TFIELDLIST(IDX)%CLONGNAME = 'latitude'
+TFIELDLIST(IDX)%CUNITS = 'degrees_north'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_latitude at mass point'
+TFIELDLIST(IDX)%NGRID = 1
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'longitude'
+TFIELDLIST(IDX)%CSTDNAME = 'longitude'
+TFIELDLIST(IDX)%CLONGNAME = 'longitude'
+TFIELDLIST(IDX)%CUNITS = 'degrees_east'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_longitude at mass point'
+TFIELDLIST(IDX)%NGRID = 1
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'latitude_u'
+TFIELDLIST(IDX)%CSTDNAME = 'latitude_at_u_location'
+TFIELDLIST(IDX)%CLONGNAME = 'latitude at u location'
+TFIELDLIST(IDX)%CUNITS = 'degrees_north'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_latitude at u point'
+TFIELDLIST(IDX)%NGRID = 2
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'longitude_u'
+TFIELDLIST(IDX)%CSTDNAME = 'longitude_at_u_location'
+TFIELDLIST(IDX)%CLONGNAME = 'longitude at u location'
+TFIELDLIST(IDX)%CUNITS = 'degrees_east'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_longitude at u point'
+TFIELDLIST(IDX)%NGRID = 2
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'latitude_v'
+TFIELDLIST(IDX)%CSTDNAME = 'latitude_at_v_location'
+TFIELDLIST(IDX)%CLONGNAME = 'latitude at v location'
+TFIELDLIST(IDX)%CUNITS = 'degrees_north'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_latitude at v point'
+TFIELDLIST(IDX)%NGRID = 3
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'longitude_v'
+TFIELDLIST(IDX)%CSTDNAME = 'longitude_at_v_location'
+TFIELDLIST(IDX)%CLONGNAME = 'longitude at v location'
+TFIELDLIST(IDX)%CUNITS = 'degrees_east'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_longitude at v point'
+TFIELDLIST(IDX)%NGRID = 3
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'latitude_f'
+TFIELDLIST(IDX)%CSTDNAME = 'latitude_at_f_location'
+TFIELDLIST(IDX)%CLONGNAME = 'latitude at f location'
+TFIELDLIST(IDX)%CUNITS = 'degrees_north'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_latitude at f point'
+TFIELDLIST(IDX)%NGRID = 5
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'longitude_f'
+TFIELDLIST(IDX)%CSTDNAME = 'longitude_at_f_location'
+TFIELDLIST(IDX)%CLONGNAME = 'longitude at f location'
+TFIELDLIST(IDX)%CUNITS = 'degrees_east'
+TFIELDLIST(IDX)%CDIR = 'XY'
+TFIELDLIST(IDX)%CCOMMENT = 'X_Y_longitude at f point'
+TFIELDLIST(IDX)%NGRID = 5
+TFIELDLIST(IDX)%NTYPE = TYPEREAL
+TFIELDLIST(IDX)%NDIMS = 2
+ALLOCATE(TFIELDLIST(IDX)%TFIELD_X2D(IMODEL))
+IDX = IDX+1
+!
+IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
+TFIELDLIST(IDX)%CMNHNAME = 'LAT'
+! TFIELDLIST(IDX)%CSTDNAME = 'latitude'
+TFIELDLIST(IDX)%CSTDNAME = ''
TFIELDLIST(IDX)%CLONGNAME = 'MesoNH: LAT'
TFIELDLIST(IDX)%CUNITS = 'degrees_north'
TFIELDLIST(IDX)%CDIR = 'XY'
@@ -657,7 +762,8 @@ IDX = IDX+1
!
IF(IDX>MAXFIELDS) CALL ERR_INI_FIELD_LIST()
TFIELDLIST(IDX)%CMNHNAME = 'LON'
-TFIELDLIST(IDX)%CSTDNAME = 'longitude'
+! TFIELDLIST(IDX)%CSTDNAME = 'longitude'
+TFIELDLIST(IDX)%CSTDNAME = ''
TFIELDLIST(IDX)%CLONGNAME = 'MesoNH: LON'
TFIELDLIST(IDX)%CUNITS = 'degrees_east'
TFIELDLIST(IDX)%CDIR = 'XY'
diff --git a/src/LIB/SURCOUCHE/src/mode_netcdf.f90 b/src/LIB/SURCOUCHE/src/mode_netcdf.f90
index bda97e58d..59e962a95 100644
--- a/src/LIB/SURCOUCHE/src/mode_netcdf.f90
+++ b/src/LIB/SURCOUCHE/src/mode_netcdf.f90
@@ -179,16 +179,23 @@ END SUBROUTINE IO_SET_KNOWNDIMS_NC4
SUBROUTINE IO_WRITE_COORDVAR_NC4(TPFILE)
USE MODD_CONF, ONLY: CPROGRAM, LCARTESIAN
+USE MODD_CONF_n, ONLY: CSTORAGE_TYPE
+USE MODD_GRID, ONLY: XLATORI, XLONORI
USE MODD_GRID_n, ONLY: LSLEVE, XXHAT, XYHAT, XZHAT
USE MODD_PARAMETERS, ONLY: JPHEXT, JPVEXT
+USE MODE_FIELD, ONLY: TFIELDLIST,FIND_FIELD_ID_FROM_MNHNAME
+USE MODE_GRIDPROJ
+
TYPE(TFILEDATA),INTENT(IN) :: TPFILE
-CHARACTER(LEN=:),ALLOCATABLE :: YSTDNAMEPREFIX
-INTEGER :: IIU, IJU, IKU
-INTEGER(KIND=IDCDF_KIND) :: INCID
-REAL,DIMENSION(:),ALLOCATABLE :: ZXHATM, ZYHATM,ZZHATM !Coordinates at mass points in the transformed space
-TYPE(IOCDF), POINTER :: PIOCDF
+CHARACTER(LEN=:),ALLOCATABLE :: YSTDNAMEPREFIX
+INTEGER :: IIU, IJU, IKU
+INTEGER :: ID, IRESP
+INTEGER(KIND=IDCDF_KIND) :: INCID
+REAL,DIMENSION(:),ALLOCATABLE :: ZXHATM, ZYHATM,ZZHATM !Coordinates at mass points in the transformed space
+REAL,DIMENSION(:,:),ALLOCATABLE :: ZLAT, ZLON
+TYPE(IOCDF), POINTER :: PIOCDF
CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_WRITE_COORDVAR_NC4','called for'//TRIM(TPFILE%CNAME))
@@ -199,7 +206,7 @@ INCID = TPFILE%NNCID
IF (TPFILE%LMASTER) THEN
IF ( TRIM(CPROGRAM)/='PGD' &
- .AND. .NOT.(TRIM(CPROGRAM)=='REAL'.AND. .NOT.ASSOCIATED(XXHAT)) ) THEN !condition to detect PREP_SURFEX
+ .AND. .NOT.(TRIM(CPROGRAM)=='REAL' .AND. CSTORAGE_TYPE=='SU') ) THEN !condition to detect PREP_SURFEX
IIU = SIZE(XXHAT)
IJU = SIZE(XYHAT)
ALLOCATE(ZXHATM(IIU),ZYHATM(IJU))
@@ -224,6 +231,51 @@ IF (TPFILE%LMASTER) THEN
CALL WRITE_HOR_COORD(PIOCDF%DIM_NJ_V,'y-dimension of the grid at v location', &
TRIM(YSTDNAMEPREFIX)//'_y_coordinate_at_v_location','Y',-0.5,JPHEXT,0, XYHAT)
!
+ IF (.NOT.LCARTESIAN) THEN
+ ALLOCATE(ZLAT(IIU,IJU),ZLON(IIU,IJU))
+ !
+ !Compute latitude/longitude for the Arakawa points
+ !
+ ! Mass point
+ CALL SM_LATLON(XLATORI,XLONORI, &
+ SPREAD(SOURCE=ZXHATM,DIM=2,NCOPIES=IJU), &
+ SPREAD(SOURCE=ZYHATM,DIM=1,NCOPIES=IIU), &
+ ZLAT,ZLON)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('latitude',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLAT,IRESP,OISCOORD=.TRUE.)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('longitude',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLON,IRESP,OISCOORD=.TRUE.)
+ ! u point
+ CALL SM_LATLON(XLATORI,XLONORI, &
+ SPREAD(SOURCE=XXHAT, DIM=2,NCOPIES=IJU), &
+ SPREAD(SOURCE=ZYHATM,DIM=1,NCOPIES=IIU), &
+ ZLAT,ZLON)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('latitude_u',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLAT,IRESP,OISCOORD=.TRUE.)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('longitude_u',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLON,IRESP,OISCOORD=.TRUE.)
+ ! v point
+ CALL SM_LATLON(XLATORI,XLONORI, &
+ SPREAD(SOURCE=ZXHATM,DIM=2,NCOPIES=IJU), &
+ SPREAD(SOURCE=XYHAT, DIM=1,NCOPIES=IIU), &
+ ZLAT,ZLON)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('latitude_v',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLAT,IRESP,OISCOORD=.TRUE.)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('longitude_v',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLON,IRESP,OISCOORD=.TRUE.)
+ ! xi vorticity point (=f point =uv point)
+ CALL SM_LATLON(XLATORI,XLONORI, &
+ SPREAD(SOURCE=XXHAT,DIM=2,NCOPIES=IJU), &
+ SPREAD(SOURCE=XYHAT,DIM=1,NCOPIES=IIU), &
+ ZLAT,ZLON)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('latitude_f',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLAT,IRESP,OISCOORD=.TRUE.)
+ CALL FIND_FIELD_ID_FROM_MNHNAME('longitude_f',ID,IRESP)
+ CALL IO_WRITE_FIELD_NC4_X2(TPFILE,TFIELDLIST(ID),ZLON,IRESP,OISCOORD=.TRUE.)
+ !
+ DEALLOCATE(ZLAT,ZLON)
+ END IF
+ !
DEALLOCATE(ZXHATM,ZYHATM)
ELSE
CALL WRITE_HOR_COORD(PIOCDF%DIM_NI,'x-dimension of the grid','x_grid_index','X',0.,JPHEXT,JPHEXT)
@@ -237,7 +289,7 @@ IF (TPFILE%LMASTER) THEN
ENDIF
!
IF (TRIM(CPROGRAM)/='PGD' .AND. TRIM(CPROGRAM)/='NESPGD' &
- .AND. .NOT.(TRIM(CPROGRAM)=='REAL'.AND. .NOT.ASSOCIATED(XXHAT)) ) THEN !condition to detect PREP_SURFEX
+ .AND. .NOT.(TRIM(CPROGRAM)=='REAL' .AND. CSTORAGE_TYPE=='SU') ) THEN !condition to detect PREP_SURFEX
!
IKU = SIZE(XZHAT)
ALLOCATE(ZZHATM(IKU))
@@ -478,14 +530,21 @@ END IF
END SUBROUTINE IO_WRITE_HEADER_NC4
-SUBROUTINE IO_WRITE_FIELD_ATTR_NC4(TPFIELD,KNCID,KVARID,HCALENDAR)
+SUBROUTINE IO_WRITE_FIELD_ATTR_NC4(TPFIELD,KNCID,KVARID,KSHAPE,HCALENDAR,OISCOORD)
!
-TYPE(TFIELDDATA), INTENT(IN) :: TPFIELD
-INTEGER(KIND=IDCDF_KIND), INTENT(IN) :: KNCID
-INTEGER(KIND=IDCDF_KIND), INTENT(IN) :: KVARID
-CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: HCALENDAR
+USE MODD_CONF, ONLY: CPROGRAM, LCARTESIAN
+USE MODD_CONF_n, ONLY: CSTORAGE_TYPE
!
-INTEGER(KIND=IDCDF_KIND) :: STATUS
+TYPE(TFIELDDATA), INTENT(IN) :: TPFIELD
+INTEGER(KIND=IDCDF_KIND), INTENT(IN) :: KNCID
+INTEGER(KIND=IDCDF_KIND), INTENT(IN) :: KVARID
+INTEGER(KIND=IDCDF_KIND),DIMENSION(:),OPTIONAL,INTENT(IN) :: KSHAPE
+CHARACTER(LEN=*), OPTIONAL,INTENT(IN) :: HCALENDAR
+LOGICAL, OPTIONAL,INTENT(IN) :: OISCOORD ! Is a coordinate variable (->do not write coordinates attribute)
+!
+INTEGER(KIND=IDCDF_KIND) :: STATUS
+CHARACTER(LEN=:),ALLOCATABLE :: YCOORDS
+LOGICAL :: GISCOORD
!
CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_WRITE_FIELD_ATTR_NC4','called for field '//TRIM(TPFIELD%CMNHNAME))
!
@@ -494,6 +553,15 @@ IF(LEN_TRIM(TPFIELD%CSTDNAME)==0 .AND. LEN_TRIM(TPFIELD%CLONGNAME)==0) THEN
&to respect CF-convention for variable '//TRIM(TPFIELD%CMNHNAME))
ENDIF
!
+IF (TPFIELD%NDIMS>1 .AND. .NOT.PRESENT(KSHAPE)) &
+ CALL PRINT_MSG(NVERB_FATAL,'IO','IO_WRITE_FIELD_ATTR_NC4','KSHAPE not provided for '//TRIM(TPFIELD%CMNHNAME))
+!
+IF (PRESENT(OISCOORD)) THEN
+ GISCOORD = OISCOORD
+ELSE
+ GISCOORD = .FALSE.
+END IF
+!
! GRID attribute definition
IF(TPFIELD%NGRID<0) THEN
CALL PRINT_MSG(NVERB_DEBUG,'IO','IO_WRITE_FIELD_ATTR_NC4','TPFIELD%NGRID not set for variable '//TRIM(TPFIELD%CMNHNAME))
@@ -542,6 +610,51 @@ IF(PRESENT(HCALENDAR)) THEN
IF (STATUS /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_ATTR_NC4 [NF90_PUT_ATT]')
ENDIF
!
+! Coordinates (CF convention)
+IF (.NOT.GISCOORD &
+ .AND. TRIM(CPROGRAM)/='PGD' &
+ .AND. .NOT.(TRIM(CPROGRAM)=='REAL' .AND. CSTORAGE_TYPE=='SU' ) ) THEN !condition to detect PREP_SURFEX
+ !0D: nothing to do
+ !1D: no direct correspondance with latitude(_x)/longitude(_x) 2D variables => nothing to do
+ IF (.NOT.LCARTESIAN .AND. TPFIELD%NDIMS>1 .AND. TPFIELD%NGRID/=0) THEN
+ IF (TPFIELD%CDIR=='XY') THEN
+ IF (KSHAPE(1)==NCOORDID(1,TPFIELD%NGRID)%LEN .AND. KSHAPE(2)==NCOORDID(2,TPFIELD%NGRID)%LEN ) THEN
+ SELECT CASE(TPFIELD%NGRID)
+ CASE (0) !Not on Arakawa grid
+ !Nothing to do
+ CASE (1) !Mass point
+ YCOORDS='latitude longitude'
+ CASE (2) !u point
+ YCOORDS='latitude_u longitude_u'
+ CASE (3) !v point
+ YCOORDS='latitude_v longitude_v'
+ CASE (4) !w point
+ YCOORDS='latitude longitude'
+ CASE (5) !xi vorticity point (=f point =uv point)
+ YCOORDS='latitude_f longitude_f'
+ CASE (6) !eta vorticity point (=uw point)
+ YCOORDS='latitude_u longitude_u'
+ CASE (7) !zeta vorticity point (=vw point)
+ YCOORDS='latitude_v longitude_v'
+ CASE (8) !fw point (=uvw point)
+ YCOORDS='latitude_f longitude_f'
+ CASE DEFAULT
+ CALL PRINT_MSG(NVERB_ERROR,'IO','IO_WRITE_FIELD_ATTR_NC4','invalid NGRID for variable '//TRIM(TPFIELD%CMNHNAME))
+ END SELECT
+ !
+ STATUS = NF90_PUT_ATT(KNCID, KVARID,'coordinates',YCOORDS)
+ IF (STATUS /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_ATTR_NC4 [NF90_PUT_ATT]')
+ DEALLOCATE(YCOORDS)
+ ELSE
+ CALL PRINT_MSG(NVERB_WARNING,'IO','IO_WRITE_FIELD_ATTR_NC4','coordinates not implemented for variable ' &
+ //TRIM(TPFIELD%CMNHNAME))
+ END IF
+ ELSE
+ !No YCOORDS for CDIR/='XY'
+ END IF
+ END IF
+ENDIF
+!
END SUBROUTINE IO_WRITE_FIELD_ATTR_NC4
FUNCTION GETDIMCDF(TPFILE, KLEN, HDIMNAME)
@@ -830,7 +943,7 @@ KRESP = IRESP
END SUBROUTINE IO_WRITE_FIELD_NC4_X1
-SUBROUTINE IO_WRITE_FIELD_NC4_X2(TPFILE,TPFIELD,PFIELD,KRESP,KVERTLEVEL,KZFILE)
+SUBROUTINE IO_WRITE_FIELD_NC4_X2(TPFILE,TPFIELD,PFIELD,KRESP,KVERTLEVEL,KZFILE,OISCOORD)
!
TYPE(TFILEDATA),TARGET,INTENT(IN) :: TPFILE
TYPE(TFIELDDATA), INTENT(IN) :: TPFIELD
@@ -838,6 +951,7 @@ REAL,DIMENSION(:,:), INTENT(IN) :: PFIELD ! array containing the data field
INTEGER, INTENT(OUT):: KRESP
INTEGER,OPTIONAL, INTENT(IN) :: KVERTLEVEL ! Number of the vertical level (needed for Z-level splitted files)
INTEGER,OPTIONAL, INTENT(IN) :: KZFILE ! Number of the Z-level splitted file
+LOGICAL,OPTIONAL, INTENT(IN) :: OISCOORD ! Is a coordinate variable (->do not write coordinates attribute)
!
INTEGER(KIND=IDCDF_KIND) :: STATUS
INTEGER(KIND=IDCDF_KIND) :: INCID
@@ -903,7 +1017,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(PFIELD),KIND=IDCDF_KIND),OISCOORD=OISCOORD)
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, PFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_X2[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
@@ -964,7 +1078,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(PFIELD),KIND=IDCDF_KIND))
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, PFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_X3[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
@@ -1025,7 +1139,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(PFIELD),KIND=IDCDF_KIND))
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, PFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_X4[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
@@ -1086,7 +1200,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(PFIELD),KIND=IDCDF_KIND))
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, PFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_X5[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
@@ -1147,7 +1261,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(PFIELD),KIND=IDCDF_KIND))
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, PFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_X6[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
@@ -1337,7 +1451,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(KFIELD),KIND=IDCDF_KIND))
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, KFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_N2[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
@@ -1397,7 +1511,7 @@ ELSE
END IF
! Write metadata
-CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID)
+CALL IO_WRITE_FIELD_ATTR_NC4(TPFIELD,INCID,IVARID,KSHAPE=INT(SHAPE(KFIELD),KIND=IDCDF_KIND))
! Write the data
STATUS = NF90_PUT_VAR(INCID, IVARID, KFIELD)
IF (status /= NF90_NOERR) CALL HANDLE_ERR(status,__LINE__,'IO_WRITE_FIELD_NC4_N3[NF90_PUT_VAR] '//TRIM(TPFIELD%CMNHNAME),IRESP)
--
GitLab