convert transient array data in UZF package to Transient2d

autotest/t007_test.py

@@ -1,4 +1,6 @@
 # Test export module
+import sys
+sys.path.append('..')
 import copy
 import glob
 import os
@@ -138,6 +140,49 @@ def test_export_output():
     arr_mask = arr.mask[0]
     assert np.array_equal(ibound_mask, arr_mask)
 
+def test_write_shapefile():
+    from flopy.utils.reference import SpatialReference
+    from flopy.export.shapefile_utils import shp2recarray
+    from flopy.export.shapefile_utils import write_grid_shapefile, write_grid_shapefile2
+
+    sr = SpatialReference(delr=np.ones(10) *1.1,  # cell spacing along model rows
+                          delc=np.ones(10) *1.1,  # cell spacing along model columns
+                          epsg=26715,
+                          lenuni=1  # MODFLOW length units
+                          )
+    outshp1 = os.path.join(tpth, 'junk.shp')
+    outshp2 = os.path.join(tpth, 'junk2.shp')
+    write_grid_shapefile(outshp1, sr, array_dict={})
+    write_grid_shapefile2(outshp2, sr, array_dict={})
+
+    for outshp in [outshp1, outshp2]:
+        # check that pyshp reads integers
+        # this only check that row/column were recorded as "N"
+        # not how they will be cast by python or numpy
+        import shapefile as sf
+        sfobj = sf.Reader(outshp)
+        for f in sfobj.fields:
+            if f[0] == 'row' or f[0] == 'column':
+                assert f[1] == 'N'
+        recs = list(sfobj.records())
+        for r in recs[0]:
+            assert isinstance(r, int)
+
+        # check that row and column appear as integers in recarray
+        ra = shp2recarray(outshp)
+        assert np.issubdtype(ra.dtype['row'], np.integer)
+        assert np.issubdtype(ra.dtype['column'], np.integer)
+
+        try: # check that fiona reads integers
+            import fiona
+            with fiona.open(outshp) as src:
+                meta = src.meta
+                assert 'int' in meta['schema']['properties']['row']
+                assert 'int' in meta['schema']['properties']['column']
+        except:
+            pass
+
+
 def test_export_array():
 
     try:
@@ -847,7 +892,8 @@ def build_sfr_netcdf():
     # for namfile in ["fhb.nam"]:
     # export_netcdf(namfile)
     #test_freyberg_export()
-    test_export_array()
+    #test_export_array()
+    test_write_shapefile()
     #test_wkt_parse()
     #test_get_rc_from_node_coordinates()
     pass

autotest/t009_test.py

@@ -49,14 +49,34 @@
              }
 
 def create_sfr_data():
-
-    dtype = np.dtype([('iseg', int), ('ireach', int)])
+    dtype = np.dtype([('k', int),
+             ('i', int),
+             ('j', int),
+             ('iseg', int),
+             ('ireach', int)])
     r = create_empty_recarray(27, dtype=dtype)
+    r['i'] = [3, 4, 5,
+              7, 8, 9,
+              0, 1, 2,
+              4, 4, 5,
+              0, 0, 0,
+              3, 4, 5,
+              0, 1, 2,
+              4, 5, 6,
+              2, 2, 2]
+    r['j'] = [0, 1, 2,
+              6, 6, 6,
+              6, 6, 6,
+              3, 4, 5,
+              9, 8, 7,
+              6, 6, 6,
+              0, 0, 0,
+              6, 6, 6,
+              9, 8, 7]
     r['iseg'] = sorted(list(range(1, 10)) * 3)
     r['ireach'] = [1, 2, 3] * 9
 
-    dtype = np.dtype([('nseg', int), ('outseg', int)])
-    d = create_empty_recarray(9, dtype=dtype)
+    d = create_empty_recarray(9, dtype=np.dtype([('nseg', int), ('outseg', int)]))
     d['nseg'] = range(1, 10)
     d['outseg'] = [4, 0, 6, 8, 3, 8, 1, 2, 8]
     return r, d
@@ -251,6 +271,38 @@ def test_const():
     assert sfr.const == 1.486 * 86400.
     assert True
 
+def test_export():
+    fm = flopy.modflow
+    m = fm.Modflow()
+    dis = fm.ModflowDis(m, 1, 10, 10, lenuni=2, itmuni=4)
+    m.sr = SpatialReference(delr=m.dis.delr.array, delc=m.dis.delc.array)
+    m.sr.write_shapefile(os.path.join(outpath, 'grid.shp'))
+    r, d = create_sfr_data()
+    sfr = flopy.modflow.ModflowSfr2(m, reach_data=r, segment_data={0: d})
+    sfr.segment_data[0]['flow'][-1] = 1e4
+    sfr.stress_period_data.export(os.path.join(outpath, 'sfr.shp'), sparse=True)
+    sfr.export_linkages(os.path.join(outpath, 'linkages.shp'))
+    sfr.export_outlets(os.path.join(outpath, 'outlets.shp'))
+    sfr.export_transient_variable(os.path.join(outpath, 'inlets.shp'),
+                                  'flow')
+
+    from flopy.export.shapefile_utils import shp2recarray
+    ra = shp2recarray(os.path.join(outpath, 'inlets.shp'))
+    assert ra.flow0[0] == 1e4
+    ra = shp2recarray(os.path.join(outpath, 'outlets.shp'))
+    assert ra.iseg[0] + ra.ireach[0] == 5
+    ra = shp2recarray(os.path.join(outpath, 'linkages.shp'))
+    crds = np.array(list(ra.geometry[2].coords))
+    assert np.array_equal(crds, np.array([[2.5,  4.5], [3.5,  5.5]]))
+    ra = shp2recarray(os.path.join(outpath, 'sfr.shp'))
+    assert ra.iseg0.sum() == sfr.reach_data.iseg.sum()
+    assert ra.ireach0.sum() == sfr.reach_data.ireach.sum()
+    y = np.concatenate([np.array(g.exterior)[:, 1] for g in ra.geometry])
+    x = np.concatenate([np.array(g.exterior)[:, 0] for g in ra.geometry])
+    assert (x.min(), y.min(), x.max(), y.max()) == m.sr.bounds
+    assert ra[(ra.iseg0 == 2) & (ra.ireach0 == 1)]['geometry'][0].bounds \
+        == (6.0, 2.0, 7.0, 3.0)
+
 def test_example():
     m = flopy.modflow.Modflow.load('test1ss.nam', version='mf2005',
                                    exe_name='mf2005.exe',
@@ -397,12 +449,13 @@ def test_sfr_plot():
     pass
 
 if __name__ == '__main__':
-    #test_sfr()
+    test_sfr()
     test_sfr_renumbering()
-    #test_example()
-    #test_transient_example()
-    #test_sfr_plot()
-    #test_assign_layers()
-    #test_SfrFile()
-    #test_const()
+    test_example()
+    test_export()
+    test_transient_example()
+    test_sfr_plot()
+    test_assign_layers()
+    test_SfrFile()
+    test_const()
     pass

autotest/t024_test.py

@@ -104,7 +104,8 @@ def test_properties_check():
 
 
 if __name__ == '__main__':
-    # for mfnam in testmodels:
-    #    checker_on_load(mfnam)
+    print('numpy version: {}'.format(np.__version__))
+    for mfnam in testmodels:
+        checker_on_load(mfnam)
     test_bcs_check()
-    #test_properties_check()
+    test_properties_check()

autotest/t031_test.py

@@ -22,6 +22,11 @@
 for f in mffiles:
     shutil.copy(f, os.path.join(path, os.path.split(f)[1]))
 
+def ra_slice(ra, cols):
+    raslice = np.column_stack([ra[c] for c in cols])
+    dtype = [(str(d[0]), d[1]) for d in ra.dtype.descr if d[0] in cols]
+    return np.array([tuple(r) for r in raslice],
+                    dtype=dtype).view(np.recarray)
 
 def test_mpsim():
     model_ws = path

autotest/t034_test.py

@@ -45,10 +45,33 @@ def test_create():
 
     finf = np.loadtxt(os.path.join(datpth, 'finf.dat'))
     finf = np.reshape(finf, (m.nper, m.nrow, m.ncol))
+    finf = {i: finf[i] for i in range(finf.shape[0])}
 
     extwc = np.loadtxt(os.path.join(datpth, 'extwc.dat'))
 
-    uzgag = {-68: [-68],
+    uzgag = {-68: [],
+              65: [3, 6, 1],
+              66: [6, 3, 2],
+              67: [10, 5, 3]}
+    uzf = flopy.modflow.ModflowUzf1(m,
+                                    nuztop=1, iuzfopt=1, irunflg=1, ietflg=1,
+                                    ipakcb=0,
+                                    iuzfcb2=61,
+                                    # binary output of recharge and groundwater discharge
+                                    ntrail2=25, nsets=20,
+                                    surfdep=1.0, uzgag=uzgag,
+                                    iuzfbnd=m.bas6.ibound.array,
+                                    irunbnd=irunbnd,
+                                    vks=vks,
+                                    finf=finf,
+                                    eps=3.5,
+                                    thts=0.3,
+                                    pet=5.000000E-08,
+                                    extdp=15.,
+                                    extwc=extwc
+                                    )
+    assert uzf.uzgag == uzgag
+    uzgag2 = {-68: [-68],
              65: [3, 6, 65, 1],
              66: [6, 3, 66, 2],
              67: [10, 5, 67, 3]}
@@ -58,8 +81,8 @@ def test_create():
                                     ipakcb=0,
                                     iuzfcb2=61,
                                     # binary output of recharge and groundwater discharge
-                                    ntrail2=25, nsets=20, nuzgag=4,
-                                    surfdep=1.0, uzgag=uzgag,
+                                    ntrail2=25, nsets=20,
+                                    surfdep=1.0, uzgag=uzgag2,
                                     iuzfbnd=m.bas6.ibound.array,
                                     irunbnd=irunbnd,
                                     vks=vks,
@@ -70,13 +93,36 @@ def test_create():
                                     extdp=15.,
                                     extwc=extwc
                                     )
+    assert uzf.uzgag == uzgag
+    uzgaglist = [[-68],
+                 [3, 6, 65, 1],
+                 [6, 3, 66, 2],
+                 [10, 5, 67, 3]]
+    uzf = flopy.modflow.ModflowUzf1(m,
+                                    nuztop=1, iuzfopt=1, irunflg=1, ietflg=1,
+                                    ipakcb=0,
+                                    iuzfcb2=61,
+                                    # binary output of recharge and groundwater discharge
+                                    ntrail2=25, nsets=20,
+                                    surfdep=1.0, uzgag=uzgaglist,
+                                    iuzfbnd=m.bas6.ibound.array,
+                                    irunbnd=irunbnd,
+                                    vks=vks,
+                                    finf=finf,
+                                    eps=3.5,
+                                    thts=0.3,
+                                    pet=5.000000E-08,
+                                    extdp=15.,
+                                    extwc=extwc
+                                    )
+    assert uzf.uzgag == uzgag
     m.write_input()
-
+    uzf2 = flopy.modflow.ModflowUzf1.load(os.path.join(cpth, uzf.file_name[0]), m)
+    assert uzf2.uzgag == uzgag
     m2 = flopy.modflow.Modflow.load('UZFtest2.nam', version='mf2005',
                                     exe_name='mf2005',
                                     verbose=True,
-                                    model_ws=os.path.split(gpth)[0])
-
+                                    model_ws=os.path.split(gpth)[0], forgive=False)
     # verify that all of the arrays in the created UZF package are the same
     # as those in the loaded example
     attrs = [attr for attr in dir(uzf)

autotest/t049_test.py

@@ -133,6 +133,8 @@ def test_pathline_plot():
         run = False
         lpth = pth
 
+    nampath = os.path.join(lpth, mfnam)
+    assert os.path.exists(nampath), "namefile {} doesn't exist.".format(nampath)
     # load the modflow files for model map
     m = flopy.modflow.Modflow.load(mfnam, model_ws=lpth, verbose=True, forgive=False,
                                    exe_name=mf2005_exe)

docs/version_changes.md

@@ -7,6 +7,8 @@ FloPy Changes
 * MODFLOW 6 identifiers are now zero based
 * Added remove_package method in MFSimulation and MFModel that removes MODFLOW 6 packages from the existing simulation/model
 * Flopy code for MODFLOW 6 generalized to support different model types
+* Added `ModflowSfr2.export_transient_variable` method to export shapefiles of segment data variables, with stress period data as attributes
+* Added support for UZF package gages
 
 * Bug fixes:
     * Fixed issue with default settings for MODFLOW 5 SUB package `dp` dataset.
@@ -17,6 +19,10 @@ FloPy Changes
 	* IMS file name conflicts now automatically resolved
 	* Fixed issue with passing in numpy ndarrays arrays as layered data
 	* Doc string formatting for MODFLOW 6 packages fixed to make doc strings easier to read
+	* UZF package: fixed issues with handling of finf, pet, extdp and extwc arrays.
+	* SFR package: fixed issue with reading stress period data where not all segments are listed for periods > 0.
+	* `SpatialReference.write_gridSpec` was not converting the model origin coordinates to model length units.
+	* shorted integer field lengths written to shapefiles to 18 characters; some readers may misinterpret longer field lengths as float dtypes.
 
 ### Version 3.2.8
 * Added `has_package(name)` method to see if a package exists. This feature goes nicely with `get_package(name)` method.
@@ -25,6 +31,11 @@ FloPy Changes
     * Installation: Added dfn files required by MODFLOW 6 functionality to MANIFEST.in so that they are included in the distribution.
     * SFR2 package: Fixed issue reading transient data when `ISFOPT` is 4 or 5 for the first stress period.
 
+### Version 3.2.7 - develop
+* Added SFR2 package functionarlity
+	* `export_inlets` method to write shapefile showing locations where external flows are entering the stream network.  
+
+
 ### Version 3.2.7
 * Added beta support for MODFLOW 6 See [here](./mf6.md) for more information.
 * Added support for retrieving time series from binary cell-by-cell files. Cell-by-cell time series are accessed in the same way they are accessed for heads and concentrations but a text string is required.