纯Go语言判断点是否在多边形内

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
package main

import (
	"encoding/json"
	"fmt"
	"log"

	"github.com/tidwall/geojson/geometry"
)

type F struct {
	Polys []*geometry.Poly
}

func (f *F) Insert(poly *geometry.Poly) {
	f.Polys = append(f.Polys, poly)
}

func (f *F) Contains(x float64, y float64) bool {
	p := geometry.Point{X: x, Y: y}
	for _, poly := range f.Polys {
		if poly.ContainsPoint(p) {
			return true
		}
	}
	return false
}

func main() {
	// 多边形坐标集字符串,geojson.io可以生成
	polygonString := `[[
        [84.06629085202707,61.09246473319979],
        [72.49198357244117,54.2478482682468],
        [76.04595361631505,49.145991557702814],
        [113.30166443088888,51.79338737652475],
        [114.95645443380147,55.102615742572425],
        [109.48319438253571,59.55098721892466],
        [104.69660340886469,60.03279304645392],
        [95.37143615702166,53.7249665598695],
        [105.9540342347463,55.274455891821276],
        [92.46849779684248,51.413269843450934],
        [90.8684973410422,58.41496021967512],
        [84.06629085202707,61.09246473319979]
        ]]`
	// 解析 JSON 字符串以获取多边形坐标
	var coords [][][2]float64
	if err := json.Unmarshal([]byte(polygonString), &coords); err != nil {
		log.Fatalf("Error parsing JSON: %v", err)
	}

	// 将坐标转换为 geometry.Point 类型
	var points []geometry.Point
	for _, coord := range coords[0] {
		points = append(points, geometry.Point{X: coord[0], Y: coord[1]})
	}

	// 创建多边形
	poly := geometry.NewPoly(points, nil, geometry.DefaultIndexOptions)

	// 创建一个查找器实例
	finder := F{}

	// 将多边形插入查找器
	finder.Insert(poly)

	// 使用提供的坐标检查点是否在多边形内
	pointX, pointY :=  87.24702218937767, 54.08457905522715 // 在圈内
	contains := finder.Contains(pointX, pointY)
	fmt.Println("Point is in polygon:", contains)
	pointX, pointY =      96.92668165010213, 47.6979436832691 // 不在圈内
	contains = finder.Contains(pointX, pointY)
	fmt.Println("Point is in polygon:", contains)
}

上述代码中提到的多边形和点的geojson文件(可在geojson.io中粘贴复现图形):

1
{"type":"FeatureCollection","features":[{"type":"Feature","properties":{},"geometry":{"coordinates":[[[84.06629085202707,61.09246473319979],[72.49198357244117,54.2478482682468],[76.04595361631505,49.145991557702814],[113.30166443088888,51.79338737652475],[114.95645443380147,55.102615742572425],[109.48319438253571,59.55098721892466],[104.69660340886469,60.03279304645392],[95.37143615702166,53.7249665598695],[105.9540342347463,55.274455891821276],[92.46849779684248,51.413269843450934],[90.8684973410422,58.41496021967512],[84.06629085202707,61.09246473319979]]],"type":"Polygon"}},{"type":"Feature","properties":{"marker-color":"#ed1d1d","marker-size":"medium","marker-symbol":"circle"},"geometry":{"coordinates":[96.92668165010213,47.6979436832691],"type":"Point"},"id":1},{"type":"Feature","properties":{},"geometry":{"coordinates":[87.24702218937767,54.08457905522715],"type":"Point"}}]}

world.geojson(大小:23+MB)生成自https://geojson-maps.ash.ms,选择High resolutionAll

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
import json,time
from shapely.geometry import Point, shape
start = time.time()
with open('world.geojson','r') as f:
    geojson_data = json.load(f)
# Given point
point = [92.5976562500021, 42.58220940617247]

# Function to check if a point is inside a polygon
def is_point_in_polygon(point, polygon):
    point = Point(point)
    polygon = shape(polygon)
    return polygon.contains(point)

# Check which feature contains the point
# feature_containing_point = None
for feature in geojson_data['features']:
    if is_point_in_polygon(point, feature['geometry']):
        # feature_containing_point = feature
        print(feature['properties']['name_zh'])
        break
print(time.time()- start)
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
package main

import (
	"encoding/json"
	"fmt"
	"os"
	"time"

	"github.com/tidwall/geojson"
	"github.com/tidwall/geojson/geometry"
)

func main() {
	// 读取 GeoJSON 文件
	start := time.Now()
	file, err := os.ReadFile("world.geojson")
	if err != nil {
		fmt.Println(err)
		os.Exit(1)
	}
	// 解析 GeoJSON
	obj, err := geojson.Parse(string(file), geojson.DefaultParseOptions)
	if err != nil {
		fmt.Println(err)
		os.Exit(1)
	}

	// 你的点坐标
	point := geometry.Point{X: 92.5976562500021, Y: 42.58220940617247}
	geojson_point := geojson.NewPoint(point)
	// 检查点是否在某个 feature 内
	fc, ok := obj.(*geojson.FeatureCollection)
	if !ok {
		fmt.Println("Error: not a FeatureCollection")
		os.Exit(1)
	}

	fc.ForEach(func(geom geojson.Object) bool {
		feature, ok := geom.(*geojson.Feature)
		// fmt.Println(feature.Contains(geojson_point))
		if !ok {
			return true // 继续迭代
		}

		if feature.Contains(geojson_point) {
			// fmt.Println("Feature within point:", point)
			jsonBytes, err := feature.MarshalJSON()
			if err != nil {
				fmt.Println("Error marshaling feature to JSON:", err)
				return true
			}

			// 解析 JSON 以获取属性
			var featureJSON map[string]interface{}
			if err := json.Unmarshal(jsonBytes, &featureJSON); err != nil {
				fmt.Println("Error unmarshaling feature JSON:", err)
				return true
			}

			// 提取 name 属性
			if properties, ok := featureJSON["properties"].(map[string]interface{}); ok {
				if name, ok := properties["name_zh"].(string); ok {
					fmt.Println("Feature Name:", name)
					return false // 停止迭代
				}
			}
		}

		return true // 继续迭代
	})
	fmt.Println("Time:", time.Since(start))
}

https://coloraven.github.io/2023/12/geopandas读取geojson文件并写入Spatialite数据库