使用形状的两个椭圆的交点面积

问题描述 投票:0回答:1

我正在整形地使用以找到两个椭圆的重叠区域。

这是椭圆函数:

import torch
from shapely.geometry.point import Point
from shapely import affinity


def create_ellipse(center, lengths, angle=0):
    """
    create a shapely ellipse. adapted from
    https://gis.stackexchange.com/a/243462
    """
    circ = Point(center).buffer(1)
    ell = affinity.scale(circ, int(lengths[0]), int(lengths[1]))
    ellr = affinity.rotate(ell, angle)
    return ellr

有2个pytorch张量:

wh1 = torch.tensor([3.6250, 2.8125])
w1 = wh1[0]
h1 = wh1[1]
wh2 = torch.tensor([[ 1.9215,  8.2504,  0.2698,  2.1242,  0.3143,  2.1111,  1.9065,  1.3590,
          1.9489,  1.3957,  0.7099,  1.2371,  0.9294,  0.3765,  0.3079,  0.2556,
          0.2194,  0.2841,  0.3109,  0.8383,  0.3096,  0.2565,  0.3197,  0.5158,
          0.3186,  0.3041,  0.5445,  0.2659,  0.5199,  7.7250,  6.2023,  7.0652,
          0.6808,  0.8072,  1.4054,  0.7802,  2.3664,  2.8730,  0.4354,  1.4630,
          4.7688,  4.6711,  1.5979,  2.4067,  1.7955,  1.8868,  1.8641,  1.5332,
          1.5486],
        [ 1.7919,  4.1976,  0.3816,  1.9401,  0.4061,  2.7557,  2.5397,  1.7599,
          1.8994,  1.5508,  0.4886,  1.1070,  0.9757,  0.7873,  0.4764,  0.3846,
          0.4159,  0.4954,  0.3313,  0.8940,  0.5923,  0.3516,  0.5340,  0.3600,
          0.4161,  0.2797,  0.8646,  0.5811,  0.4164,  0.8625,  5.2044, 10.5977,
          0.8569,  0.6457,  2.4552,  0.8479,  3.8021,  2.6354,  1.2007,  1.7143,
          3.6156,  7.8706,  2.1638,  7.8473,  1.9153,  2.0610,  2.2354,  1.9746,
          1.9436]])
w2 = wh2[0]
h2 = wh2[1]

w1和w2是椭圆的宽度,h1和h2是椭圆的高度。我想找到第一个椭圆(w1和h1)以及每个带有w2和h2的椭圆的重叠区域。输出的大小应与w2或h2相同。我运行以下命令:

ell1 = create_ellipse((0, 0), (w1.item() / 2.0, h1.item() / 2.0), 0)
for w2_t, h2_t in zip(w2, h2):
  print('w2: ', w2_t.item(), 'h2: ', h2_t.item())
  ell2 = create_ellipse((0, 0), (w2_t.item() / 2.0, h2_t.item() / 2.0), 0)
  inter_area = ell1.intersection(ell2).area
  print('inter area: ', inter_area)

然后显示以下错误:

TopologyException: Input geom 1 is invalid: Self-intersection at or near point -1 0 at -1 0
w2:  1.9214999675750732 h2:  1.7919000387191772
inter area:  0.0
w2:  8.250399589538574 h2:  4.1975998878479
inter area:  3.1365484905459384
w2:  0.26980000734329224 h2:  0.3815999925136566
inter area:  0.0
w2:  2.1242001056671143 h2:  1.9400999546051025
---------------------------------------------------------------------------
TopologicalError                          Traceback (most recent call last)
<ipython-input-143-68e4119a1c74> in <module>()
      3   print('w2: ', w2_t.item(), 'h2: ', h2_t.item())
      4   ell2 = create_ellipse((0, 0), (w2_t.item() / 2.0, h2_t.item() / 2.0), 0)
----> 5   inter_area = ell1.intersection(ell2).area
      6   print('inter area: ', inter_area)

2 frames
/usr/local/lib/python3.6/dist-packages/shapely/geometry/base.py in intersection(self, other)
    647     def intersection(self, other):
    648         """Returns the intersection of the geometries"""
--> 649         return geom_factory(self.impl['intersection'](self, other))
    650 
    651     def symmetric_difference(self, other):

/usr/local/lib/python3.6/dist-packages/shapely/topology.py in __call__(self, this, other, *args)
     68             err = TopologicalError(
     69                     "This operation could not be performed. Reason: unknown")
---> 70             self._check_topology(err, this, other)
     71         return product
     72 

/usr/local/lib/python3.6/dist-packages/shapely/topology.py in _check_topology(self, err, *geoms)
     36                     "The operation '%s' could not be performed. "
     37                     "Likely cause is invalidity of the geometry %s" % (
---> 38                         self.fn.__name__, repr(geom)))
     39         raise err
     40 

TopologicalError: The operation 'GEOSIntersection_r' could not be performed. Likely cause is invalidity of the geometry <shapely.geometry.polygon.Polygon object at 0x7f7bc3100940>
python intersection ellipse shapely
1个回答
0
投票

我发现了问题。函数:

def create_ellipse(center, lengths, angle=0):
    """
    create a shapely ellipse. adapted from
    https://gis.stackexchange.com/a/243462
    """
    circ = Point(center).buffer(1)
    ell = affinity.scale(circ, int(lengths[0]), int(lengths[1]))
    ellr = affinity.rotate(ell, angle)
    return ellr

应更改为:

def create_ellipse(center, lengths, angle=0):
    """
    create a shapely ellipse. adapted from
    https://gis.stackexchange.com/a/243462
    """
    circ = Point(center).buffer(1)
    ell = affinity.scale(circ, lengths[0], lengths[1])
    ellr = affinity.rotate(ell, angle)
    return ellr
© www.soinside.com 2019 - 2024. All rights reserved.