Menemukan Segmen Garis Terdekat ke Titik menggunakan rupawan?

Latar Belakang

Dari Titik yang diketahui, saya perlu membuat "batas sekeliling" terdekat yang berdekatan dengan tabel MultiLineStrings, seperti yang ditunjukkan pada diagram.

Saya telah mencari situs ini dengan sejumlah istilah (mis. Tepi minimum, batas minimum, tetangga terdekat, klip, yang mengandung poligon, visibilitas, jepret, potong node, jejak sinar, pengisian banjir, batas dalam, rute, cekung lambung kapal) tetapi tidak dapat menemukan pertanyaan sebelumnya yang sepertinya cocok dengan skenario ini.

Diagram

• Lingkaran hijau adalah titik yang dikenal.
• Garis hitam adalah MultiLineStrings yang dikenal.
• Garis abu-abu adalah indikasi sapuan radial dari Titik yang diketahui.
• Titik merah adalah persimpangan terdekat sapuan radial dan MultiLineStrings.

Parameter

• Intinya tidak akan pernah memotong MultiLineStrings.
• Point akan selalu dipusatkan secara nominal di dalam MultiLineStrings.
• MultiLineStrings tidak akan pernah sepenuhnya menyertakan Point, oleh karena itu perimeter akan menjadi MultiLineString.
• Akan ada tabel berisi sekitar 1.000 MultiLineStrings (biasanya berisi satu baris sekitar 100 poin).

Metodologi Dianggap

• Lakukan penyapuan radial dengan membangun serangkaian garis dari Point yang diketahui (pada, katakanlah, kenaikan 1 derajat).
• Tetapkan titik persimpangan terdekat dari setiap garis sapuan radial dengan MultiLineStrings.
• Ketika salah satu garis sapuan radial tidak bersinggungan dengan MultiLineStrings, ini akan menunjukkan celah dalam perimeter yang akan ditampung dalam konstruksi MultiLineString perimeter.

Ringkasan

Sementara teknik ini akan menemukan persimpangan terdekat, itu tidak akan selalu menemukan semua titik perimeter simpul terdekat, tergantung pada resolusi sapuan radial. Adakah yang bisa merekomendasikan metode alternatif untuk menetapkan semua titik perimeter atau melengkapi teknik sapuan radial dengan beberapa bentuk buffering, sektoring atau offset?

Perangkat lunak

Preferensi saya adalah menggunakan SpatiaLite dan / atau Shapely untuk solusi tetapi akan menerima saran yang dapat diimplementasikan menggunakan perangkat lunak sumber terbuka.

Sunting: Solusi yang Berfungsi (berdasarkan jawaban oleh @ gen)

from shapely.geometry import Point, LineString, mapping, shape
from shapely.ops import cascaded_union
from shapely import affinity
import fiona

sweep_res = 10  # sweep resolution (degrees)
focal_pt = Point(0, 0)  # radial sweep centre point
sweep_radius = 100.0  # sweep radius

# create the radial sweep lines
line = LineString([(focal_pt.x,focal_pt.y), \
                   (focal_pt.x, focal_pt.y + sweep_radius)])

sweep_lines = [affinity.rotate(line, i, (focal_pt.x, focal_pt.y)) \
               for i in range(0, 360, sweep_res)]

radial_sweep = cascaded_union(sweep_lines)

# load the input lines and combine them into one geometry
input_lines = fiona.open("input_lines.shp")
input_shapes = [shape(f['geometry']) for f in input_lines]
all_input_lines = cascaded_union(input_shapes)

perimeter = []
# traverse each radial sweep line and check for intersection with input lines
for radial_line in radial_sweep:
    inter = radial_line.intersection(all_input_lines)

    if inter.type == "MultiPoint":
       # radial line intersects at multiple points
       inter_dict = {}
       for inter_pt in inter:
           inter_dict[focal_pt.distance(inter_pt)] = inter_pt
       # save the nearest intersected point to the sweep centre point
       perimeter.append(inter_dict[min(inter_dict.keys())])

    if inter.type == "Point":
       # radial line intersects at one point only
       perimeter.append(inter)

    if inter.type == "GeometryCollection":
       # radial line doesn't intersect, so skip
       pass

# combine the nearest perimeter points into one geometry
solution = cascaded_union(perimeter)

# save the perimeter geometry
schema = {'geometry': 'MultiPoint', 'properties': {'test': 'int'}}
with fiona.open('perimeter.shp', 'w', 'ESRI Shapefile', schema) as e:
     e.write({'geometry':mapping(solution), 'properties':{'test':1}})

Saya telah mereproduksi contoh Anda dengan shapefile.

Anda dapat menggunakan Shapely dan Fiona untuk menyelesaikan masalah Anda.

1) Masalah Anda (dengan rupawan Point):

2) dimulai dengan garis arbitrer (dengan panjang yang memadai):

from shapely.geometry import Point, LineString
line = LineString([(point.x,point.y),(final_pt.x,final_pt.y)])

3) menggunakan shapely.affinity.rotate untuk membuat jari-jari (memutar garis dari titik, lihat juga jawaban Mike Toews di Python, perpustakaan pustaka: apakah mungkin untuk melakukan operasi affine pada poligon bentuk? ):

from shapely import affinity
# Rotate i degrees CCW from origin at point (step 10°)
radii= [affinity.rotate(line, i, (point.x,point.y)) for i in range(0,360,10)]

4) sekarang, menggunakan shapely: cascaded_union (atau shapely: unary_union ) untuk mendapatkan MultiLineString:

from shapely.ops import cascaded_union
mergedradii = cascaded_union(radii)
print mergedradii.type
MultiLineString

5) sama dengan garis asli (shapefile)

import fiona
from shapely.geometry import shape
orlines = fiona.open("orlines.shp")
shapes = [shape(f['geometry']) for f in orlines]
mergedlines = cascaded_union(shapes)
print mergedlines.type
MultiLineString

6) persimpangan antara dua multigeometry dihitung dan hasilnya disimpan ke shapefile:

 points =  mergedlines.intersection(mergedradii)
 print points.type
 MultiPoint
 from shapely.geometry import mapping
 schema = {'geometry': 'MultiPoint','properties': {'test': 'int'}}
 with fiona.open('intersect.shp','w','ESRI Shapefile', schema) as e:
      e.write({'geometry':mapping(points), 'properties':{'test':1}})

Hasil:

7) tetapi, masalah, jika Anda menggunakan radius yang lebih panjang, hasilnya berbeda:

8) Dan jika Anda ingin mendapatkan hasil, Anda hanya perlu memilih titik dengan jarak terpendek dari titik asli pada radius:

points_ok = []
for line in mergeradii:
   if line.intersects(mergedlines):
       if line.intersection(mergedlines).type == "MultiPoint":
          # multiple points: select the point with the minimum distance
          a = {}
          for pt in line.intersection(merged):
              a[point.distance(pt)] = pt
          points_ok.append(a[min(a.keys())])
       if line.intersection(mergedlines).type == "Point":
          # ok, only one intersection
          points_ok.append(line.intersection(mergedlines))
solution = cascaded_union(points_ok)
schema = {'geometry': 'MultiPoint','properties': {'test': 'int'}}
with fiona.open('intersect3.shp','w','ESRI Shapefile', schema) as e:
     e.write({'geometry':mapping(solution), 'properties':{'test':1}})

Hasil akhir:

Saya harap itu yang Anda inginkan.