Python2.7：射线法匹配坐标是否在范围坐标内
判断一个坐标是否在一个坐标范围内，可以示用射线法的方式来判断
因为这我这个是栅格匹配景区，里面会出现scenic_id（景区编码） 以及 grid_id（栅格编码）等字眼
首先先说一下射线法
就是以一个点位基准，像某一方向发射线，判断这根线与图形边缘的交点有几个，然后判断是否在图形范围内，其中如果交点位基数，则判断在图形内，交点位偶数，则判断在图形外。
如图所示，A，C点在图形内，与图形边缘的交点分别是1个和3个，而B点在图形外，与图形边缘交点是4个。在地图上判断坐标点也同样适用。
Python代码
因为是Python初学者，所以很多说法可能用的不是很对，以及会有更加简便的语法，如果你看了这边文章，可以提出来，我会修改的，谢谢。
栅格(grid)类 以及 点(point)类
因为是判断栅格是否在景区内，首先将栅格看作一个类，里面包含了3个属性，分别是grid_id(栅格编码)、lon(栅格中心经度)、lat(栅格中心维度)
class Grid:
grid_id =''
lon = ''
lat = ''
def __init__(self, grid_id,lon, lat):
self.grid_id = grid_id
self.lon = lon
self.lat = lat
之后创建景区的合围点对应的类，为 点类，里面包含了2个属性，分别是lon(栅格中心经度)、lat(栅格中心维度)
class Point:
lon = ''
lat = ''
def __init__(self, lon, lat):
self.lon = lon
self.lat = lat
具体方法
这里的使用的方法借鉴了网上的方法，只是把网上的方法理解一下，处理方式大同小异
方法1：get_polygon_bounds： 获取不规则景区的 四至边界下的矩形范围
方法2：is_point_in_rect：判断点是否在矩形范围外，在外就直接排除，如果直接采用射线法，会浪费大量资源，先简单的排除大部分数据，虽然下面获取数据源有点重复，这个是在获取数据源的时候进一步缩小的范围。
方法3：isRayIntersectsSegment：射线法，判断射线和景区不规范范围边的关系
方法4：is_point_in_polygon：主要是对上述3中方法的整体调用，以及判断交点有几个。得出是否在景区防卫内
其中方法3中的下面的代码，如图解释
xseg=point_end.lon-(point_end.lon-point_start.lon)*(point_end.lat-point.lat)/(point_end.lat-point_start.lat)
其中交点的 纵坐标一定是 py
那么主要计算的就是横坐标：ax-(ax-cx)*(ay-py)/(ay-cy)
class Utils:
# 求外包矩形
def get_polygon_bounds(self,points):
length = len(points)
top = down = left = right = points[0]
for i in range(1, length):
if points[i].lon > top.lon:
top = points[i]
elif points[i].lon < down.lon:
down = points[i]
else:
pass
if points[i].lat > right.lat:
right = points[i]
elif points[i].lat < left.lat:
left = points[i]
else:
pass
top_left = Point(top.lon, left.lat)
top_right = Point(top.lon, right.lat)
down_right = Point(down.lon, right.lat)
down_left = Point(down.lon, left.lat)
return [top_left, top_right, down_right, down_left]
# 判断点是否在外包矩形外
def is_point_in_rect(self,point, polygon_bounds):
top_left = polygon_bounds[0]
top_right = polygon_bounds[1]
down_right = polygon_bounds[2]
down_left = polygon_bounds[3]
return (down_left.lon <= point.lon <= top_right.lon
and top_left.lat <= point.lat <= down_right.lat)
#判断两点和射线的关系
def isRayIntersectsSegment(self,point,point_start,point_end): #[x,y] [lng,lat]
#输入：判断点，边起点，边终点，都是[lng,lat]格式数组
if point_start.lat==point_end.lat: #排除与射线平行、重合，线段首尾端点重合的情况
return False
if point_start.lat>point.lat and point_end.lat>point.lat: #线段在射线上边
return False
if point_start.lat<point.lat and point_end.lat<point.lat: #线段在射线下边
return False
if point_start.lat==point.lat and point_end.lat>point.lat: #交点为下端点，对应point_start
return False
if point_end.lat==point.lat and point_start.lat>point.lat: #交点为下端点，对应point_end
return False
if point_start.lon<point.lon and point_end.lat<point.lat: #线段在射线左边
return False
xseg=point_end.lon-(point_end.lon-point_start.lon)*(point_end.lat-point.lat)/(point_end.lat-point_start.lat) #求交
if xseg<point.lon: #交点在射线起点的左侧
return False
return True #排除上述情况之后
def is_point_in_polygon(self,scenic_id,grid, points):
polygon_bounds = self.get_polygon_bounds(points)
if not self.is_point_in_rect(grid, polygon_bounds):
# print '不在矩形范围内，排除'
return False
length = len(points)
sinsc=0
for i in range(0, length):
j=i+1
if j==length:
j=0
point_start = points[i]
point_end = points[j]
if self.isRayIntersectsSegment(grid,point_start,point_end):
sinsc=sinsc+1
#将交点为基数的左边和栅格编码取出
if sinsc%2==1:
# print scenic_id + '\t' + grid.grid_id
return [scenic_id,grid.grid_id]
else :
return False
数据源
数据来源数据库，因此，需要连接数据库，这里连接的是Gbase数据库，使用了Python中的pymysql包
下面的代码中包含了4个方法。 1.获取景区合围、2.获取景区合围的四至边界，3.获取栅格点（因数据量特别大，点的数量接近1000万，因此限制在景区四至边界点，加快程序运行）、4.将栅格匹配到的景区关系插入数据库。
class Gbase_Conn:
#获取景区合围数据
def get_encircle(self,scenic_id):
print '开始连接数据库'
db = pymysql.connect(host='', port=, user='',
password='', db='', charset='') #获取数据库连接
cursor = db.cursor()
sql="select 字段 from 表名 " \
"where 条件字段= '"+scenic_id+"';" #后面的scenic_id是传入的参数
try:
cursor.execute(sql)
enc_infos=cursor.fetchall()
print '获取数据完成'
except:
print '获取数据失败'
return enc_infos
#获取景区的四至边界
def get_around(self,scenic_id):
print '开始连接数据库'
db = pymysql.connect(host='', port=, user='',
password='', db='', charset='') #获取数据库连接
cursor = db.cursor()
sql_max_lon="select ceil(max(经度字段)*10)/10 from 表名 where 条件字段 = '"+scenic_id+"';" #后面的scenic_id是传入的参数，向上获取最大值，确保不会遗漏点
sql_min_lon="select floor(min(经度字段)*10)/10 from 表名 where 条件字段 = '"+scenic_id+"';" #后面的scenic_id是传入的参数，向下获取最小值，确保不会遗漏点
sql_max_lat="select ceil(max(维度字段)*10)/10 from 表名 where 条件字段 = '"+scenic_id+"';" #后面的scenic_id是传入的参数，向上获取最大值，确保不会遗漏点
sql_min_lat="select floor(min(维度字段*10)/10 from 表名 where 条件字段 = '"+scenic_id+"';" #后面的scenic_id是传入的参数，向下获取最小值，确保不会遗漏点
try:
cursor.execute(sql_max_lon)
max_lon=cursor.fetchall()
cursor.execute(sql_min_lon)
min_lon=cursor.fetchall()
cursor.execute(sql_max_lat)
max_lat=cursor.fetchall()
cursor.execute(sql_min_lat)
min_lat=cursor.fetchall()
print '获取数据完成'
except:
print '获取数据失败'
return [max_lon[0],min_lon[0],max_lat[0],min_lat[0]]
#获取所有景区四至边界内所有栅格
def get_grid_info(self,max_lon,min_lon,max_lat,min_lat):
print '开始连接数据库'
db = pymysql.connect(host='', port=, user='',
password='', db='', charset='') #获取数据库连接
cursor = db.cursor()
sql="select 字段 from 表名 " \
"where 经度字段 <='"+ max_lon +"' and 经度字段>= '"+ min_lon + "'" \
"and 维度字段<= '"+ max_lat +"' and 维度字段>= '"+ min_lat +"';"
try:
cursor.execute(sql)
grid_infos=cursor.fetchall()
print '获取数据完成'
except:
print '获取数据失败'
return grid_infos
#将匹配到的数据插入数据库中
def ins_table(self,scenic_id,mapping):
print '开始连接数据库'
db = pymysql.connect(host='', port=, user='',
password='', db='', charset='') #获取数据库连接
cursor = db.cursor()
# sql_tru="delete from 表名 where 字段 = '"+scenic_id+"';"
# try:
# cursor.execute(sql_tru)
# db.commit()
# except:
# db.rollback()
for map in mapping:
sql="insert into 表名 " \
"values ('%s','%s') " % \
(map[0],map[1]) #我这里的表只存对应关系。
try:
cursor.execute(sql) #将本次爬虫爬取的结果集插入数据库表中
db.commit()
except:
db.rollback()
return
最后把主函数贴出来 主要是最上面的方法调用
class Grid_ascription(object):
def __init__(self):
self.database_connect = Gbase_Conn()
self.utils = Utils()
def main(self,scenic_id):
#获取合围点
enc_infos=self.database_connect.get_encircle(scenic_id)
print enc_infos
points=[]
for enc_info in enc_infos:
lon=float(str(enc_info[2]).strip('0'))
lat=float(str(enc_info[3]).strip('0'))
point=Point(lon,lat)
points.append(point)
#找出景区范围内 最大 最小，最左，最右边界
max_round=self.database_connect.get_around(scenic_id)
max_lon=str(float(str(max_round[0][0]).strip('0')))
min_lon=str(float(str(max_round[1][0]).strip('0')))
max_lat=str(float(str(max_round[2][0]).strip('0')))
min_lat=str(float(str(max_round[3][0]).strip('0')))
# 获取需要判断是否在合围内的点
grid_infos=self.database_connect.get_grid_info(max_lon,min_lon,max_lat,min_lat)
print grid_infos
mapping=[]
for grid_info in grid_infos:
grid_id=grid_info[0]
grid_lon=grid_info[1]
grid_lat=grid_info[2]
grid=Grid(grid_id,grid_lon,grid_lat)
if self.utils.is_point_in_polygon(scenic_id,grid,points) :
map = self.utils.is_point_in_polygon(scenic_id,grid,points)
# print map
mapping.append(map)
self.database_connect.ins_table(scenic_id,mapping)
if __name__ == "__main__":
grid_main = Grid_ascription()
grid_main.main('') #传参景区编码
这个程序亲测可以，但是有一个问题，就是每次只能匹配一个坐标范围，如果需要匹配的坐标范围较多，还是没有很方便。
如果有好方法，可以告知一下，谢谢
浙公网安备 33010602011771号