python - Testing point with in/out of a vector shapefile -


here question.

1. intro

  • a shapefile in polygon type represent study area

http://i8.tietuku.com/08fdccbb7e11c0a9.png

  • some point located in whole rectangle map

http://i8.tietuku.com/877f87022bf817b8.png

i want test whether each point located within/out polygon , further operation(for example, sum grid point amount within study area)

2. idea

i have 2 methods information on stack overflow.

2.1 idea

rasterize shapefile raster file , test.

i haven't done yet, have asked 1 question here , answer.

2.2 idea b

i have tried using poly.contain() test scatter point's location, result wasn't match reality.

3. code based on idea b:

for example:

  • original data represent pt(a pandas dataframe) contain 1000 grids x,y.
  • shapefile shown study area, want filter original data leaving point within area.
3.1 preparation 
# map 4 boundaries xc1,xc2,yc1,yc2 = 113.49805889531724,115.5030664238035,37.39995194888143,38.789235929357105 # grid definition lon_grid  = np.linspace(x_map1,x_map2,38) lat_grid  = np.linspace(y_map1,y_map2,32)
3.1 preparation 
# generate (lon,lat)    xx = lon_grid[pt.x.iloc[:].as_matrix()] yy = lat_grid[pt.y.iloc[:].as_matrix()]  sh = (len(xx),2) data = np.zeros(len(xx)*2).reshape(*sh) in range(0,len(xx),1):     data[i] = np.array([xx[i],yy[i]])  # reading shapefile  map = basemap(llcrnrlon=x_map1,llcrnrlat=y_map1,urcrnrlon=x_map2,\               urcrnrlat=y_map2) map.readshapefile('/xx,'xx')
3.2 test 
patches=[] info, shape in zip(map.xxx_info, map.xxx):     x,y=zip(*shape)     patches.append(polygon(np.array(shape), true) ) poly in patches:      mask = np.array([poly.contains_point(xy) xy in data])
  • then, have numpy array mask value 0,1 represent within/out.
  • combine mask pt ==> pt = pt[[pt.mask == 1]], can filter points

but problem using poly,contains_point(xy), couldn't results match attempt.

an example idea 2

sum value 0,1:

unique, counts = np.unique(mask, return_counts=true)       print np.asarray((unique, counts)).t #result:   > [[0 7]     [1 3]]

http://i4.tietuku.com/7d156db62c564a30.png

from unique value, there must 3 point within shapefile area, result shows 1 point besides.

another test 40 points

http://i4.tietuku.com/5fc12514265b5a50.png

4. problem

the result wrong, , haven't figured out.
think problem may happen 2 reasons:

  • the polygon shapefile wrong(a simple polygon don't think problem remains here).
  • using poly.contains_point(xy) incorrect.

add 2016-01-16

thanks answer, reason found out shapefile itself.
when change shapely.polygon, works well.

here code , result 

c =    fiona.open("xxx.shp") pol = c.next() geom = shape(pol['geometry']) poly_data = pol["geometry"]["coordinates"][0] poly = polygon(poly_data) ax.add_patch(plt.polygon(poly_data))  xx = lon_grid[pt_select.x.iloc[:].as_matrix()] yy = lat_grid[pt_select.y.iloc[:].as_matrix()]  sh = (len(xx),2) points = np.zeros(len(xx)*2).reshape(*sh) in range(0,len(xx),1):     points[i] = np.array([xx[i],yy[i]]) mask = np.array([poly.contains(point(x, y)) x, y in points])  ax.plot(points[:, 0], points[:, 1], "rx") ax.plot(points[mask, 0], points[mask, 1], "ro")

http://i4.tietuku.com/8d895efd3d9d29ff.png

you can use shapely:

import numpy np shapely.geometry import polygon, point  poly_data = [[0, 0], [0, 1], [1, 0], [0.2, 0.5]] poly = polygon(poly_data)  points = np.random.rand(100, 2)  mask = np.array([poly.contains(point(x, y)) x, y in points])

and here plot code:

import pylab pl

fig, ax = pl.subplots() ax.add_patch(pl.polygon(poly_data)) ax.plot(points[:, 0], points[:, 1], "rx") ax.plot(points[mask, 0], points[mask, 1], "ro")

the output:

enter image description here

you can use multipoint speed calculation:

from shapely.geometry import polygon, multipoint  poly_data = [[0, 0], [0, 1], [1, 0], [0.2, 0.5]] poly = polygon(poly_data) points = np.random.rand(100, 2) inside_points = np.array(multipoint(points).intersection(poly))

you can use polygon.contains_point() in matplotlib:

poly = pl.polygon(poly_data) mask = [poly.contains_point(p) p in points]

Comments

Post a Comment

Popular posts from this blog

get url and add instance to a model with prefilled foreign key :django admin -

Image
i have page in each contest shown in row. each row contains contest id , link other model names have foreign key contest. news model: class news(models.model): title = models.charfield(max_length=50, null=true) description = models.textfield() time = models.timefield( default= datetime.datetime.now().time()) contest = models.foreignkey(contest) i want use admin adding , changing these models. want when user clicks on 'news' has foreign key contest id=1, shows list of news value foreign key. use {% instance in instances %} <tr class="{% cycle 'row1' 'row2' %}"> <td>{{ instance.title }}</td> <td>{{ instance.starttime }}</td> <td>{{ instance.timespan }}</td> <td>{{ instance.date }}</td> <td><a href="/adm...
Read more

css - Make div keyboard-scrollable in jQuery Mobile? -

from can tell, although jquery-mobile-powered pages can contain divs overflow set scroll or auto, , these divs can scrolled one-screen bar or mouse wheel, cannot scrolled using arrow keys, page-up/page-down, or home/end. instead, official "page" div (with data-role="page") absorbs input. perhaps other divs can't acquire focus, i'm not sure. is there way around this? edit: jsfiddle of simple example: https://jsfiddle.net/qogz0shx/ <script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.4/jquery.min.js"></script> <script src="https://ajax.googleapis.com/ajax/libs/jquerymobile/1.4.5/jquery.mobile.min.js"> </script> <link rel="stylesheet" href="https://code.jquery.com/mobile/1.4.5/jquery.mobile-1.4.5.css"/> <style> #outer { overflow:scroll; height: 50vh; width: 50vw; } #inner { height: 500vh; width: 500vw; } ...
Read more

ruby on rails - Seeing duplicate requests handled with Unicorn -

recently upgraded our ec2 instances c3.large in order boost more unicorn workers handle increased traffic our site. (6xworkers per machine - 2 ec2 instances). when did this, started seeing duplicate records being created! looks somehow, maybe 2 unicorn workers attempting process same request? in nginx logs see 1 [post] request particular rails controller - yet 2 records being created in database. i'm hitting kind of race condition issue - have no idea how debug this. doesn't happen time, when it's frustrating well. noticed 2 database records within 5 seconds of each other. (and yes disable button when user clicks submit button). unicorn 4.8.2 ruby 1.9.3 rails 3.2.14 any great! thanks! although may not fix problem, upgrade ruby 2.0.0. there blog mentioning idea upgrade ruby 2.0.0 if using unicorn. for example: https://www.digitalocean.com/community/articles/how-to-optimize-unicorn-workers-in-a-ruby-on-rails-app specifically: if using ruby 1.9, ...
Read more