<SYSTEM_TASK:>
Given the following text description, write Python code to implement the functionality described below step by step
<END_TASK>
<USER_TASK:>
Description:
 determine_region
A notebook to determine what region (e.g. neighborhood, ward, census district) the issue is referring to
Step1:  Read in neighborhood shapefiles
Step3:  Now plot the shapefiles
Step4:  Read in issues and determine the region
In this section I will read determine the region each issue occurs in. First to read in the issues
Step5:  Remove issues that do not have correct coordinates
<ASSISTANT_TASK:>
Python Code:
import fiona
from shapely.geometry import shape
import nhrc2
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
from collections import defaultdict
import numpy as np
from matplotlib.patches import Polygon
from shapely.geometry import Point
%matplotlib inline
#the project root directory:
nhrc2dir = ('/').join(nhrc2.__file__.split('/')[:-1])+'/'
Explanation: determine_region
A notebook to determine what region (e.g. neighborhood, ward, census district) the issue is referring to
End of explanation
c = fiona.open(nhrc2dir+'data/nh_neighborhoods/nh_neighborhoods.shp')
pol = c.next()
geom = shape(pol['geometry'])
c.crs
for i in c.items():
    print(i[1])
len(c)
for i in c:
    pol = i
    geom = shape(pol['geometry'])
    geom
geom
Explanation: Read in neighborhood shapefiles:
End of explanation
#Based on code from Kelly Jordhal: 
#http://nbviewer.ipython.org/github/mqlaql/geospatial-data/blob/master/Geospatial-Data-with-Python.ipynb
def plot_polygon(ax, poly):
    a = np.asarray(poly.exterior)
    ax.add_patch(Polygon(a, facecolor='#46959E', alpha=0.3))
    ax.plot(a[:, 0], a[:, 1], color='black')
def plot_multipolygon(ax, geom):
     Can safely call with either Polygon or Multipolygon geometry
    
    if geom.type == 'Polygon':
        plot_polygon(ax, geom)
    elif geom.type == 'MultiPolygon':
        for poly in geom.geoms:
            plot_polygon(ax, poly)
nhv_geom = defaultdict()
#colors = ['red', 'green', 'orange', 'brown', 'purple']
fig, ax = plt.subplots(figsize=(12,12))
for rec in c:
    #print(rec['geometry']['type'])
    hood = rec['properties']['name']
    nhv_geom[hood] = shape(rec['geometry'])
    plot_multipolygon(ax, nhv_geom[hood])
labels = ax.get_xticklabels() 
for label in labels: 
    label.set_rotation(90) 
ax.set_xlabel('Longitude')
ax.set_ylabel('Latitude')
ax.plot(scf_df.loc[0, 'lng'], scf_df.loc[0, 'lat'], 'o')
Explanation: Now plot the shapefiles
End of explanation
import nhrc2.backend.read_seeclickfix_api_to_csv as rscf
scf_cats = rscf.read_categories(readfile=True)
scf_df = rscf.read_issues(scf_cats, readfile=True)
len(scf_cats)
scf_df.head(3)
len(scf_df)
Explanation: Read in issues and determine the region
In this section I will read determine the region each issue occurs in. First to read in the issues:
End of explanation
scf_df = scf_df[((scf_df['lat'] < 41.36) & (scf_df['lat'] > 41.24) & (scf_df['lng']>=-73.00) & (scf_df['lng'] <= -72.86))]
print(len(scf_df))
scf_df.loc[0, 'lat']
grid_point = Point(scf_df.loc[0, 'lng'], scf_df.loc[0, 'lat']) 
for idx in range(5):
    grid_point = Point(scf_df.loc[idx, 'lng'], scf_df.loc[idx, 'lat']) 
    print('Point {} at {}'.format(idx, scf_df.loc[idx, 'address']))
    print('Downtown: {}'.format(grid_point.within(nhv_geom['Downtown'])))
    print('East Rock: {}'.format(grid_point.within(nhv_geom['East Rock'])))
    print('Fair Haven Heights: {}'.format(grid_point.within(nhv_geom['Fair Haven Heights'])))
print('Number of neighborhoods: {}'.format(len(nhv_geom.keys())))
for hood in nhv_geom.keys():
    print(hood)
def get_neighborhoods(scf_df, neighborhoods):
    hoods = []
    for idx in scf_df.index:
        grid_point = Point(scf_df.loc[idx, 'lng'], scf_df.loc[idx, 'lat']) 
        for hoodnum, hood in enumerate(nhv_geom.keys()):
            if grid_point.within(nhv_geom[hood]):
                hoods.append(hood)
                break
            if hoodnum == 19:
                #There are 20 neighborhoods. If you are the 20th (element 19 in 
                #zero-based indexing) and have not continued out of the iteration
                #set the neighborhood name to "Other":
                hoods.append('Other')
    return hoods
            
%time nbrhoods = get_neighborhoods(scf_df, nhv_geom)
print(len(scf_df))
print(len(nbrhoods))
Explanation: Remove issues that do not have correct coordinates:
End of explanation
<END_TASK>