17

我在互联网上四处寻找，找不到解决这个特定问题的完美算法：

我们的客户有一组点和重量数据以及每个点，如下图所示：

加权点 http://chakrit.net/files/stackoverflow/so_heightmap_points.png

其中，我们有一个 GIS 程序，可以从这些点及其权重值生成“高度图”或一种地形数据，但由于我们有近一千个数据点并且这些数据会随着时间而变化，我们希望创建我们自己的工具来自动生成这些高度图。

到目前为止，我已经尝试计算每个像素从其到最近数据点的距离的Sqrt((x1 - x2) ^ 2 + (y1 - y2) ^ 2)权重，并将权重和距离因子应用于数据点的颜色以生成该特定像素的结果渐变颜色：

高度图结果 http://chakrit.net/files/stackoverflow/so_heightmap_result.png

您可以看到某些数据点的配置仍然存在问题，并且当有很多数据点时，该算法有时会产生一个相当多边形的图像。理想的结果应该看起来更像一个省略号，而不是一个多边形。

这是来自维基百科关于梯度上升的文章的一个示例图像，它展示了我想要的结果：

山 http://chakrit.net/files/stackoverflow/so_gradient_descent.png

梯度上升算法不是我感兴趣的。我感兴趣的东西；是首先计算该图片中的原始函数的算法，提供具有权重的数据点。

我没有上过任何拓扑数学课，但我可以做一些微积分。我想我可能遗漏了一些东西，并且对我应该在那个谷歌搜索框中输入什么感到迷茫。

我需要一些指示。

谢谢！

4

9 回答 9

7

您正在寻找的是表面插值。

有些产品可以做到这一点（这里有一个）

然后可以以所需的分辨率询问生成的函数/样条/其他数学结构以提供高度图。

你的插值函数

Sqrt((x1 - x2) ^ 2 + (y1 - y2) ^ 2)

类似于反距离加权 方法，只是您应用了任意过滤器并丢弃了许多其他数据点。

这些技术中的大多数都依赖于合理数量的样本和支撑这些值的“类似地形”的行为。

我建议使用权重作为高度样本，并在第二个链接中尝试简单的 Shepard 方法（不要过滤任何像素开始），方法是在可以混合的插值点处获取样本点对整体高度值的贡献比例这些比率中的样本颜色也可以为点着色。使用强度（粗略地说是简单 RGB 空间中的灰度）来显示高度或像示例图像一样添加黑色轮廓线。

于 2009-02-17T16:14:52.950 回答
5

这个问题并不像表面上看起来那么容易。您的问题是两个区域的边界两侧需要具有相同的高度，也就是说，给定像素的高度由多个最近邻确定。

如果我理解正确，您至少需要两种算法（和第三条行话）。

要正确执行此操作，您需要将平面分解为Voronoi 镶嵌。

您可能会想要使用kd-tree来帮助您找到最近的邻居。而不是采用 O(n^2)，而是将其降低到 O(n log(n))（额外的好处是您的 Voronoi 区域生成阶段在开发中将足够快，可以在高度计算阶段工作）。

既然您有一个将每个点索引到其最近邻居 i 的二维地图，您需要遍历地图上的每个 x,y 点并计算其高度。

要对给定的点 x,y 执行此操作，首先获取其最近的邻居 i 并将其粘贴到列表中，然后收集 Voronoi 图上的所有连续区域。一个简单的方法是使用洪水填充找到该区域中的所有点，然后查看边界并收集其他身份。

使用所有最近邻居的列表，您现在可以正确插值！（有关插值方案，请参见其他答案）。

于 2009-02-17T16:17:20.293 回答
4

您已询问有关不规则数据的二维插值算法的信息，这是一个相当复杂的领域。既然你说你有 ArcGIS，我强烈建议你在 ArcGIS 中使用其内置功能进行自动插值以进行自动计算。我相信这比编写自己的插值算法要容易得多。我已经完成了 ArcGIS 的一些自动化操作，它相当简单。

如果你确实编写了自己的插值代码——我建议你不要——首先要选择合适的算法，因为有几个算法都有自己的优缺点。这是从优秀的插值工具Surfer的帮助中摘录的一些建议（顺便说一句，它也可以很容易地自动化）。还有更多算法，这些只是我尝试过的。

  • 克里金法是更灵活的方法之一，可用于对几乎任何类型的数据集进行网格化。对于大多数数据集，使用默认线性变差函数的克里金法非常有效。一般来说，我们最常推荐这种方法。克里金法是默认的网格化方法，因为它可以为大多数数据集生成良好的地图。对于较大的数据集，克里金法可能会相当慢。克里金法可以推断超出数据 Z 范围的网格值。
  • 逆距离加权速度很快，但倾向于在数据点周围生成同心轮廓的“靶心”模式。幂的反距离不会外推超出数据范围的 Z 值。一个简单的逆距离加权算法很容易实现，但速度较慢。
  • 使用线性插值进行三角剖分很快。当您使用小型数据集时，带有线性插值的三角剖分会在数据点之间生成不同的三角形面。使用线性插值进行三角剖分不会将 Z 值外推到数据范围之外。
  • Shephard 方法类似于 Inverse Distance to a Power，但不会产生“靶心”模式，尤其是在使用平滑因子时。Shepard 方法可以推断超出数据 Z 范围的值。

要实现算法：您可以尝试谷歌搜索，或点击其他一些答案中的链接。有一些包含插值的开源 GIS 包，所以如果你喜欢通过 C++ 挖坑，也许你可以从中提取算法。或者大卫沃森的这本书显然被认为是经典之作，虽然它是一个棘手的阅读和示例代码是意大利面条基本！但是，据我所知，这是最好的。如果 Stack Overflow 上的其他人知道得更好，请纠正我，因为我也不敢相信。

于 2009-02-21T20:42:36.053 回答
3

克里金法是执行此操作的重量级方法之一，尤其是在 GIS 领域。它有几个很好的数学特性 - 缺点是它可能会很慢，具体取决于您的变异函数。

如果你想要更简单的东西，有很多插值程序可以很好地处理这个问题。如果你能找到一份《数值配方》，第 3 章将专门解释许多插值的变体，并包含代码示例及其功能属性的描述。

于 2009-02-17T17:02:36.067 回答
2

是首先计算该图片中的原始函数的算法，提供具有权重的数据点。

这是可能的。如果您从单点开始，您将始终以圆形结束，但如果您对数据点进行加权并考虑到这一点，您可以将圆形挤压成椭圆形，如图所示。

您最终使用多边形的原因是您在计算中使用了离散函数 - 首先找到最接近的颜色，然后确定颜色。

相反，您应该研究梯度算法，该算法根据将点包围在三角形中的三个数据点的距离和权重为点分配颜色。

梯度算法

这取决于您要显示的内容。一个简单的算法是：

对于每个像素：

  • 找到形成围绕该像素的最小三角形的三个点

  • 将此点设置为受权重和到每个数据点的距离影响的颜色（HSV 颜色系统）：

    pixel.color = datapoint[1].weight * distance(pixel, datapoint[1]) * datapoint[1].color + datapoint[2].weight * distance(pixel, datapoint[2]) * datapoint[2].color + datapoint[3].weight * distance(pixel, datapoint[3]) * datapoint[3].color

我在这里使用 +，但您需要确定适合您的应用程序的“平均”算法。

-亚当

于 2009-02-17T16:13:24.920 回答
1

表面插值似乎是一个困难的数学问题。另一种更便宜的方法是：

For each pixel:
For each point:
pixel.addWeight(weight(point, pixel))

def addWeight(w):
totalweight += w
numberofweights += 1
weight = totalweight / numberofweights

示例权重函数：

def weight(point, pixel):
return point.weight * 1/(1 + sqrt((point.x - pixel.x)^2 + (point.y - pixel.y)^2))

这是一种蛮力的方法，但它很简单。

于 2009-02-17T16:46:57.680 回答
1

不久前，我在 Winamp AVS 中实现了类似的东西。它使用“metaballs”类型的方法计算每个数据点的平方反比距离（以避免速度的 sqrt），将其限制（例如到 1.0），并为 2D 网格上的每个点计算这些距离的总和。这将给出一个平滑变化的颜色/高度图。

如果您想查看代码，它位于我的J10 AVS 包中的“Glowy”预设中。

编辑：只是看着它，我添加了一些其他爵士乐以使其看起来更漂亮，最重要的部分是：

d1=s/(sqr(px1-rx)+sqr(py1-ry));
d2=s/(sqr(px2-rx)+sqr(py2-ry));
d3=s/(sqr(px3-rx)+sqr(py3-ry));
d4=s/(sqr(px4-rx)+sqr(py4-ry));
d5=s/(sqr(px5-rx)+sqr(py5-ry));
d6=s/(sqr(px6-rx)+sqr(py6-ry));
d=d1+d2+d3+d4+d5+d6;

这需要 6 点的总和。对红色、绿色和蓝色输出值所做的一切都是为了让它看起来更漂亮。6 分不算多，但请记住，当它是新机器时（它以 ~20fps 的速度运行），我试图在 400MHz 机器上的 320x200 网格上实时运行它。:)

用 red = d 替换 red =、green = 和 blue = ... 行；等等……看看我的意思。所有的美丽都消失了，你只剩下一个灰度图像，数据点周围有平滑变化的斑点。

另一个编辑：我忘了说“s”是所有点的共享权重，为每个点更改它会赋予每个点单独的权重，例如 d1 = 2/(...) 和 d2 = 1/(.. .) d1 的中心高度是 d2 的两倍。您可能还想用 d1 = 2/max(..., 1.0) 之类的东西来限制底部的表达式，以平滑点的顶部，这样它们就不会在中间的无穷大处达到峰值。:)

对不起答案的混乱......我认为发布代码示例就足够了，但在检查时我的代码令人困惑且难以阅读。:(

于 2009-05-25T10:05:49.793 回答
1

我知道这是一个相当古老的问题，但我在尝试解决类似问题时偶然发现了它。

有一个名为Surfit的开源项目正是实现了这种类型的功能。

于 2010-02-19T05:17:18.133 回答
0

您正在寻找Blender称为“元球”的东西（带有链接的维基百科文章，示例）。这样想：

你的物体是从地面伸出的圆锥体。它们都是抛物线，重量表明它们从地面伸出多远。或者，使它们都具有相同的高度并相应地调整抛物线的“平面度”，因此较大的重量会使圆锥体非常宽，而较低的重量会使圆锥体锋利。甚至在一定程度上两者兼而有之。

我建议您实现它并查看它的外观。

接下来，您需要在结果上挂一块布或橡胶布。布料会拉伸一定量，通常会由于重力而下垂。锥体保持它。

只要你靠近一个圆锥的中心，Z坐标就是圆锥表面上的位置。当你离开锥体中心时，重力开始下降，其他锥体的影响增加。

于 2009-02-17T16:34:22.310 回答