修复线剪切算法的代码

Posted 2023-03-23

【中文标题】修复线剪切算法的代码

【英文标题】：Fixing the code for Line Clipping Algorithm

【发布时间】：2015-10-04 12:03:02

【问题描述】：

Midpoint subdivision algorithm [page-93(104)]的工作原理是将一条线分成更小的段，并测试每个段以发现它们是否在剪辑区域的可见边界内。

在二分搜索算法中，我们找到中间元素，然后选择右手边或左手边。

但是，在这里，如下图所示，在第一次分割之后，我们发现这两个小节实际上都是有争议的。因此，它们都是进一步细分的候选者。所以，我们不能像二分查找那样进行。

我正在使用迭代方法。但是，下面的代码不起作用：

    Line2d GetClippedLine()
    
        Line2d clippingCandidate = this->line;

        std::vector<Line2d> lines = clippingCandidate.GetMidpointSubLines();

        while(lines[0] != lines[1])
        
            lines = clippingCandidate.GetMidpointSubLines();

            Line2d one = lines[0];
            Line2d two = lines[1]; 

            if(one.IsClippingCandidate(rectangle))
            
                clippingCandidate = one;
            
            if(two.IsClippingCandidate(rectangle))
                       
                clippingCandidate = two;
            

            if(one.IsVisible(rectangle))
            
                Coordinates2d::Draw(one, Yellow);
            
            if(two.IsVisible(rectangle))
            
                Coordinates2d::Draw(two, Yellow);
            

            clippingCandidate.Show();
            //std::cout<<"++";
            //two.Show();
            std::cout<<"\n";
        

        return clippingCandidate;
    

【问题讨论】：

细分算法很容易通过递归实现。你接受递归解决方案吗？或者你只想要迭代版本？另一个问题：请描述剪辑区域的属性。是凸的吗？它总是一个矩形吗？对于矩形，细分方法没有意义，最好用显式公式实现裁剪。

我只是在学习CG。我需要一个迭代来更好地理解它。至于目前的问题：裁剪区域是一个凸矩形。如果需要凹面或其他多边形，我稍后会考虑。目前，我需要开始做一些事情。 “对于矩形，细分方法没有意义” - 可能是。但是，我需要它来学习算法。

【参考方案1】：

您只需细分任何不完全位于矩形内部或外部的线段。

我不明白这个算法有什么好处。至少看看你的代码，只计算与矩形的交点会快得多。

【讨论】：

这称为中点细分线裁剪算法。我正在阅读这本书并尝试实施它。

【参考方案2】：

从事几何建模多年，我一直以递归方式实现所有细分算法。它的编码和理解要简单得多，而且我很确定它不会比迭代解决方案慢。

一般情况

如果剪裁是凹的，那么答案可能由任意多个片段组成。在这种情况下，您需要一个堆栈（从递归中替换代码堆栈）。下面是一个实现的例子：

vector<Line> result;  //ordered array of line segments = clipping result
stack<Line> pieces;  //stack of pieces that are not processed yet
pieces.push(inputLine);

while (!pieces.empty()) 
    auto curr = pieces.top();
    pieces.pop();
    if (GetLength(curr) < eps)  //hard stop criterion
      continue;
    auto relative = GetSegmentPositionRelativeToRegion(curr, clippingRegion);
    if (relative == FULLY_INSIDE)
        result.push_back(curr);
    else if (relative == INTERSECTS)  //not fully inside, not fully outside
        Line left, right;   //halves of the curr line segment
        SplitAtMiddle(curr, left, right);
        pieces.push(left);
        pieces.push(right);
    

凸案例

如果你的剪切区域保证是凸的，那么结果中就只有一个线段。在这种情况下，细分过程具有非常简单的结构。假设输入段不是完全在剪辑区域之内或之外（即它相交）。

首先，您将片段细分，使一半在外部，而另一半与剪辑相交。这意味着您始终只能记住一个工作线段。

在某些时候这条规则被打破了，然后有两个细分分支。在左分支可能有两种分段细分的情况：(intersects + inside) 或 (outside + intersects) >）。在右侧反映：(inside + intersects) 或 (intersects + outside)。在所有情况下，只有一个子段相交，因此只剩下一个工作段。您可以为循环中的每个分支编写单独的代码。

附：当然，这种描述忽略了剪裁区域恰好通过中间点的特殊情况。至于我，实现这个算法是不值得的。

【参考方案3】：

你问的完全正确。 Explanations of midpoint subdivision have arisen that are very sloppy or just wrong。您的代码似乎基于这些不良来源之一。

M-S 仅适用于查找交点当您已经知道线段跨越剪切边界（每边一个端点）时，并且它通常使用整数实现。它最初被用作 Cohen 和 Sutherland 的完整裁剪算法的变体中的子程序。

如果您不熟悉 C-S，请参阅 the Wikipedia article。 “输出代码”引导对包含视口边界的无限线进行连续裁剪。在那里的伪代码中，您将用 M-S 替换浮点数学。

假设您在 x=C 处对左边界进行裁剪，跨越它的线段是 P0(x0,y0)---P1(x1,y1)。还说x0<C<=x1，所以P0 已知在边界之外。那么M-S算法就是：

tx1 = x1; // don't modify P1; it's inside the boundary
ty1 = y1;
while (x0 < C) 
  xm = (x0 + tx1 + 1) >> 1;
  ym = (y0 + ty1 + 1) >> 1;
  if (xm <= C)  // the midpoint is on or outside the boundary
    x0 = xm;     //   move P0
    y0 = ym;
   else        // the midpoint is strictly inside
    tx1 = xm;    //   move P1 
    ty1 = ym;
  

// The clipped segment is (x0,x1)--(y0,y1).

对于其他 3 个边界，您需要 3 个其他较小的变化。

终止条件很棘手。 + 1s 是必要的，以避免在 x0 = C-1 和 tx1 = C: (C + C - 1 + 1) >> 1 == C 的情况下永远循环，因此下一次迭代将终止。

话虽如此，中点细分已经过时了。它对只有整数算术的处理器很有用（至少在 90 年代中期之前都是这种情况；我在 1984 年用 8088 汇编语言实现了它）。找到中点只需要除以 2，这是整数右移，因此对于最大大小为 n 的坐标，可以使用不超过上限（log_2 n）次快速迭代进行剪辑。如今，浮点单元以 gigaflop 速率运行，使用浮点进行剪辑可能更快，当然也更容易。

加法 只是为了好玩，用 C 实现：

#include <stdio.h>
#include <stdlib.h>

typedef unsigned OUTCODE;
typedef int COORD;
typedef int BOOL;
#define TRUE 1
#define FALSE 0

#define XMIN 0
#define YMIN 0
#define XMAX 5000
#define YMAX 3000

// Not strictly portable, but usually fine.
#define SIGN_BIT (~(~0u >> 1))

#define LEFT SIGN_BIT
#define TOP (LEFT >> 1)
#define RIGHT (TOP >> 1)
#define BOTTOM (RIGHT >> 1)
#define ALL (LEFT | BOTTOM | RIGHT | TOP)

// Mask the sign bit.
#define M(X) ((X) & SIGN_BIT)

// Shift previous value and mask in the new sign bit.
#define SM(Prev, New) (((OUTCODE)(Prev) >> 1) | M(New))

__inline OUTCODE outcode(COORD x, COORD y) 
  return SM(SM(SM(M(YMAX - y), XMAX - x), y - YMIN), x - XMIN);


// In the S-T coordinate system, pO is outside boundary C and will be moved
// to the boundary while pI doesn't move. I is the termination correction.
#define MOVE_TO_BOUNDARY(SO, TO, SI, TI, C, I, IS_OUTSIDE) do  \
  COORD tsi = SI, tti = TI; \
  while (SO IS_OUTSIDE C)  \
    COORD sm = (tsi + SO + I) >> 1; \
    COORD tm = (tti + TO + I) >> 1; \
    if (sm IS_OUTSIDE ## = C)  \
      SO = sm; \
      TO = tm; \
     else  \
      tsi = sm; \
      tti = tm; \
     \
   \
 while (0)

BOOL clip(COORD *x0p, COORD *y0p, COORD *x1p, COORD *y1p) 
  COORD x0 = *x0p, y0 = *y0p, x1 = *x1p, y1 = *y1p;
  OUTCODE code0 = outcode(x0, y0);
  OUTCODE code1 = outcode(x1, y1);
  for (;;) 
    if ((code0 | code1) == 0) 
      *x0p = x0; *y0p = y0; *x1p = x1; *y1p = y1;
      return TRUE;
     else if (code0 & code1) 
      return FALSE;
     else if (code0) 
      if (code0 & BOTTOM)     MOVE_TO_BOUNDARY(y0, x0, y1, x1, YMAX, 0, >);
      else if (code0 & RIGHT) MOVE_TO_BOUNDARY(x0, y0, x1, y1, XMAX, 0, >);
      else if (code0 & TOP)   MOVE_TO_BOUNDARY(y0, x0, y1, x1, YMIN, 1, <);
      else /* LEFT */         MOVE_TO_BOUNDARY(x0, y0, x1, y1, XMIN, 1, <);
      code0 = outcode(x0, y0);
     else 
      if (code1 & BOTTOM)     MOVE_TO_BOUNDARY(y1, x1, y0, x0, YMAX, 0, >);
      else if (code1 & RIGHT) MOVE_TO_BOUNDARY(x1, y1, x0, y0, XMAX, 0, >);
      else if (code1 & TOP)   MOVE_TO_BOUNDARY(y1, x1, y0, x0, YMIN, 1, <);
      else /* LEFT */         MOVE_TO_BOUNDARY(x1, y1, x0, y0, XMIN, 1, <);
      code1 = outcode(x1, y1);
    
  


int main(void) 
  int n = 0, margin = 2000;
  for (;;) 
    // Generate some random points around the viewport.
    int x0 = rand() % (2 * margin + XMAX - XMIN) - margin;
    int y0 = rand() % (2 * margin + YMAX - YMIN) - margin;
    int x1 = rand() % (2 * margin + XMAX - XMIN) - margin;
    int y1 = rand() % (2 * margin + YMAX - YMIN) - margin;
    printf("a(%d, %d)--(%d, %d) %x--%x\n", x0, y0, x1, y1,
      outcode(x0,y0) >> 28, outcode(x1,y1) >> 28);
    BOOL r = clip(&x0, &y0, &x1, &y1);
    printf("a(%d, %d)--(%d, %d): %d\n", x0, y0, x1, y1, r);
  
  return 0;

在我的 MacBook 上，它在 90 秒内剪辑了十亿个片段。看看这与浮点 C-S 相比如何会很有趣。

【讨论】：

老实说，在我看来，这与通常的二分搜索一模一样 =)

我们正在计算 2D 中的线段与框相交（得到线段作为结果）。如果我们将 2D 点存储为双精度对，那么我们现在可以在 SSE2 中轻松地做到这一点，例如 here。大多数时间将花费在计算坐标增量的倒数上。在 Sandy Bridge 上需要 22 个周期，所有其他操作每个只需要 1 个周期。所以我觉得35个周期左右是真的，也就是说单核可以在10秒内处理10^9个段。

但是，我们可以以单精度存储坐标，一次处理两个段。用approximation + NR 可以在 4 个周期内计算 4 个倒数。我想有可能在 3-4 秒内完成 =)

@stgatilov 非常酷。我有点惊讶差异甚至更大。但是 MSD 的时间取决于输入的大小，因为这设置了 # 次迭代的上限。我的测试仅使用 14 个有效位。这在当时是一个不错的算法。

【参考方案4】：

我认为这是硬件实现的原因。不要将其视为剪辑。将其视为如何将线（或后来投影的三角形）排序为四叉树的练习。最后，您可以继续裁剪到半节点细分，直到每个节点都只有一个像素宽。不要将您的外部部件视为特殊部件。它们只是非常大的四边形。