diff --git a/src/LIB/SURCOUCHE/src/modd_io.f90 b/src/LIB/SURCOUCHE/src/modd_io.f90 index 1602c603008a5cfb69271edc35f7f3c762cd30e8..100ea3c8ed6ea5c0b7765f64175e4a1ea607e876 100644 --- a/src/LIB/SURCOUCHE/src/modd_io.f90 +++ b/src/LIB/SURCOUCHE/src/modd_io.f90 @@ -9,7 +9,7 @@ MODULE MODD_IO_ll ! -USE MODD_NETCDF, ONLY: IDCDF_KIND, IOCDF +USE MODD_NETCDF, ONLY: IDCDF_KIND, IOCDF, TPTR2DIMCDF USE MODD_PARAMETERS, ONLY: NDIRNAMELGTMAX, NFILENAMELGTMAX ! IMPLICIT NONE @@ -100,6 +100,7 @@ TYPE TFILEDATA LOGICAL :: LNCCOMPRESS = .FALSE. ! Do compression on fields INTEGER(KIND=IDCDF_KIND) :: NNCCOMPRESS_LEVEL = 0 ! Compression level TYPE(IOCDF),POINTER :: TNCDIMS => NULL() ! Structure containing netCDF dimensions + TYPE(TPTR2DIMCDF),DIMENSION(:,:),ALLOCATABLE :: TNCCOORDS ! Structure pointing to the coordinates variables ! !Fields for other files INTEGER :: NLU = -1 !Logical unit number diff --git a/src/LIB/SURCOUCHE/src/modd_netcdf.f90 b/src/LIB/SURCOUCHE/src/modd_netcdf.f90 index c8c5c8a1f7185a92f72500373b1da6e89f538b9a..24541396d350d2544517e86aaef83918b1bc81f0 100644 --- a/src/LIB/SURCOUCHE/src/modd_netcdf.f90 +++ b/src/LIB/SURCOUCHE/src/modd_netcdf.f90 @@ -32,7 +32,10 @@ TYPE DIMCDF TYPE(DIMCDF), POINTER :: NEXT => NULL() END TYPE DIMCDF -TYPE(DIMCDF),DIMENSION(3,0:8),TARGET :: NCOORDID !X,Y,Z coordinates for the Arakawa points - !0 2nd-dimension is to treat NGRID=0 case without crash +TYPE TPTR2DIMCDF + TYPE(DIMCDF),POINTER :: TDIM +END TYPE TPTR2DIMCDF + +TYPE(DIMCDF),TARGET :: TDIM_DUMMY !Dummy dimension END MODULE MODD_NETCDF diff --git a/src/LIB/SURCOUCHE/src/mode_netcdf.f90 b/src/LIB/SURCOUCHE/src/mode_netcdf.f90 index 14442683d1ec466bef0b6dfe97196c8f15ecf60e..49d0b1499f5b2feac88d18d7723366e6c84f5f0d 100644 --- a/src/LIB/SURCOUCHE/src/mode_netcdf.f90 +++ b/src/LIB/SURCOUCHE/src/mode_netcdf.f90 @@ -125,7 +125,7 @@ USE MODD_CONF_n, ONLY: CSTORAGE_TYPE USE MODD_DIM_n, ONLY: NIMAX_ll, NJMAX_ll, NKMAX USE MODD_PARAMETERS_ll, ONLY: JPHEXT, JPVEXT -TYPE(TFILEDATA),INTENT(IN) :: TPFILE +TYPE(TFILEDATA),INTENT(INOUT) :: TPFILE CHARACTER(LEN=*),OPTIONAL,INTENT(IN) :: HPROGRAM_ORIG !To emulate a file coming from this program CHARACTER(LEN=:),ALLOCATABLE :: YPROGRAM @@ -166,39 +166,46 @@ ELSE IF (.NOT. ASSOCIATED(PIOCDF%DIMTIME)) ALLOCATE(PIOCDF%DIMTIME) END IF -! Store X,Y,Z coordinates for the Arakawa points +!Store X,Y,Z coordinates for the Arakawa points +!0 2nd-dimension is to treat NGRID=0 case without crash +IF (.NOT.ALLOCATED(TPFILE%TNCCOORDS)) ALLOCATE(TPFILE%TNCCOORDS(3,0:8)) +!Dummy point +TPFILE%TNCCOORDS(1,0)%TDIM => TDIM_DUMMY +TPFILE%TNCCOORDS(2,0)%TDIM => TDIM_DUMMY +TPFILE%TNCCOORDS(3,0)%TDIM => TDIM_DUMMY ! Mass point -NCOORDID(1,1) = PIOCDF%DIM_NI -NCOORDID(2,1) = PIOCDF%DIM_NJ -NCOORDID(3,1) = PIOCDF%DIM_LEVEL +TPFILE%TNCCOORDS(1,1)%TDIM => PIOCDF%DIM_NI +TPFILE%TNCCOORDS(2,1)%TDIM => PIOCDF%DIM_NJ +TPFILE%TNCCOORDS(3,1)%TDIM => PIOCDF%DIM_LEVEL ! u point -NCOORDID(1,2) = PIOCDF%DIM_NI_U -NCOORDID(2,2) = PIOCDF%DIM_NJ_U -NCOORDID(3,2) = PIOCDF%DIM_LEVEL +TPFILE%TNCCOORDS(1,2)%TDIM => PIOCDF%DIM_NI_U +TPFILE%TNCCOORDS(2,2)%TDIM => PIOCDF%DIM_NJ_U +TPFILE%TNCCOORDS(3,2)%TDIM => PIOCDF%DIM_LEVEL ! v point -NCOORDID(1,3) = PIOCDF%DIM_NI_V -NCOORDID(2,3) = PIOCDF%DIM_NJ_V -NCOORDID(3,3) = PIOCDF%DIM_LEVEL +TPFILE%TNCCOORDS(1,3)%TDIM => PIOCDF%DIM_NI_V +TPFILE%TNCCOORDS(2,3)%TDIM => PIOCDF%DIM_NJ_V +TPFILE%TNCCOORDS(3,3)%TDIM => PIOCDF%DIM_LEVEL ! w point -NCOORDID(1,4) = PIOCDF%DIM_NI -NCOORDID(2,4) = PIOCDF%DIM_NJ -NCOORDID(3,4) = PIOCDF%DIM_LEVEL_W +TPFILE%TNCCOORDS(1,4)%TDIM => PIOCDF%DIM_NI +TPFILE%TNCCOORDS(2,4)%TDIM => PIOCDF%DIM_NJ +TPFILE%TNCCOORDS(3,4)%TDIM => PIOCDF%DIM_LEVEL_W ! xi vorticity point (=f point =uv point) -NCOORDID(1,5) = PIOCDF%DIM_NI_U -NCOORDID(2,5) = PIOCDF%DIM_NJ_V -NCOORDID(3,5) = PIOCDF%DIM_LEVEL +TPFILE%TNCCOORDS(1,5)%TDIM => PIOCDF%DIM_NI_U +TPFILE%TNCCOORDS(2,5)%TDIM => PIOCDF%DIM_NJ_V +TPFILE%TNCCOORDS(3,5)%TDIM => PIOCDF%DIM_LEVEL ! eta vorticity point (=uw point) -NCOORDID(1,6) = PIOCDF%DIM_NI_U -NCOORDID(2,6) = PIOCDF%DIM_NJ_U -NCOORDID(3,6) = PIOCDF%DIM_LEVEL_W +TPFILE%TNCCOORDS(1,6)%TDIM => PIOCDF%DIM_NI_U +TPFILE%TNCCOORDS(2,6)%TDIM => PIOCDF%DIM_NJ_U +TPFILE%TNCCOORDS(3,6)%TDIM => PIOCDF%DIM_LEVEL_W ! zeta vorticity point (=vw point) -NCOORDID(1,7) = PIOCDF%DIM_NI_V -NCOORDID(2,7) = PIOCDF%DIM_NJ_V -NCOORDID(3,7) = PIOCDF%DIM_LEVEL_W +TPFILE%TNCCOORDS(1,7)%TDIM => PIOCDF%DIM_NI_V +TPFILE%TNCCOORDS(2,7)%TDIM => PIOCDF%DIM_NJ_V +TPFILE%TNCCOORDS(3,7)%TDIM => PIOCDF%DIM_LEVEL_W ! fw point (=uvw point) -NCOORDID(1,8) = PIOCDF%DIM_NI_U -NCOORDID(2,8) = PIOCDF%DIM_NJ_V -NCOORDID(3,8) = PIOCDF%DIM_LEVEL_W +TPFILE%TNCCOORDS(1,8)%TDIM => PIOCDF%DIM_NI_U +TPFILE%TNCCOORDS(2,8)%TDIM => PIOCDF%DIM_NJ_V +TPFILE%TNCCOORDS(3,8)%TDIM => PIOCDF%DIM_LEVEL_W + END SUBROUTINE IO_SET_KNOWNDIMS_NC4 @@ -795,7 +802,7 @@ IF (.NOT.GISCOORD) THEN !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 + IF (KSHAPE(1)==TPFILE%TNCCOORDS(1,TPFIELD%NGRID)%TDIM%LEN .AND. KSHAPE(2)==TPFILE%TNCCOORDS(2,TPFIELD%NGRID)%TDIM%LEN ) THEN SELECT CASE(TPFIELD%NGRID) CASE (0) !Not on Arakawa grid !Nothing to do @@ -1042,24 +1049,24 @@ END IF ! DO JI=1,SIZE(KSHAPE) IF (JI == 1) THEN - IF ( (YDIR == 'XX' .OR. YDIR == 'XY') .AND. KSHAPE(1)==NCOORDID(1,IGRID)%LEN) THEN - KVDIMS(1) = NCOORDID(1,IGRID)%ID - ELSE IF ( YDIR == 'YY' .AND. KSHAPE(1)==NCOORDID(2,IGRID)%LEN) THEN - KVDIMS(1) = NCOORDID(2,IGRID)%ID - ELSE IF ( YDIR == 'ZZ' .AND. KSHAPE(1)==NCOORDID(3,IGRID)%LEN) THEN - KVDIMS(1) = NCOORDID(3,IGRID)%ID + IF ( (YDIR == 'XX' .OR. YDIR == 'XY') .AND. KSHAPE(1)==TPFILE%TNCCOORDS(1,IGRID)%TDIM%LEN) THEN + KVDIMS(1) = TPFILE%TNCCOORDS(1,IGRID)%TDIM%ID + ELSE IF ( YDIR == 'YY' .AND. KSHAPE(1)==TPFILE%TNCCOORDS(2,IGRID)%TDIM%LEN) THEN + KVDIMS(1) = TPFILE%TNCCOORDS(2,IGRID)%TDIM%ID + ELSE IF ( YDIR == 'ZZ' .AND. KSHAPE(1)==TPFILE%TNCCOORDS(3,IGRID)%TDIM%LEN) THEN + KVDIMS(1) = TPFILE%TNCCOORDS(3,IGRID)%TDIM%ID ELSE PTDIM => GETDIMCDF(TPFILE, KSHAPE(1)); KVDIMS(1) = PTDIM%ID END IF ELSE IF (JI == 2) THEN - IF ( YDIR == 'XY' .AND. KSHAPE(2)==NCOORDID(2,IGRID)%LEN) THEN - KVDIMS(2) = NCOORDID(2,IGRID)%ID + IF ( YDIR == 'XY' .AND. KSHAPE(2)==TPFILE%TNCCOORDS(2,IGRID)%TDIM%LEN) THEN + KVDIMS(2) = TPFILE%TNCCOORDS(2,IGRID)%TDIM%ID ELSE PTDIM => GETDIMCDF(TPFILE, KSHAPE(2)); KVDIMS(2) = PTDIM%ID END IF ELSE IF (JI == 3) THEN - IF ( YDIR == 'XY' .AND. KSHAPE(3)==NCOORDID(3,IGRID)%LEN) THEN - KVDIMS(3) = NCOORDID(3,IGRID)%ID + IF ( YDIR == 'XY' .AND. KSHAPE(3)==TPFILE%TNCCOORDS(3,IGRID)%TDIM%LEN) THEN + KVDIMS(3) = TPFILE%TNCCOORDS(3,IGRID)%TDIM%ID ELSE PTDIM => GETDIMCDF(TPFILE, KSHAPE(3)); KVDIMS(3) = PTDIM%ID END IF @@ -1131,9 +1138,9 @@ ELSE IF (TPFIELD%CLBTYPE/='NONE') THEN END IF ! IF (TPFIELD%CLBTYPE=='LBX' .OR. TPFIELD%CLBTYPE=='LBXU') THEN - PTDIM => NCOORDID(2,IGRID) + PTDIM => TPFILE%TNCCOORDS(2,IGRID)%TDIM TPDIMS(2) = PTDIM - PTDIM => NCOORDID(3,IGRID) + PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM TPDIMS(3) = PTDIM ILEN = TPDIMS(2)%LEN * TPDIMS(3)%LEN ISIZE = KLEN/ILEN @@ -1143,9 +1150,9 @@ ELSE IF (TPFIELD%CLBTYPE/='NONE') THEN TPDIMS(1) = PTDIM ILEN = ILEN * PTDIM%LEN ELSE IF (TPFIELD%CLBTYPE=='LBY' .OR. TPFIELD%CLBTYPE=='LBYV') THEN - PTDIM => NCOORDID(1,IGRID) + PTDIM => TPFILE%TNCCOORDS(1,IGRID)%TDIM TPDIMS(1) = PTDIM - PTDIM => NCOORDID(3,IGRID) + PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM TPDIMS(3) = PTDIM ILEN = TPDIMS(1)%LEN * TPDIMS(3)%LEN ISIZE = KLEN/ILEN @@ -1164,11 +1171,11 @@ ELSE DO JI=1,TPFIELD%NDIMS IF (JI == 1) THEN IF ( (YDIR == 'XX' .OR. YDIR == 'XY') ) THEN - PTDIM => NCOORDID(1,IGRID) + PTDIM => TPFILE%TNCCOORDS(1,IGRID)%TDIM ELSE IF ( YDIR == 'YY' ) THEN - PTDIM => NCOORDID(2,IGRID) + PTDIM => TPFILE%TNCCOORDS(2,IGRID)%TDIM ELSE IF ( YDIR == 'ZZ' ) THEN - PTDIM => NCOORDID(3,IGRID) + PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM ELSE CALL PRINT_MSG(NVERB_WARNING,'IO','IO_GUESS_DIMIDS_NC4','can not guess 1st dimension for field '//TRIM(TPFIELD%CMNHNAME)) END IF @@ -1176,7 +1183,7 @@ ELSE TPDIMS(JI) = PTDIM ELSE IF (JI == 2) THEN IF ( YDIR == 'XY') THEN - PTDIM => NCOORDID(2,IGRID) + PTDIM => TPFILE%TNCCOORDS(2,IGRID)%TDIM ELSE IF (JI==TPFIELD%NDIMS) THEN !Guess last dimension ISIZE = KLEN/ILEN IF (MOD(KLEN,ILEN)/=0) CALL PRINT_MSG(NVERB_WARNING,'IO','IO_GUESS_DIMIDS_NC4', & @@ -1189,7 +1196,7 @@ ELSE TPDIMS(JI) = PTDIM ELSE IF (JI == 3) THEN IF ( YDIR == 'XY' ) THEN - PTDIM => NCOORDID(3,IGRID) + PTDIM => TPFILE%TNCCOORDS(3,IGRID)%TDIM ELSE IF (JI==TPFIELD%NDIMS) THEN !Guess last dimension ISIZE = KLEN/ILEN IF (MOD(KLEN,ILEN)/=0) CALL PRINT_MSG(NVERB_WARNING,'IO','IO_GUESS_DIMIDS_NC4', &