Переместите географические точки, которые лежат слишком близко друг к другу (Python / Shapely)

Question

У меня есть набор данных о нефтеперерабатывающих заводах в Техасе (GeoJSON здесь - https://pastebin.com/R0D9fif9):

Name,Latitude,Longitude
Marathon Petroleum,29.374722,-94.933611
Marathon Petroleum,29.368733,-94.903253
Valero,29.367617,-94.909515
LyondellBasell,29.71584,-95.234814
Valero,29.722213,-95.255198
Exxon,29.743865,-95.009208
Shell,29.720425,-95.12495
Petrobras,29.722466,-95.208807

Я хотел бы создать печатную карту из этих точек. Но они лежат слишком близко друг к другу в данной резолюции.

Поскольку каждый легочный завод должен быть упомянут в легенде, я не могу кластеризоваться. Так что я бы хотел

• Получить центроид - это было легко

  import json
  import csv
  from shapely.geometry import shape, Point, MultiPoint
  
  with open('refineries.csv', 'rU') as infile:
      reader = csv.DictReader(infile)
      data = {}
      for row in reader:
          for header, value in row.items():
              try:
                  data[header].append(value)
              except KeyError:
                  data[header] = [value]
  
   listo = list(zip(data['Longitude'], data['Latitude']))
   points1 = MultiPoint(points=listo)
  
   points = MultiPoint([(-94.933611, 29.374722), (-94.903253, 29.368733), (-94.909515, 29.367617), (-95.234814, 29.71584), (-95.255198, 29.722213), (-95.009208, 29.743865), (-95.12495, 29.720425), (-95.208807, 29.722466)])
  
   print(points.centroid)
  

• Смещайте все точки от центра тяжести, пока не будет достигнуто минимальное расстояние между всеми

Можете ли вы помочь мне здесь? Заранее спасибо!

Answer 1

Это зависит от того, как именно вы хотите сдвинуть точки от центра тяжести. Один из способов состоит в том, чтобы рассчитать для каждой точки ее расстояние и азимут по большому кругу относительно центроида и перемасштабировать все расстояния, чтобы гарантировать, что расстояние между двумя ближайшими точками больше указанного порогового значения. В приведенном ниже примере pyproj используется для расчета азимутов и расстояний.

import json
import csv
import sys
from shapely.geometry import shape, Point, MultiPoint
from pyproj import Geod

with open('refineries.csv', 'rU') as infile:
    reader = csv.DictReader(infile)
    data = {}
    for row in reader:
        for header, value in row.items():
            if not header in data:
                data[header] = []
            data[header].append(value)

listo = list(zip(map(float, data['Longitude']), map(float, data['Latitude'])))

def scale_coords(coords, required_dist = 1000.):
    g = Geod(ellps = 'WGS84')

    num_of_points = len(coords)

    #calculate centroid
    C = MultiPoint(coords).centroid

    #determine the minimum distance among points
    dist_min, dist_max = float('inf'), float('-inf')
    for i in range(num_of_points):
        lon_i, lat_i = coords[i]
        for j in range(i+1, num_of_points):
            lon_j, lat_j = coords[j]
            _,_,dist = g.inv(lon_i, lat_i, lon_j, lat_j)
            dist_min = min(dist_min, dist)
            dist_max = max(dist_max, dist)

    if dist_min > required_dist:
        return coords

    coords_scaled = [None]*num_of_points
    scaling = required_dist / dist_min
    for i, (lon_i, lat_i) in enumerate(coords):
        az,_,dist = g.inv(C.x, C.y, lon_i, lat_i)
        lon_f,lat_f,_ = g.fwd(C.x, C.y, az, dist*scaling)
        coords_scaled[i] = (lon_f, lat_f)

    return coords_scaled

В качестве альтернативы, это может сочетаться с подходом, при котором вы также ослабляете азимуты. Это в принципе привело бы к меньшему коэффициенту масштабирования для "радиальных" расстояний. Однако это также немного исказит "визуальное распределение" точек. Кроме того, метод, представленный выше, может быть "улучшен" за счет игнорирования любых выбросов в пересчете, то есть точек, которые уже достаточно далеко от центроида и не имеют соседних соседей.