Получить разницу двух файлов формы

Question

У меня есть два файла формы, которые содержат полигоны. Я пытаюсь найти дельту из этого. Попытка сделать это с помощью следующего кода, но не работает так, как я ожидал. Ниже приведены два файла формы, один из которых синий, это файл формы буфера, мне нужно удалить ту область буфера, которая пересекается с синим буфером.

import fiona
from shapely.geometry import shape, mapping, Polygon

green = fiona.open(
    "/home/gulve/manual_geo_ingestion/probe-data/images/r/shape_out/dissolved.shp")
blue = fiona.open(
    "/home/gulve/manual_geo_ingestion/probe-data/images/g/shape/shape.shp")

print([not shape(i['geometry']).difference(shape(j['geometry'])).is_empty for i, j in zip(list(blue), list(green))])

schema = {'geometry': 'Polygon',
          'properties': {}}
crs = {'init': u'epsg:3857'}

with fiona.open(
        '/home/gulve/manual_geo_ingestion/probe-data/images/r/shape_out/diff.shp', 'w',
        driver='ESRI Shapefile', crs=crs, schema=schema
) as write_shape:
    for geom in [shape(i['geometry']).difference(shape(j['geometry'])) for i, j in zip(list(blue), list(green))]:
        if not geom.empty:
            write_shape.write({'geometry': mapping((shape(geom))), 'properties': {}})

Ожидаемый результат:

Answer 1

После того, как вы импортировали шейп-файлы в PostgreSQL, просто выполните этот запрос:

CREATE TABLE not_intersects AS 
SELECT * FROM shape 
WHERE id NOT IN (SELECT DISTINCT shape.id
         FROM buffer,shape 
         WHERE ST_Intersects(buffer.geom,shape.geom)); 

Этот запрос создаст третью таблицу (называется здесь not_intersects) содержит многоугольники, которые не пересекаются между двумя таблицами (файлы форм).

Здесь результат представлен желтым цветом.

Answer 2

Я могу решить это с помощью стройных функций

Вот мой код..

import fiona
from shapely.geometry import shape, mapping, Polygon
from shapely.ops import unary_union

buffered_shape = fiona.open(
    "dissolved.shp", 'r', encoding='UTF-8')
color_shape = fiona.open(
    "shape.shp", 'r', encoding='UTF-8')

print([not shape(i['geometry']).difference(shape(j['geometry'])).is_empty for i, j in
       zip(list(color_shape), list(buffered_shape))])

outmulti = []
for pol in color_shape:
    green = shape(pol['geometry'])
    for pol2 in buffered_shape:
        red = shape(pol2['geometry'])
        if red.intersects(green):
            # If they intersect, create a new polygon that is
            # essentially pol minus the intersection
            intersect = green.intersection(red)
            nonoverlap = green.symmetric_difference(intersect)
            outmulti.append(nonoverlap)
        else:
            outmulti.append(green)

finalpol = unary_union(outmulti)

schema = {'geometry': 'MultiPolygon',
          'properties': {}}
crs = {'init': u'epsg:4326'}

with fiona.open(
        'shape_out/diff.shp', 'w',
        driver='ESRI Shapefile', crs=crs, schema=schema
) as write_shape:
    for geom in finalpol:
        # if not geom.empty:
        write_shape.write({'geometry': mapping(Polygon(shape(geom))), 'properties': {}})