使用 LINQ 指向多边形或多边形上的点

Question

正如之前的问题How to Zip enumerable with itself 所述，我正在研究一些基于点列表的数学算法。我目前正在研究多边形中的点。我有如何做到这一点的代码，并且在 SO 上找到了几个很好的参考资料，例如这个链接 Hit test。所以，我可以弄清楚一个点是否在多边形中。作为确定的一部分，我想确定该点是否实际上在多边形上。这个我也可以。如果我能做到这一切，你可能会问我什么问题？

我可以使用 LINQ 有效地做到这一点吗？我已经可以执行以下操作（假设我之前的问题以及我的问题/答案链接到的链接中描述的 Pairwise 扩展方法，并假设具有 X 和 Y 成员的 Position 类型）。我没有进行太多测试，所以 lambda 可能不是 100% 正确的。此外，它不会考虑非常小的差异。

编辑 请注意，此代码的最新编辑已从 TakeWhile 更改为 TakeWhileInclusive。请参阅本文末尾的扩展方法。

//
// Extend a ray from pt towards the left.  Does this ray intersect the segment p1->p2?
// By definition, the ray extending from pt cannot intersect a horizontal segment, so our
// first check is to see whether or not the segment is horizontal.  
// If it is, there cannot be an intersection.
// If it is not, there could be an intersection.
//
public static bool PointIntersectSegment(Position p1, Position p2, Position pt)
{
  bool intersect = false;
  if (p1.Y != p2.Y)
  {
    // Is pt between (vertically) p1 and p2?  
    // If so, the ray from pt might intersect.  
    // If not, the ray from pt cannot intersect.
    if ((p1.Y >= pt.Y && p2.Y < pt.Y) || (p1.Y < pt.Y && p2.Y >= pt.Y))
    {
      if (p1.X < pt.X && p2.X < pt.X)
      {
        // If the segment is to the left of pt, then the ray extending leftwards from pt will intersect it.
        intersect = true;
      }
      else
        if ((p1.X < pt.X || p2.X < pt.X))
        {
          // If either end of the segment is to the left of pt, calculate intersection (x only) of the 
          // ray from pt and the segment.  If the intersection (x only) is to the left of pt, then
          // the ray extending leftwards from pt will intersect it.
          double inverseSlope = (p1.X - p2.X) / (p1.Y - p2.Y);
          double intersectionX = (pt.Y - p1.Y) * inverseSlope + p1.X;
          if (intersectionX < pt.X)
          {
            intersect = true;
          }
        }
    }
  }
  return intersect;
}

public static bool PointOnSegment(Position p1, Position p2, Position pt)
{
  // Obviously, this is not really going to tell us if pt is on the segment p1->p2.  I am still
  // working on that.  For testing the PointInPolygon algorithm, I can still simulate the "on"
  // case by passing in a pt that is equal to one of the points on the polygon.
  return (pt == p1 || pt == p2);
}

public static PointInPolygonLocation PointInPolygon(IEnumerable<Position> pts, GTPosition pt)
{
  //
  // Implemention of the Jordan Curve theorem to determine if a point is in a polygon.  In essence,
  // this algorithm extends a ray infinitely from pt.  In my implementation I am extending to the left.
  // If the point is inside the polygon, then the ray will intersect the polygon a odd number of times.
  // If the point is outside the polygon, then the ray will intersect the polygon an even number of times.
  // Ideally, we would be able to report not only if the point is inside or outside of the polygon, but
  // also if it is "on" the polygon (i.e. it lies on one of the segments).  Note that "on" and 
  // inside/outside are exlusive.  The point is either on the polygon or it is inside or outside, it
  // cannot be "inside and on" or "outside and on".
  // So, the algorithm is as follows:
  // 1.  Consider the points of the polygon as pairs making up the segments (p1->p2, p2->p3, etc).
  // 2.  For each segment, perform two calculations:  
  //       Does the ray extending left from pt intersect the segment?
  //       Is pt on the segment p[i], p[i+1]?
  // 3.  Count the total number of intersections.  If odd, point is inside.  If even, point is outside.
  // 4.  Here is the tricky part, if the point is on any segment, then we can stop.  If it is on the 
  //     first segment, then there is no need to go through the on/off calculation and the intersection
  //     calculation for the rest of the segments.
  //
  var result = pts.Pairwise( (p1, p2) =>
                         {
                           int intersect = 0;
                           int on = 0;

                           on = PointOnSegment(p1, p2, pt) ? 1 : 0;
                           if (on == 0)
                           {
                             //Don't really need to determine intersection if we already know that pt is on p1->p2.
                             intersect = PointIntersectSegment(p1, p2, pt) ? 1 : 0;
                           }
                           return new { on, intersect };
                         })
                         .TakeWhileInclusive((item) => item.on == 0) //Only consider segments until (or if) pt is on a segment.
                         .Aggregate((curr, next) => new     //Keep a running total of how many intersections.
                                                    { 
                                                      on = curr.on + next.on, 
                                                      intersect = curr.intersect + next.intersect 
                                                    });
  if (result.on != 0)
  {
    return PointInPolygonLocation.On;
  }
  else
  {
    return result.intersect % 2 == 0 ? PointInPolygonLocation.Outside : PointInPolygonLocation.Inside;
  }
}

PointInPolygon 这个函数接受输入位置 pt，迭代输入的位置值序列，并使用 Jordan Curve 方法确定从 pt 向左延伸的光线与多边形相交的次数。 lambda 表达式将在“压缩”列表中生成 1 用于交叉的每个段，而 0 用于其余部分。这些值的总和决定了 pt 是在多边形内部还是外部（奇数 == 内部，偶数 == 外部）。到目前为止，一切顺利。

现在，对于序列中任何连续的位置值对（即在 lambda 的任何执行中），我们还可以确定 pt 是否在段 p1、p2 上。如果是这样，我们可以停止计算，因为我们有答案。

最终，我的问题是：我可以执行这个计算（也许使用聚合？）这样我们将只对序列进行不超过 1 次的迭代，并且如果我们遇到 pt 为 ON 的段，我们是否可以停止迭代?换句话说，如果 pt 在第一个段上，则无需检查其余段，因为我们有答案。

很可能此操作（尤其是可能提前停止迭代的要求/愿望）并不适合 LINQ 方法。

我突然想到，也许 lambda 表达式可以产生一个元组、交集值（1 或 0 或者可能是真或假）和“on”值（真或假）。也许那时我可以使用 TakeWhile(anontype.PointOnPolygon == false)。如果我对元组求和并且如果 ON == 1，则该点位于多边形上。否则，元组其他部分之和的奇数或偶数表明该点是在内部还是在外部。

编辑 清理代码格式，将 lambda 移动到独立函数，如果点在段上，则添加一个新函数。新函数实际上是一个伪函数，它只是将输入点与段的起点和终点进行比较，如果两者相等则返回 true。因此，它可能应该称为 PointIsOnEndpointOfSegment。这里已经有很多代码了，我不想用更多的玩弄 x 和 y 来进行“真正的”ON 计算来混淆这个问题。

使用现在的代码，我可以执行以下操作： 1. 将点序列视为一对（段）序列（通过 Pairwise）。 2. 在 Pairwise 结果选择器中，如果从 pt 向左延伸的射线与线段相交，我计算一个匿名类型，该类型包含值 intersection=1，如果点在线段“上”，则值 on=1。 3. 将匿名类型的交集表示为一个序列。

我认为我可以使用 TakeWhile 来获取所有匿名类型，直到（或如果）其中一种匿名类型指示 pt 在段上。然后，使用 Aggregate 对交集和 on 的值求和。如果 on != 0 的总和，则该点位于其中一个段上。否则，交点的总和将告诉我们该点是在内部（奇数）还是在外部（偶数 - 或零）。

不幸的是，我遇到了与 TakeWhile 类似的问题，如下所述：

How to TakeWhile PLUS one more element.

My TakeWhile 从序列中获取所有项目，但不包括指示交叉点的项目（如果有）。因此，当我汇总结果时，相交项不存在。

看起来人们之前处理这种情况的一种方法是编写一个类似于 TakeWhile 的 TakeUntil 扩展，但它包含第一个不符合谓词的项目。这样的 TakeUntil 扩展会强制评估整个序列吗？

这主要是一个练习，看看是否可以使用 Linq 实现此算法，但需要满足以下要求： 1. 仅当输入点不在任何段上时，才会迭代整个输入点序列。 2. 如果输入点在其中一个段上，则对输入点的迭代应该停止。如果输入点位于 1000 点多边形的第一段，则无需对其余段进行交集和“开”计算。

我们之前在 C++ 中使用 for 循环在点上实现了这个算法，如果点曾经在一个段上，就会中断。我当然可以实现相同的方式，我只是想看看我是否可以 LINQify 而不产生比“旧”for 循环方式更多的交互。

我可能会尝试 TakeUntil 方法，看看会发生什么。

编辑

使用此代码（TakeWhileInclusive IEnumerable 几乎与关于 TakeWhile PLUS 另一个元素的链接中描述的 TakeUntil 相同）：

// Mimic TakeWhile's overload which takes an index as a parameter to the predicate.
public static IEnumerable<T> TakeWhileInclusive<T>(this IEnumerable<T> seq, Func<T, int, bool> predicate)
{
  int i = 0;
  foreach( T e in seq)
  {
    if (!predicate(e,i))
    {
      // If here, then this is first item to fail predicate.  Yield this item and then break.
      yield return e;
      yield break;
    }
    // yield each item from the input sequence until end of sequence or first failure (see above).
    yield return e;
    i++;
  }
}

//
// First saw this here: http://blog.foxxtrot.net/2009/06/inclusively-take-elements-using-linq-and-custom-extensions.html
// TakeWhileInclusive - IEnumerable extension to mimic TakeWhile, but also to return the first
// item that failed the predicate.
// e.g.  seq = 1 2 3 4
// TakeWhile(p => p != 3) will yield 1 2
// TakeWhileInclusive(p => p != 3) will yield 1 2 3
// e.g.  seq = 0 0 0 1 0
// TakeWhile(p => p == 0) will yield 0 0 0
// TakeWhileInclusive(p => p == 0) will yield 0 0 0 1
// Similar to TakeUntil from here: https://stackoverflow.com/questions/2242318/how-could-i-take-1-more-item-from-linqs-takewhile
//
public static IEnumerable<T> TakeWhileInclusive<T>(this IEnumerable<T> seq, Func<T, bool> predicate)
{
  return seq.Select((x, i) => new { Item = x, Index = i })
            .TakeWhileInclusive((x, i) => predicate(x.Item))
            .Select(x => x.Item);
}

并且我的原始代码中的相应更改（将 TakeWhile 替换为 TakeWhileInclusive），我可以得到答案（点在多边形上或者它在内部或外部）并且迭代停止（通过 TakeWhileInclusive）在它命中后该点为 ON 的多边形段。我会再次指出，我的“PointOnSegment”代码是伪造的，但只要我正在测试一个等于多边形点之一的点，就可以进行测试。

我将暂时保留此内容，以防其他人想发表评论。现在我倾向于接受我自己的答案，因为它符合我的意图。这种特殊的方法最终是否是一个好的方法，我还不知道。

Answer 1

我不太了解算法本身，但我想我可以提供帮助：
假设我了解很多，您基本上需要的是检查是否有任何元组适合某个谓词。然后，您可以将每对 Positions 放入一个元组中，并将其作为 IEnumerable<Tuple<Position, Position>> 传递，您可以在其上执行 .Any(predicate goes here) 以查看它们中的任何一个是否符合描述（或者，如果您需要实际谓词为真的值，使用 .FirstOrDefault() 并检查 null)
所以你的代码看起来像这样：

   public static PointInPolygonLocation PointInPolygon(IEnumerable<Position> pts, Position pt)
    {
        bool isIn = pts.Pairwise((p1, p2) => Tuple.Create(p1, p2)).Any(tuple => tuple.Item1.Y != tuple.Item2.Y &&
        ((tuple.Item1.Y >= pt.Y && tuple.Item2.Y < pt.Y) || (tuple.Item1.Y < pt.Y && tuple.Item2.Y >= pt.Y))
        && ((tuple.Item1.X < tuple.Item1.X && tuple.Item2.X < pt.X) || ((pt.Y - tuple.Item1.Y) * ((tuple.Item1.X - tuple.Item2.X) / (tuple.Item1.Y - tuple.Item2.Y)) * tuple.Item1.X) < pt.X));
        return isIn ? PointInPolygonLocation.Inside : PointInPolygonLocation.Outside;
    }

我会认真考虑在这个可怕的谓词中添加一些 cmets。我把它拉上了拉链。嵌套的 ifs 变成了 &&s 等等。当然，你应该检查一下，但我 90% 确定它和你所做的一样。引用 knuth 的话：Beware of bugs in the above code; I have only proved it correct, not tried it.

编辑：现在可读性更强了，你说呢？

   public static PointInPolygonLocation PointInPolygon(IEnumerable<Position> pts, Position pt)
    {
        bool isIn = pts.Pairwise((p1, p2) => Tuple.Create(p1, p2)).Any(tuple => ReadablePredicate(tuple.Item1, tuple.Item2, pt));
        return isIn ? PointInPolygonLocation.Inside : PointInPolygonLocation.Outside;
    }
    public static bool ReadablePredicate(Position p1, Position p2, Position pt)
    {
        if (p1.Y == p2.Y)
            return false;
        if (!((p1.Y >= pt.Y && p2.Y < pt.Y) || (p1.Y < pt.Y && p2.Y >= pt.Y)))
            return false;
        if (p1.X < pt.X && p2.X < pt.X) // Originally was (p1.X < p1.X && ...). that makes no sense, so I assumed you meant p1.X < pt.X
            return true;
        if ((p1.X < pt.X || p2.X < pt.X) && ((pt.Y - p1.Y) * ((p1.X - p2.X) / (p1.Y - p2.Y)) * p1.X) < pt.X)
            return true;
        return false;
    }

Answer 2

我已经经历了这个问题的几次迭代（哈哈）。最终有两个关键：

1. 我的结果选择器函数计算多边形的每个段的两个值。从该点向左延伸的射线是否与线段相交？ （1 或 0 便于稍后通过 Aggregate 求和）。点在段上吗？ （1 或 0 便于稍后通过 Aggregate 求和）。如果段上的点值为 1（点在段上），则迭代应该停止。

2. 由于我想迭代尽可能少的项目，因此我需要在点位于段上的情况下限制序列。如果该点不在任何段上，我别无选择，只能遍历所有段。如果点在段上，我想退出，但我还需要聚合最后一个值（这导致“迭代直到点在段上”条件失败）。我认为 TakeWhile 会起作用，因为我会使用 { intersect, on } 对直到 on == 1（即直到该点位于段上）。不幸的是，TakeWhile 检测到这种情况并停止迭代，但我的输出（提供给 Aggregate）不包含“失败”条件，因此 Aggregate to 的结果没有反映该点位于段上。

我认为，如果我要通过 LINQ 完成此操作，我将不得不遍历整个序列，然后检查（聚合的）结果以确定该点是否落在了我刚刚处理过的任何段上迭代。使用从上面链接的帖子和代码注释中注明的非 SO 网站启发/复制的 TakeWhileInclusive 扩展方法，我可以获得所需的 TakeWhile + 1 行为。

也许有更好的方法来做到这一点？我不知道。我真的才刚刚开始。 Pairwise 模式适用于许多涉及点序列的几何计算。给定一个点序列和一个用点对表示的算法（p1->p2、p2->p3、p3->p4 等），可以“简单”地用 lambda 表示算法（简单的，如距离）或独立的函数（如本文中的，质心等）。

如果其他人想发表评论或提出建议或提出其他答案，我将暂时保留此内容。最终，我可能会接受这个答案（如原始问题中的附加代码和 cmets 所示。