<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 This notebook shows some lower-level functionality in flopy for working with shapefiles
including
Step1:  write a numpy record array to a shapefile
in this case, we want to visualize output from the checker
first make a toy model
Step2:  set the spatial reference
note that the grid spacings have to be converted explicitly (in this case from feet to meters)
Step3:  make geometry objects for the cells with errors
geometry objects allow the shapefile writer to be simpler and agnostic about the kind of geometry
Step4:  write the shapefile
the projection (.prj) file can be written using an epsg code  
or copied from an existing .prj file
Step5:  read it back in
flopy geometry objects representing the shapes are stored in the 'geometry' field
Step6:  How the epsg feature works
requires an internet connection the first time to get the prj text from spatialreference.org using requests  
if it doesn't exist, epsgref.py is created in the site-packages folder  
the prj text is then stashed in this file as an entry in a dictionary prj  
this dictionary is imported and checked for the prjtext each time flopy.utils.reference.getprj is called
Step7:  working with the flopy.utils.reference.epsgRef handler
Step8:  remove an entry
Step9:  start over with a new file
Step10:  Other geometry types
Linestring
create geometry objects for pathlines from a MODPATH simulation
plot the paths using the built in plotting method
Step11:  Points
<ASSISTANT_TASK:>
Python Code:
%matplotlib inline
try:
    from importlib import reload # python3
except:
    pass # python2 (reload in default namespace)
import sys
import shutil
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
import flopy
from flopy.utils.geometry import Polygon, LineString, Point
from flopy.utils.reference import SpatialReference
from flopy.export.shapefile_utils import recarray2shp, shp2recarray
from flopy.utils.modpathfile import PathlineFile, EndpointFile
from flopy.utils.reference import epsgRef
ep = epsgRef()
ep.reset()
print(sys.version)
print('numpy version: {}'.format(np.__version__))
print('matplotlib version: {}'.format(mpl.__version__))
print('flopy version: {}'.format(flopy.__version__))
Explanation: This notebook shows some lower-level functionality in flopy for working with shapefiles
including:
* recarray2shp convience function for writing a numpy record array to a shapefile
* shp2recarray convience function for quickly reading a shapefile into a numpy recarray
* utils.geometry classes for writing shapefiles of model input/output. For example, quickly writing a shapefile of model cells with errors identified by the checker
* demonstration of how the epsgRef class works for retrieving projection file information (WKT text) from spatialreference.org, and caching that information locally for when an internet connection isn't available
* how to reset epsgRef if it becomes corrupted
* examples of how the Point and LineString classes can be used to quickly plot pathlines and endpoints from MODPATH (these are also used by the PathlineFile and EndpointFile classes to write shapefiles of this output)
End of explanation
m = flopy.modflow.Modflow('toy_model', model_ws='temp')
botm = np.zeros((2, 10, 10))
botm[0, :, :] = 1.5
botm[1, 5, 5] = 4 # negative layer thickness!
botm[1, 6, 6] = 4
dis = flopy.modflow.ModflowDis(nrow=10, ncol=10, 
                               nlay=2, delr=100, delc=100,
                               top=3, botm=botm, model=m)
Explanation: write a numpy record array to a shapefile
in this case, we want to visualize output from the checker
first make a toy model
End of explanation
m.sr = SpatialReference(delr=m.dis.delr * .3048, delc=m.dis.delc * .3048, xul=600000, yul=5170000, 
                        proj4_str='EPSG:26715', rotation=45)
chk = dis.check()
chk.summary_array
Explanation: set the spatial reference
note that the grid spacings have to be converted explicitly (in this case from feet to meters)
End of explanation
get_vertices = m.sr.get_vertices # function to get the referenced vertices for a model cell
geoms = [Polygon(get_vertices(i, j)) for i, j in chk.summary_array[['i', 'j']]]
geoms[0].type
geoms[0].exterior
geoms[0].bounds
geoms[0].plot() # this feature requires descartes
Explanation: make geometry objects for the cells with errors
geometry objects allow the shapefile writer to be simpler and agnostic about the kind of geometry
End of explanation
recarray2shp(chk.summary_array, geoms, 'temp/test.shp', epsg=26715)
shutil.copy('temp/test.prj', 'temp/26715.prj')
recarray2shp(chk.summary_array, geoms, 'temp/test.shp', prj='temp/26715.prj')
Explanation: write the shapefile
the projection (.prj) file can be written using an epsg code  
or copied from an existing .prj file
End of explanation
ra = shp2recarray('temp/test.shp')
ra
ra.geometry[0].plot()
Explanation: read it back in
flopy geometry objects representing the shapes are stored in the 'geometry' field
End of explanation
import epsgref
reload(epsgref)
from epsgref import prj
prj
from flopy.utils.reference import getprj, epsgRef
getprj(4326)
reload(epsgref)
from epsgref import prj
for k, v in prj.items():
    print('{}:\n{}\n'.format(k, v))
Explanation: How the epsg feature works
requires an internet connection the first time to get the prj text from spatialreference.org using requests  
if it doesn't exist, epsgref.py is created in the site-packages folder  
the prj text is then stashed in this file as an entry in a dictionary prj  
this dictionary is imported and checked for the prjtext each time flopy.utils.reference.getprj is called
End of explanation
ep = epsgRef()
ep.add(9999, 'junk')
epsgRef.show()
Explanation: working with the flopy.utils.reference.epsgRef handler
End of explanation
ep.remove(9999)
epsgRef.show()
Explanation: remove an entry
End of explanation
ep.reset()
reload(epsgref)
from epsgref import prj
prj
len(prj.keys())
Explanation: start over with a new file
End of explanation
pthfile = PathlineFile('../data/mp6/EXAMPLE-3.pathline')
pthdata = pthfile._data.view(np.recarray)
length_mult = 1. # multiplier to convert coordinates from model to real world
rot = 0 # grid rotation
particles = np.unique(pthdata.particleid)
geoms = []
for pid in particles:
    ra = pthdata[pthdata.particleid == pid]
    x, y = SpatialReference.rotate(ra.x * length_mult,
                                   ra.y * length_mult,
                                   theta=rot)
    z = ra.z
    geoms.append(LineString(list(zip(x, y, z))))
geoms[0]
geoms[0].plot()
fig, ax = plt.subplots()
for g in geoms:
    g.plot(ax=ax)
ax.autoscale()
ax.set_aspect(1)
Explanation: Other geometry types
Linestring
create geometry objects for pathlines from a MODPATH simulation
plot the paths using the built in plotting method
End of explanation
eptfile = EndpointFile('../data/mp6/EXAMPLE-3.endpoint')
eptdata = eptfile.get_alldata()
x, y = SpatialReference.rotate(eptdata['x0'] * length_mult,
                               eptdata['y0'] * length_mult,
                               theta=rot)
z = eptdata['z0']
geoms = [Point(x[i], y[i], z[i]) for i in range(len(eptdata))]
fig, ax = plt.subplots()
for g in geoms:
    g.plot(ax=ax)
ax.autoscale()
ax.set_aspect(2e-6)
Explanation: Points
End of explanation
<END_TASK>