SQL使用PostgreSQL
GeoDjango
- 判断点在边界内
1. point 和 boundary 都存在sql中
The proper way to check whether a Point is
contained by a MultiPolygon is to use
point.intersects(multipolygon).
>>> Rental.objects.filter(location__intersects=preferences.locations)
[<Rental: Rental object>, <Rental: Rental object>]
2. boundary 存在sql中
from django.contrib.gis.geos import Point
from dqchina.dt_crawler.models import BusinessCircle
bc = BusinessCircle.objects.filter(business_circle_id=1372)
point = Point((121.367859, 28.588330))
bc.filter(boundary__contains=point2) # 不在边界内返回空列表
3. point 存在sql中
Points.objects.filter(center_point__intersects=multipolygon)
from django.contrib.gis.geos import Polygon, MultiPolygon
s = ((120.31760343703833, 30.318206356552544),
(120.31710114001254, 30.318962214074713),
(120.31764131672955, 30.318961186851823),
(120.31760343703833, 30.318206356552544))
polygon = Polygon(s)
multipolygon = MultiPolygon(polygon, polygon...)
or
mp = MultiPolygon()
mp.append(polygon)
mp.append(polygon)
SQL_RAW = """
SELECT name, id, total_hushu, tencent_lng, tencent_lat FROM economic_dqchina_loupaninfo
WHERE st_contains(st_geomfromtext('{}', 4326) :: geometry, center_point :: geometry)""".format(mp.wkt)
- 判断两个多边形是否相交的部分
有两种判断包含的方法，一个是contains，一个是intersects
区别：
1. intersects: 表示包含和相交，没有主动和被动的关系
2. contains：表示包含，有主被动关系
例：
p1 = Polygon([[0,0], [2, 0], [2,2], [0, 2], [0,0]])
p2 = Polygon([[1,1], [3, 1], [3,3], [1, 3], [1,1]])
重合部分的多边形：
p3 = (p1 & p2) or p3 = p1.intersection(p2)
p1.intersects(p2) True
p2.intersects(p1) True
p1.intersects(p3) True
p3.intersects(p1) True
p1.contains(p2) False
p2.contains(p1) False
p1.contains(p3) True
p3.contains(p1) False
- 合并两个部分重合的多边形
example.png
import Polygon, MultiPolygon
p1 = Polygon(..)
p2 = Ploygon(..)
mp = MultiPolygon(p1, p2)
or
mp = MultiPolygon()
mp.append(p1)
mp.append(p2)
1. p1.union(p2).wkt
2. mp.unary_union.wkt
- 一定距离内的点
from django.contrib.gis.measure import D
Model.objects.filter(center_point__distance_lte=(center_point, D(km=3)))
Model.objects.filter(center_point__distance_lte=(center_point, 300))
# 注：
友情提示，这个select效率很低。。大量用的话推荐还是推荐用sql语句
- 两点间距离
from math import radians, cos, sin, asin, sqrt
def haversine(lon1, lat1, lon2, lat2, default="GCJ02"):
"""
Calculate the great circle distance between two points
on the earth (specified in decimal degrees)
origin: https://stackoverflow.com/questions/15736995/how-can-i-quickly-estimate-the-distance-between-two-latitude-longitude-points
"""
if default == "WGS84":
pass
elif default == "GCJ02":
lon1, lat1 = gcj02towgs84(lon1, lat1)
lon2, lat2 = gcj02towgs84(lon2, lat2)
pass
elif default == "BD09":
lon1, lat1 = bd09towgs84(lon1, lat1)
lon2, lat2 = bd09towgs84(lon2, lat2)
else:
pass
# convert decimal degrees to radians
lon1, lat1, lon2, lat2 = map(radians, [lon1, lat1, lon2, lat2])
# haversine formula
dlon = lon2 - lon1
dlat = lat2 - lat1
a = sin(dlat / 2) ** 2 + cos(lat1) * cos(lat2) * sin(dlon / 2) ** 2
c = 2 * asin(sqrt(a))
km = 6367 * c
return km
-
地图下方增加wkt调试窗口
wkt.png
如图
class xxxAdmin(admin.OSMGeoAdmin):
display_wkt = True
display_srid = True
modifiable = False (该边界不允许修改)
- 热力图
需要用到postgis的聚合(Cluster)
SELECT
ST_X(ST_Centroid(gc)),
ST_Y(ST_Centroid(gc)),
ST_NumGeometries(gc)
FROM (
SELECT unnest(ST_ClusterWithin(A.center_point::GEOMETRY, 1000::FLOAT / 111195)) gc
FROM economic_dqchina_xiezilou AS A WHERE city_id=310100
) f;
