From 5a0f8548b18a397a636bb5821d4ce6a325d710f1 Mon Sep 17 00:00:00 2001 From: mjreno Date: Tue, 17 Dec 2024 17:04:00 -0500 Subject: [PATCH] refactor(netcdf): utility package input adjustments (#2106) * netcdf improved error handling and utl-ncf input adjustments * escape underscores * add autotest assert for netcdf data model * add baseline testing for utl-ncf options * remove in_record attribute from chunk params * set chunk_time param optional --------- Co-authored-by: mjreno --- autotest/test_netcdf_gwe_cnd.py | 5 ++ autotest/test_netcdf_gwf_disv.py | 21 ++++- .../test_netcdf_gwf_lak_wetlakbedarea02.py | 5 ++ autotest/test_netcdf_gwf_rch01.py | 5 ++ autotest/test_netcdf_gwf_rch03.py | 5 ++ autotest/test_netcdf_gwf_sto01.py | 46 +++++++++-- autotest/test_netcdf_gwf_vsc03_sfr.py | 5 ++ autotest/test_netcdf_gwt_dsp01.py | 39 +++++++-- autotest/test_netcdf_gwt_prudic2004t2.py | 5 ++ doc/mf6io/mf6ivar/dfn/utl-ncf.dfn | 42 +++------- .../mf6ivar/examples/utl-ncf-example.dat | 5 +- src/Idm/utl-ncfidm.f90 | 70 ++++------------ src/Utilities/Export/DisNCStructured.f90 | 80 +++++++++++-------- src/Utilities/Export/MeshNCModel.f90 | 26 +++--- src/Utilities/Export/ModelExport.f90 | 15 +++- src/Utilities/Export/NCModel.f90 | 19 +++-- 16 files changed, 234 insertions(+), 159 deletions(-) diff --git a/autotest/test_netcdf_gwe_cnd.py b/autotest/test_netcdf_gwe_cnd.py index 95a47210343..829bc78ace1 100644 --- a/autotest/test_netcdf_gwe_cnd.py +++ b/autotest/test_netcdf_gwe_cnd.py @@ -24,6 +24,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -57,6 +58,10 @@ def check_output(idx, test, export, gridded_input): name = "gwe-" + test.name + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{name}.nc") as ds: + assert ds.data_model == "NETCDF4" + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{name}.nc" diff --git a/autotest/test_netcdf_gwf_disv.py b/autotest/test_netcdf_gwf_disv.py index 78ae82e393b..121090657f8 100644 --- a/autotest/test_netcdf_gwf_disv.py +++ b/autotest/test_netcdf_gwf_disv.py @@ -14,6 +14,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") wkt = ( 'PROJCS["NAD83 / UTM zone 18N", ' @@ -53,7 +54,15 @@ def build_models(idx, test, export, gridded_input): gwf.name_file.nc_mesh2d_filerecord = f"{name}.nc" # netcdf config - ncf = flopy.mf6.ModflowUtlncf(gwf.disv, ogc_wkt=wkt, filename=f"{name}.disv.ncf") + ncf = flopy.mf6.ModflowUtlncf( + gwf.disv, + deflate=9, + shuffle=True, + chunk_time=1, + chunk_face=3, + wkt=wkt, + filename=f"{name}.disv.ncf", + ) # output control oc = flopy.mf6.ModflowGwfoc( @@ -72,6 +81,16 @@ def check_output(idx, test, export, gridded_input): name = test.name + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{name}.nc") as ds: + assert ds.data_model == "NETCDF4" + cmpr = ds.variables["head_l1"].filters() + chnk = ds.variables["head_l1"].chunking() + assert cmpr["shuffle"] + assert cmpr["complevel"] == 9 + assert chnk == [1, 3] + assert ds.variables["projection"].getncattr("wkt").lower() == wkt.lower() + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{name}.nc" diff --git a/autotest/test_netcdf_gwf_lak_wetlakbedarea02.py b/autotest/test_netcdf_gwf_lak_wetlakbedarea02.py index 025cf886cd3..afc94df4abd 100644 --- a/autotest/test_netcdf_gwf_lak_wetlakbedarea02.py +++ b/autotest/test_netcdf_gwf_lak_wetlakbedarea02.py @@ -24,6 +24,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -61,6 +62,10 @@ def check_output(idx, test, export, gridded_input): name = cases[idx] gwfname = "gwf-" + name + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{gwfname}.nc") as ds: + assert ds.data_model == "NETCDF4" + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{gwfname}.nc" diff --git a/autotest/test_netcdf_gwf_rch01.py b/autotest/test_netcdf_gwf_rch01.py index c54ab869334..f89c0b4853a 100644 --- a/autotest/test_netcdf_gwf_rch01.py +++ b/autotest/test_netcdf_gwf_rch01.py @@ -24,6 +24,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -58,6 +59,10 @@ def check_output(idx, test, export, gridded_input): name = "rch" + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{name}.nc") as ds: + assert ds.data_model == "NETCDF4" + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{name}.nc" diff --git a/autotest/test_netcdf_gwf_rch03.py b/autotest/test_netcdf_gwf_rch03.py index e473fdc5133..7b9b7dd74d7 100644 --- a/autotest/test_netcdf_gwf_rch03.py +++ b/autotest/test_netcdf_gwf_rch03.py @@ -24,6 +24,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -59,6 +60,10 @@ def check_output(idx, test, export, gridded_input): name = "rch" + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{name}.nc") as ds: + assert ds.data_model == "NETCDF4" + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{name}.nc" diff --git a/autotest/test_netcdf_gwf_sto01.py b/autotest/test_netcdf_gwf_sto01.py index c987c537d40..38fa7a1e69e 100644 --- a/autotest/test_netcdf_gwf_sto01.py +++ b/autotest/test_netcdf_gwf_sto01.py @@ -12,6 +12,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") htol = [None for _ in range(len(cases))] @@ -52,15 +53,28 @@ def build_models(idx, test, export, gridded_input): if export == "ugrid": gwf.name_file.nc_mesh2d_filerecord = f"{name}.nc" + ncf = flopy.mf6.ModflowUtlncf( + gwf.dis, + deflate=5, + shuffle=True, + chunk_time=1, + chunk_face=10, + wkt=wkt, + filename=f"{name}.dis.ncf", + ) elif export == "structured": gwf.name_file.nc_structured_filerecord = f"{name}.nc" - - # netcdf config - ncf = flopy.mf6.ModflowUtlncf( - gwf.dis, - ogc_wkt=wkt, - filename=f"{name}.dis.ncf", - ) + ncf = flopy.mf6.ModflowUtlncf( + gwf.dis, + deflate=5, + shuffle=True, + chunk_time=1, + chunk_z=1, + chunk_y=5, + chunk_x=5, + wkt=wkt, + filename=f"{name}.dis.ncf", + ) return sim, dummy @@ -68,6 +82,24 @@ def build_models(idx, test, export, gridded_input): def check_output(idx, test, export, gridded_input): from test_gwf_sto01 import check_output as check + # verify format of generated netcdf file + with nc.Dataset(test.workspace / "gwf_sto01.nc") as ds: + assert ds.data_model == "NETCDF4" + if export == "structured": + cmpr = ds.variables["head"].filters() + chnk = ds.variables["head"].chunking() + assert chnk == [1, 1, 5, 5] + assert ( + ds.variables["projection"].getncattr("crs_wkt").lower() == wkt.lower() + ) + elif export == "ugrid": + cmpr = ds.variables["head_l1"].filters() + chnk = ds.variables["head_l1"].chunking() + assert chnk == [1, 10] + assert ds.variables["projection"].getncattr("wkt").lower() == wkt.lower() + assert cmpr["shuffle"] + assert cmpr["complevel"] == 5 + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = "gwf_sto01.nc" diff --git a/autotest/test_netcdf_gwf_vsc03_sfr.py b/autotest/test_netcdf_gwf_vsc03_sfr.py index 50be2d19098..97a7fcc0b13 100644 --- a/autotest/test_netcdf_gwf_vsc03_sfr.py +++ b/autotest/test_netcdf_gwf_vsc03_sfr.py @@ -24,6 +24,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -61,6 +62,10 @@ def check_output(idx, test, export, gridded_input): name = "gwf-" + test.name + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{name}.nc") as ds: + assert ds.data_model == "NETCDF4" + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{name}.nc" diff --git a/autotest/test_netcdf_gwt_dsp01.py b/autotest/test_netcdf_gwt_dsp01.py index fd37f1f8c09..621b359d3c8 100644 --- a/autotest/test_netcdf_gwt_dsp01.py +++ b/autotest/test_netcdf_gwt_dsp01.py @@ -12,6 +12,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -29,14 +30,26 @@ def build_models(idx, test, export, gridded_input): if export == "ugrid": gwt.name_file.nc_mesh2d_filerecord = f"{gwtname}.nc" + ncf = flopy.mf6.ModflowUtlncf( + gwt.dis, + deflate=3, + shuffle=False, + chunk_time=1, + chunk_face=5, + filename=f"{gwtname}.dis.ncf", + ) elif export == "structured": gwt.name_file.nc_structured_filerecord = f"{gwtname}.nc" - - # netcdf config - ncf = flopy.mf6.ModflowUtlncf( - gwt.dis, - filename=f"{gwtname}.dis.ncf", - ) + ncf = flopy.mf6.ModflowUtlncf( + gwt.dis, + deflate=3, + shuffle=False, + chunk_time=1, + chunk_z=1, + chunk_y=1, + chunk_x=20, + filename=f"{gwtname}.dis.ncf", + ) oc = flopy.mf6.ModflowGwtoc( gwt, @@ -56,6 +69,20 @@ def check_output(idx, test, export, gridded_input): name = cases[idx] gwtname = "gwt_" + name + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{gwtname}.nc") as ds: + assert ds.data_model == "NETCDF4" + if export == "structured": + cmpr = ds.variables["concentration"].filters() + chnk = ds.variables["concentration"].chunking() + assert chnk == [1, 1, 1, 20] + elif export == "ugrid": + cmpr = ds.variables["concentration_l1"].filters() + chnk = ds.variables["concentration_l1"].chunking() + assert chnk == [1, 5] + assert not cmpr["shuffle"] + assert cmpr["complevel"] == 3 + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{gwtname}.nc" diff --git a/autotest/test_netcdf_gwt_prudic2004t2.py b/autotest/test_netcdf_gwt_prudic2004t2.py index a4772c47c79..90253c97000 100644 --- a/autotest/test_netcdf_gwt_prudic2004t2.py +++ b/autotest/test_netcdf_gwt_prudic2004t2.py @@ -12,6 +12,7 @@ xa = pytest.importorskip("xarray") xu = pytest.importorskip("xugrid") +nc = pytest.importorskip("netCDF4") def build_models(idx, test, export, gridded_input): @@ -47,6 +48,10 @@ def check_output(idx, test, export, gridded_input): name = test.name gwtname = "gwt_" + name + # verify format of generated netcdf file + with nc.Dataset(test.workspace / f"{gwtname}.nc") as ds: + assert ds.data_model == "NETCDF4" + if gridded_input == "netcdf": # re-run the simulation with model netcdf input input_fname = f"{gwtname}.nc" diff --git a/doc/mf6io/mf6ivar/dfn/utl-ncf.dfn b/doc/mf6io/mf6ivar/dfn/utl-ncf.dfn index 68694fda43a..f6dc80d5035 100644 --- a/doc/mf6io/mf6ivar/dfn/utl-ncf.dfn +++ b/doc/mf6io/mf6ivar/dfn/utl-ncf.dfn @@ -4,7 +4,7 @@ # --------------------- utl ncf options --------------------- block options -name ogc_wkt +name wkt type string shape lenbigline reader urword @@ -28,67 +28,45 @@ optional true longname description is the keyword used to turn on the netcdf variable shuffle filter when the deflate option is also set. The shuffle filter has the effect of storing the first byte of all of a variable's values in a chunk contiguously, followed by all the second bytes, etc. This can be an optimization for compression with certain types of data. -block options -name chunk_record -type record chunking chunk_time chunk_face chunk_z chunk_y chunk_x -reader urword -optional true -longname netcdf export chunking record -description netcdf export chunking record - -block options -name chunking -type keyword -in_record true -reader urword -optional false -longname keyword when defining chunking parameters -description is a keyword for providing netcdf export chunk sizes. Chunking can dramatically impact data access times and optimal chunking is highly dependent on access patterns (timeseries vs spatial, for example). It can also significantly impact compressibility of the data. A valid input record specifies chunk\_time and either chunk\_face (MESH) or chunk\_z, chunk\_y, and chunk\_x (STRUCTURED). - block options name chunk_time type integer -in_record true reader urword -optional false +optional true longname chunking parameter for the time dimension -description is the keyword used to provide a netcdf export time dimension chunk size. +description is the keyword used to provide a data chunk size for the time dimension in a NETCDF\_MESH2D or NETCDF\_STRUCTURED output file. Must be used in combination with the the chunk\_face parameter (NETCDF\_MESH2D) or the chunk\_z, chunk\_y, and chunk\_x parameter set (NETCDF\_STRUCTURED) to have an effect. block options name chunk_face type integer -in_record true reader urword optional true longname chunking parameter for the mesh face dimension -description is the keyword used to provide a mesh face dimension chunk size. +description is the keyword used to provide a data chunk size for the face dimension in a NETCDF\_MESH2D output file. Must be used in combination with the the chunk\_time parameter to have an effect. block options name chunk_z type integer -in_record true reader urword optional true longname chunking parameter for structured z -description is the keyword used to provide a structured grid z dimensions chunk size. +description is the keyword used to provide a data chunk size for the z dimension in a NETCDF\_STRUCTURED output file. Must be used in combination with the the chunk\_time, chunk\_x and chunk\_y parameter set to have an effect. block options name chunk_y type integer -in_record true reader urword optional true longname chunking parameter for structured y -description is the keyword used to provide a structured grid y dimensions chunk size. +description is the keyword used to provide a data chunk size for the y dimension in a NETCDF\_STRUCTURED output file. Must be used in combination with the the chunk\_time, chunk\_x and chunk\_z parameter set to have an effect. block options name chunk_x type integer -in_record true reader urword optional true longname chunking parameter for structured x -description is the keyword used to provide a structured grid x dimensions chunk size. +description is the keyword used to provide a data chunk size for the x dimension in a NETCDF\_STRUCTURED output file. Must be used in combination with the the chunk\_time, chunk\_y and chunk\_z parameter set to have an effect. block options name modflow6_attr_off @@ -107,12 +85,12 @@ type integer optional true reader urword longname number of cells in layer -description is the number of cells in a in a projected plane layer. +description is the number of cells in a projected plane layer. # --------------------- utl ncf griddata --------------------- block griddata -name lat +name latitude type double precision shape (ncpl) optional true @@ -121,7 +99,7 @@ longname cell center latitude description cell center latitude. block griddata -name lon +name longitude type double precision shape (ncpl) optional true diff --git a/doc/mf6io/mf6ivar/examples/utl-ncf-example.dat b/doc/mf6io/mf6ivar/examples/utl-ncf-example.dat index 926069d7e9c..d23a56f68f1 100644 --- a/doc/mf6io/mf6ivar/examples/utl-ncf-example.dat +++ b/doc/mf6io/mf6ivar/examples/utl-ncf-example.dat @@ -1,8 +1,9 @@ BEGIN OPTIONS DEFLATE 5 SHUFFLE - CHUNKING UGC_TIME 10 UGC_FACE 19690 - OGC_WKT 'PROJCS["NAD_1983_UTM_Zone_14N",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID[" + CHUNK_TIME 10 + CHUNK_FACE 19690 + WKT 'PROJCS["NAD_1983_UTM_Zone_14N",GEOGCS["GCS_North_American_1983",DATUM["D_North_American_1983",SPHEROID[" GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transver se_Mercator"],PARAMETER["False_Easting",500000.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",-99. 0],PARAMETER["Scale_Factor",0.9996],PARAMETER["Latitude_Of_Origin",0.0],UNIT["Meter",1.0]]' diff --git a/src/Idm/utl-ncfidm.f90 b/src/Idm/utl-ncfidm.f90 index dd82d5a1a96..3df61da5b00 100644 --- a/src/Idm/utl-ncfidm.f90 +++ b/src/Idm/utl-ncfidm.f90 @@ -12,11 +12,9 @@ module UtlNcfInputModule public utl_ncf_subpackages type UtlNcfParamFoundType - logical :: ogc_wkt = .false. + logical :: wkt = .false. logical :: deflate = .false. logical :: shuffle = .false. - logical :: chunk_record = .false. - logical :: chunking = .false. logical :: chunk_time = .false. logical :: chunk_face = .false. logical :: chunk_z = .false. @@ -24,8 +22,8 @@ module UtlNcfInputModule logical :: chunk_x = .false. logical :: attr_off = .false. logical :: ncpl = .false. - logical :: lat = .false. - logical :: lon = .false. + logical :: latitude = .false. + logical :: longitude = .false. end type UtlNcfParamFoundType logical :: utl_ncf_multi_package = .false. @@ -37,13 +35,13 @@ module UtlNcfInputModule ] type(InputParamDefinitionType), parameter :: & - utlncf_ogc_wkt = InputParamDefinitionType & + utlncf_wkt = InputParamDefinitionType & ( & 'UTL', & ! component 'NCF', & ! subcomponent 'OPTIONS', & ! block - 'OGC_WKT', & ! tag name - 'OGC_WKT', & ! fortran variable + 'WKT', & ! tag name + 'WKT', & ! fortran variable 'STRING', & ! type 'LENBIGLINE', & ! shape 'CRS well-known text (WKT) string', & ! longname @@ -90,42 +88,6 @@ module UtlNcfInputModule .false. & ! timeseries ) - type(InputParamDefinitionType), parameter :: & - utlncf_chunk_record = InputParamDefinitionType & - ( & - 'UTL', & ! component - 'NCF', & ! subcomponent - 'OPTIONS', & ! block - 'CHUNK_RECORD', & ! tag name - 'CHUNK_RECORD', & ! fortran variable - 'RECORD CHUNKING CHUNK_TIME CHUNK_FACE CHUNK_Z CHUNK_Y CHUNK_X', & ! type - '', & ! shape - 'netcdf export chunking record', & ! longname - .false., & ! required - .false., & ! multi-record - .false., & ! preserve case - .false., & ! layered - .false. & ! timeseries - ) - - type(InputParamDefinitionType), parameter :: & - utlncf_chunking = InputParamDefinitionType & - ( & - 'UTL', & ! component - 'NCF', & ! subcomponent - 'OPTIONS', & ! block - 'CHUNKING', & ! tag name - 'CHUNKING', & ! fortran variable - 'KEYWORD', & ! type - '', & ! shape - 'keyword when defining chunking parameters', & ! longname - .true., & ! required - .true., & ! multi-record - .false., & ! preserve case - .false., & ! layered - .false. & ! timeseries - ) - type(InputParamDefinitionType), parameter :: & utlncf_chunk_time = InputParamDefinitionType & ( & @@ -253,13 +215,13 @@ module UtlNcfInputModule ) type(InputParamDefinitionType), parameter :: & - utlncf_lat = InputParamDefinitionType & + utlncf_latitude = InputParamDefinitionType & ( & 'UTL', & ! component 'NCF', & ! subcomponent 'GRIDDATA', & ! block - 'LAT', & ! tag name - 'LAT', & ! fortran variable + 'LATITUDE', & ! tag name + 'LATITUDE', & ! fortran variable 'DOUBLE1D', & ! type 'NCPL', & ! shape 'cell center latitude', & ! longname @@ -271,13 +233,13 @@ module UtlNcfInputModule ) type(InputParamDefinitionType), parameter :: & - utlncf_lon = InputParamDefinitionType & + utlncf_longitude = InputParamDefinitionType & ( & 'UTL', & ! component 'NCF', & ! subcomponent 'GRIDDATA', & ! block - 'LON', & ! tag name - 'LON', & ! fortran variable + 'LONGITUDE', & ! tag name + 'LONGITUDE', & ! fortran variable 'DOUBLE1D', & ! type 'NCPL', & ! shape 'cell center longitude', & ! longname @@ -291,11 +253,9 @@ module UtlNcfInputModule type(InputParamDefinitionType), parameter :: & utl_ncf_param_definitions(*) = & [ & - utlncf_ogc_wkt, & + utlncf_wkt, & utlncf_deflate, & utlncf_shuffle, & - utlncf_chunk_record, & - utlncf_chunking, & utlncf_chunk_time, & utlncf_chunk_face, & utlncf_chunk_z, & @@ -303,8 +263,8 @@ module UtlNcfInputModule utlncf_chunk_x, & utlncf_attr_off, & utlncf_ncpl, & - utlncf_lat, & - utlncf_lon & + utlncf_latitude, & + utlncf_longitude & ] type(InputParamDefinitionType), parameter :: & diff --git a/src/Utilities/Export/DisNCStructured.f90 b/src/Utilities/Export/DisNCStructured.f90 index 235b8a4fa57..823cfe36661 100644 --- a/src/Utilities/Export/DisNCStructured.f90 +++ b/src/Utilities/Export/DisNCStructured.f90 @@ -10,8 +10,8 @@ module DisNCStructuredModule use KindModule, only: DP, I4B, LGP use ConstantsModule, only: LINELENGTH, LENBIGLINE, LENCOMPONENTNAME, & LENMEMPATH, LENVARNAME - use SimVariablesModule, only: errmsg - use SimModule, only: store_error, store_error_filename + use SimVariablesModule, only: errmsg, warnmsg + use SimModule, only: store_error, store_warning, store_error_filename use MemoryManagerModule, only: mem_setptr use InputDefinitionModule, only: InputParamDefinitionType use CharacterStringModule, only: CharacterStringType @@ -43,8 +43,8 @@ module DisNCStructuredModule integer(I4B) :: x_bnds !< x boundaries 2D array integer(I4B) :: y_bnds !< y boundaries 2D array integer(I4B) :: z_bnds !< z boundaries 2D array - integer(I4B) :: lat !< latitude 2D array - integer(I4B) :: lon !< longitude 2D array + integer(I4B) :: latitude !< latitude 2D array + integer(I4B) :: longitude !< longitude 2D array contains end type StructuredNCVarIdType @@ -53,8 +53,8 @@ module DisNCStructuredModule type(StructuredNCVarIdType) :: var_ids !< structured variable ids type type(DisType), pointer :: dis => null() !< pointer to model dis package integer(I4B) :: nlay !< number of layers - real(DP), dimension(:), pointer, contiguous :: lat => null() !< lat input array pointer - real(DP), dimension(:), pointer, contiguous :: lon => null() !< lon input array pointer + real(DP), dimension(:), pointer, contiguous :: latitude => null() !< lat input array pointer + real(DP), dimension(:), pointer, contiguous :: longitude => null() !< lon input array pointer integer(I4B), pointer :: chunk_z !< chunking parameter for z dimension integer(I4B), pointer :: chunk_y !< chunking parameter for y dimension integer(I4B), pointer :: chunk_x !< chunking parameter for x dimension @@ -136,21 +136,31 @@ subroutine dis_export_init(this, modelname, modeltype, modelfname, nc_fname, & if (this%chunk_time > 0 .and. this%chunk_z > 0 .and. & this%chunk_y > 0 .and. this%chunk_x > 0) then this%chunking_active = .true. + else if (this%chunk_time > 0 .or. this%chunk_z > 0 .or. & + this%chunk_y > 0 .or. this%chunk_x > 0) then + this%chunk_time = -1 + this%chunk_z = -1 + this%chunk_y = -1 + this%chunk_x = -1 + write (warnmsg, '(a)') 'Ignoring user provided NetCDF chunking & + ¶meters. Define chunk_time, chunk_x, chunk_y and chunk_z input & + ¶meters to see an effect.' + call store_warning(warnmsg) end if - call get_isize('LAT', this%ncf_mempath, latsz) - call get_isize('LON', this%ncf_mempath, lonsz) + call get_isize('LATITUDE', this%ncf_mempath, latsz) + call get_isize('LONGITUDE', this%ncf_mempath, lonsz) if (latsz > 0 .and. lonsz > 0) then this%latlon = .true. - call mem_setptr(this%lat, 'LAT', this%ncf_mempath) - call mem_setptr(this%lon, 'LON', this%ncf_mempath) + call mem_setptr(this%latitude, 'LATITUDE', this%ncf_mempath) + call mem_setptr(this%longitude, 'LONGITUDE', this%ncf_mempath) end if end if ! create the netcdf file call nf_verify(nf90_create(this%nc_fname, & - IAND(NF90_CLOBBER, NF90_NETCDF4), this%ncid), & + IOR(NF90_CLOBBER, NF90_NETCDF4), this%ncid), & this%nc_fname) end subroutine dis_export_init @@ -705,7 +715,7 @@ subroutine define_dim(this) 'projection_y_coordinate'), this%nc_fname) call nf_verify(nf90_put_att(this%ncid, this%var_ids%y, 'long_name', & 'Northing'), this%nc_fname) - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then call nf_verify(nf90_put_att(this%ncid, this%var_ids%y, 'grid_mapping', & this%gridmap_name), this%nc_fname) end if @@ -728,7 +738,7 @@ subroutine define_dim(this) 'projection_x_coordinate'), this%nc_fname) call nf_verify(nf90_put_att(this%ncid, this%var_ids%x, 'long_name', & 'Easting'), this%nc_fname) - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then call nf_verify(nf90_put_att(this%ncid, this%var_ids%x, 'grid_mapping', & this%gridmap_name), this%nc_fname) end if @@ -791,12 +801,12 @@ end subroutine define_dependent subroutine define_gridmap(this) class(DisNCStructuredType), intent(inout) :: this integer(I4B) :: var_id - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then call nf_verify(nf90_redef(this%ncid), this%nc_fname) call nf_verify(nf90_def_var(this%ncid, this%gridmap_name, NF90_INT, & var_id), this%nc_fname) ! TODO: consider variants epsg_code, spatial_ref, esri_pe_string, wkt, etc - call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%ogc_wkt), & + call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%wkt), & this%nc_fname) call nf_verify(nf90_enddef(this%ncid), this%nc_fname) call nf_verify(nf90_put_var(this%ncid, var_id, 1), & @@ -808,28 +818,28 @@ end subroutine define_gridmap !< subroutine define_projection(this) class(DisNCStructuredType), intent(inout) :: this - if (this%latlon .and. this%ogc_wkt /= '') then + if (this%latlon .and. this%wkt /= '') then ! lat call nf_verify(nf90_def_var(this%ncid, 'lat', NF90_DOUBLE, & (/this%dim_ids%x, this%dim_ids%y/), & - this%var_ids%lat), this%nc_fname) - call nf_verify(nf90_put_att(this%ncid, this%var_ids%lat, 'units', & - 'degrees_north'), this%nc_fname) - call nf_verify(nf90_put_att(this%ncid, this%var_ids%lat, 'standard_name', & - 'latitude'), this%nc_fname) - call nf_verify(nf90_put_att(this%ncid, this%var_ids%lat, 'long_name', & - 'latitude'), this%nc_fname) + this%var_ids%latitude), this%nc_fname) + call nf_verify(nf90_put_att(this%ncid, this%var_ids%latitude, & + 'units', 'degrees_north'), this%nc_fname) + call nf_verify(nf90_put_att(this%ncid, this%var_ids%latitude, & + 'standard_name', 'latitude'), this%nc_fname) + call nf_verify(nf90_put_att(this%ncid, this%var_ids%latitude, & + 'long_name', 'latitude'), this%nc_fname) ! lon call nf_verify(nf90_def_var(this%ncid, 'lon', NF90_DOUBLE, & (/this%dim_ids%x, this%dim_ids%y/), & - this%var_ids%lon), this%nc_fname) - call nf_verify(nf90_put_att(this%ncid, this%var_ids%lon, 'units', & - 'degrees_east'), this%nc_fname) - call nf_verify(nf90_put_att(this%ncid, this%var_ids%lon, 'standard_name', & - 'longitude'), this%nc_fname) - call nf_verify(nf90_put_att(this%ncid, this%var_ids%lon, 'long_name', & - 'longitude'), this%nc_fname) + this%var_ids%longitude), this%nc_fname) + call nf_verify(nf90_put_att(this%ncid, this%var_ids%longitude, & + 'units', 'degrees_east'), this%nc_fname) + call nf_verify(nf90_put_att(this%ncid, this%var_ids%longitude, & + 'standard_name', 'longitude'), this%nc_fname) + call nf_verify(nf90_put_att(this%ncid, this%var_ids%longitude, & + 'long_name', 'longitude'), this%nc_fname) end if end subroutine define_projection @@ -837,16 +847,16 @@ end subroutine define_projection !< subroutine add_proj_data(this) class(DisNCStructuredType), intent(inout) :: this - if (this%latlon .and. this%ogc_wkt /= '') then + if (this%latlon .and. this%wkt /= '') then ! lat - call nf_verify(nf90_put_var(this%ncid, this%var_ids%lat, & - this%lat, start=(/1, 1/), & + call nf_verify(nf90_put_var(this%ncid, this%var_ids%latitude, & + this%latitude, start=(/1, 1/), & count=(/this%dis%ncol, this%dis%nrow/)), & this%nc_fname) ! lon - call nf_verify(nf90_put_var(this%ncid, this%var_ids%lon, & - this%lon, start=(/1, 1/), & + call nf_verify(nf90_put_var(this%ncid, this%var_ids%longitude, & + this%longitude, start=(/1, 1/), & count=(/this%dis%ncol, this%dis%nrow/)), & this%nc_fname) end if diff --git a/src/Utilities/Export/MeshNCModel.f90 b/src/Utilities/Export/MeshNCModel.f90 index a773bebfd80..0062e8e7fc6 100644 --- a/src/Utilities/Export/MeshNCModel.f90 +++ b/src/Utilities/Export/MeshNCModel.f90 @@ -9,8 +9,8 @@ module MeshModelModule use KindModule, only: DP, I4B, LGP use ConstantsModule, only: LINELENGTH, LENCOMPONENTNAME, LENMEMPATH, & DNODATA, DHNOFLO - use SimVariablesModule, only: errmsg - use SimModule, only: store_error, store_error_filename + use SimVariablesModule, only: errmsg, warnmsg + use SimModule, only: store_error, store_warning, store_error_filename use MemoryManagerModule, only: mem_setptr use InputDefinitionModule, only: InputParamDefinitionType use CharacterStringModule, only: CharacterStringType @@ -123,11 +123,17 @@ subroutine mesh_init(this, modelname, modeltype, modelfname, nc_fname, & if (this%chunk_time > 0 .and. this%chunk_face > 0) then this%chunking_active = .true. + else if (this%chunk_time > 0 .or. this%chunk_face > 0) then + this%chunk_face = -1 + this%chunk_time = -1 + write (warnmsg, '(a)') 'Ignoring user provided NetCDF chunking parameter. & + &Define chunk_time and chunk_face input parameters to see an effect.' + call store_warning(warnmsg) end if ! create the netcdf file call nf_verify(nf90_create(this%nc_fname, & - IAND(NF90_CLOBBER, NF90_NETCDF4), this%ncid), & + IOR(NF90_CLOBBER, NF90_NETCDF4), this%ncid), & this%nc_fname) end subroutine mesh_init @@ -319,16 +325,16 @@ subroutine define_gridmap(this) integer(I4B) :: var_id ! was projection info provided - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then ! create projection variable call nf_verify(nf90_redef(this%ncid), this%nc_fname) call nf_verify(nf90_def_var(this%ncid, this%gridmap_name, NF90_INT, & var_id), this%nc_fname) ! cf-conventions prefers 'crs_wkt' - !call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%ogc_wkt), & + !call nf_verify(nf90_put_att(this%ncid, var_id, 'crs_wkt', this%wkt), & ! this%nc_fname) ! QGIS recognizes 'wkt' - call nf_verify(nf90_put_att(this%ncid, var_id, 'wkt', this%ogc_wkt), & + call nf_verify(nf90_put_att(this%ncid, var_id, 'wkt', this%wkt), & this%nc_fname) call nf_verify(nf90_enddef(this%ncid), this%nc_fname) call nf_verify(nf90_put_var(this%ncid, var_id, 1), & @@ -378,7 +384,7 @@ subroutine create_mesh(this) call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, & 'long_name', 'Easting'), this%nc_fname) - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then ! associate with projection call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_x, & 'grid_mapping', this%gridmap_name), & @@ -399,7 +405,7 @@ subroutine create_mesh(this) call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, & 'long_name', 'Northing'), this%nc_fname) - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then ! associate with projection call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_node_y, & 'grid_mapping', this%gridmap_name), & @@ -421,7 +427,7 @@ subroutine create_mesh(this) 'long_name', 'Easting'), this%nc_fname) call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, 'bounds', & 'mesh_face_xbnds'), this%nc_fname) - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then ! associate with projection call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_x, & 'grid_mapping', this%gridmap_name), & @@ -451,7 +457,7 @@ subroutine create_mesh(this) call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, 'bounds', & 'mesh_face_ybnds'), this%nc_fname) - if (this%ogc_wkt /= '') then + if (this%wkt /= '') then ! associate with projection call nf_verify(nf90_put_att(this%ncid, this%var_ids%mesh_face_y, & 'grid_mapping', this%gridmap_name), & diff --git a/src/Utilities/Export/ModelExport.f90 b/src/Utilities/Export/ModelExport.f90 index 327ed135de2..1d88a4c1f0c 100644 --- a/src/Utilities/Export/ModelExport.f90 +++ b/src/Utilities/Export/ModelExport.f90 @@ -8,6 +8,8 @@ module ModelExportModule use KindModule, only: DP, I4B, LGP + use SimModule, only: store_error, store_error_filename + use SimVariablesModule, only: errmsg use ConstantsModule, only: LINELENGTH, LENMODELNAME, LENCOMPONENTNAME, & LENMEMPATH use ListModule, only: ListType @@ -79,10 +81,11 @@ subroutine modelexports_create(iout) use InputLoadTypeModule, only: GetDynamicModelFromList use SimVariablesModule, only: idm_context use NCModelExportModule, only: NETCDF_MESH2D, NETCDF_STRUCTURED + use SourceCommonModule, only: file_ext integer(I4B), intent(in) :: iout type(ModelDynamicPkgsType), pointer :: model_dynamic_input type(ExportModelType), pointer :: export_model - character(len=LENMEMPATH) :: modelnam_mempath, model_mempath + character(len=LENMEMPATH) :: modelnam_mempath, model_mempath, ext integer(I4B), pointer :: disenum integer(I4B) :: n logical(LGP) :: found @@ -119,6 +122,16 @@ subroutine modelexports_create(iout) end if end if + if (found) then + ext = file_ext(export_model%nc_fname) + if (ext /= 'nc') then + errmsg = 'NetCDF output file name must use ".nc" extension. '// & + 'Filename="'//trim(export_model%nc_fname)//'".' + call store_error(errmsg) + call store_error_filename(export_model%modelfname) + end if + end if + ! add model to list call add_export_model(export_model) end do diff --git a/src/Utilities/Export/NCModel.f90 b/src/Utilities/Export/NCModel.f90 index 6cad5115f46..cb4d83d28c7 100644 --- a/src/Utilities/Export/NCModel.f90 +++ b/src/Utilities/Export/NCModel.f90 @@ -74,7 +74,7 @@ module NCModelExportModule character(len=LINELENGTH) :: mesh_name = 'mesh' !< name of mesh container variable character(len=LENMEMPATH) :: dis_mempath !< discretization input mempath character(len=LENMEMPATH) :: ncf_mempath !< netcdf utility package input mempath - character(len=LENBIGLINE) :: ogc_wkt !< wkt user string + character(len=LENBIGLINE) :: wkt !< wkt user string character(len=LINELENGTH) :: datetime !< export file creation time character(len=LINELENGTH) :: xname !< dependent variable name type(NCExportAnnotation) :: annotation !< export file annotation @@ -293,7 +293,7 @@ subroutine export_init(this, modelname, modeltype, modelfname, nc_fname, & this%nc_fname = nc_fname this%gridmap_name = '' this%ncf_mempath = '' - this%ogc_wkt = '' + this%wkt = '' this%datetime = '' this%xname = '' this%disenum = disenum @@ -342,8 +342,8 @@ subroutine export_init(this, modelname, modeltype, modelfname, nc_fname, & found_mempath) if (found_mempath) then - call mem_set_value(this%ogc_wkt, 'OGC_WKT', this%ncf_mempath, & - ncf_found%ogc_wkt) + call mem_set_value(this%wkt, 'WKT', this%ncf_mempath, & + ncf_found%wkt) call mem_set_value(this%deflate, 'DEFLATE', this%ncf_mempath, & ncf_found%deflate) call mem_set_value(this%shuffle, 'SHUFFLE', this%ncf_mempath, & @@ -354,7 +354,7 @@ subroutine export_init(this, modelname, modeltype, modelfname, nc_fname, & ncf_found%chunk_time) end if - if (ncf_found%ogc_wkt) then + if (ncf_found%wkt) then this%gridmap_name = 'projection' end if @@ -364,12 +364,11 @@ subroutine export_init(this, modelname, modeltype, modelfname, nc_fname, & end if ! set datetime string - if (isim_mode /= MVALIDATE .and. datetime0 == '') then - errmsg = 'TDIS parameter START_DATE_TIME required for NetCDF export.' - call store_error(errmsg) - call store_error_filename(modelfname) - else + if (datetime0 /= '') then this%datetime = 'days since '//trim(datetime0) + else + ! January 1, 1970 at 00:00:00 UTC + this%datetime = 'days since 1970-01-01T00:00:00' end if ! set total nstp