让我们问一下如何对多边形进行平滑处理以将所有顶点保留在其精确位置吗？
实际上，轮廓是在网格数据的顶部进行的，您可能会建议随后对网格数据进行平滑处理，因此生成的轮廓将更加平滑。请注意，这不符合我的期望，因为诸如高斯滤波器之类的平滑功能将删除少量数据包，并将更改第三个变量的范围，例如在我的应用程序中不允许的高度。
实际上，我正在寻找一段代码（最好在Python中），该代码可以对2D多边形（任何类型：凸形，凹形，自相交等）进行平滑处理（而无需担心代码）和准确。
FYI，ArcGIS中有一个功能可以完美地做到这一点，但是对于这个问题，我不是选择使用第三方商业应用程序。
1）
Scipy.interpolate：
/>如您所见，样条曲线（红色）不令人满意！
2）
这里是使用此处给出的代码得到的结果。运行不正常！
3）
对我来说，最好的解决方案应该是如下图所示，其中仅通过更改正方形即可逐渐平滑正方形一个价值。我希望有类似的概念可以平滑任何形式的多边形。
满足样条通过点的条件：
4）
这是我在他的数据中逐行在Python中实现“ whuber的想法”的方法。由于结果不好，可能存在一些错误。
K = 2是一场灾难，因此对于k> = 4。
5）
我在有问题的位置上删除了一个点，现在所得到的样条曲线与垫圈的样条曲线相同。但是仍然有一个问题，为什么该方法不适用于所有情况？
6）
很好的平滑wuber数据的方法如下（由矢量图形软件绘制），其中加号具有已顺利添加（与更新
4比较）：
7）
查看Python的结果一些标志性形状的Whuber代码的版本：
请注意，该方法似乎不适用于折线。对于角折线（轮廓），绿色是我想要的，但变成红色。这需要解决，因为等高线图总是多段线，尽管闭合的多段线可以像我的示例中那样视为多边形。同样不是尚未解决更新4中出现的问题。8）[我的最后]
这是最终解决方案（不完美！）：
请记住，您将不得不对恒星所指向的区域执行某些操作。我的代码中可能存在错误，或者建议的方法需要进一步开发以考虑所有情况并提供所需的输出。
#1 楼
样条数字序列的大多数方法将样条多边形。诀窍是使样条线在端点处平滑“闭合”。为此，将“顶点”环绕两端。然后分别对x和y坐标进行样条化。这是
R中的一个有效示例。它使用基本统计信息包中提供的默认立方spline过程。为了获得更多控制，请替换几乎所有您喜欢的过程：只需确保它能在数字上花键（即对它们进行插值），而不是仅将它们用作“控制点”。#
# Splining a polygon.
#
# The rows of 'xy' give coordinates of the boundary vertices, in order.
# 'vertices' is the number of spline vertices to create.
# (Not all are used: some are clipped from the ends.)
# 'k' is the number of points to wrap around the ends to obtain
# a smooth periodic spline.
#
# Returns an array of points.
#
spline.poly <- function(xy, vertices, k=3, ...) {
# Assert: xy is an n by 2 matrix with n >= k.
# Wrap k vertices around each end.
n <- dim(xy)[1]
if (k >= 1) {
data <- rbind(xy[(n-k+1):n,], xy, xy[1:k, ])
} else {
data <- xy
}
# Spline the x and y coordinates.
data.spline <- spline(1:(n+2*k), data[,1], n=vertices, ...)
x <- data.spline$x
x1 <- data.spline$y
x2 <- spline(1:(n+2*k), data[,2], n=vertices, ...)$y
# Retain only the middle part.
cbind(x1, x2)[k < x & x <= n+k, ]
}
为了说明其用法，让我们创建一个小（但很复杂）的多边形。
#
# Example polygon, randomly generated.
#
set.seed(17)
n.vertices <- 10
theta <- (runif(n.vertices) + 1:n.vertices - 1) * 2 * pi / n.vertices
r <- rgamma(n.vertices, shape=3)
xy <- cbind(cos(theta) * r, sin(theta) * r)
使用前面的代码对它进行样条化。若要使样条曲线更平滑，请从100增加顶点数；否则，增加顶点数。要使其不那么平滑，请减少顶点数量。
s <- spline.poly(xy, 100, k=3)
为查看结果，我们绘制了（a）以红色虚线表示的原始多边形，显示了它们之间的间隙第一个和最后一个顶点（即不闭合其边界折线）； （b）花键为灰色，再次显示其间隙。 （由于间隙太小，其端点用蓝点突出显示。）
plot(s, type="l", lwd=2, col="Gray")
lines(xy, col="Red", lty=2, lwd=2)
points(xy, col="Red", pch=19)
points(s, cex=0.8)
points(s[c(1,dim(s)[1]),], col="Blue", pch=19)
评论
好答案。有什么方法可以确保轮廓不会由于平滑而最终交叉？
– Kirk Kuykendall
2012年5月7日20:50
这是个好问题，@ Kirk。我不知道有任何方法可以保证这种平滑形式不交叉。 （实际上，我什至没有看到如何保证平滑多段线是非自相交的。但是，对于大多数轮廓来说，这并不是一个大问题。）为此，您需要返回原始轮廓DEM，而是首先使用更好的方法来计算轮廓。 （有更好的方法-人们已经知道了很长时间-但AFAIK某些最受欢迎的GIS并未使用它们。）
– hu
2012年5月7日在21:01
首先，我仍在努力以Python来实现您的答案，但仍未成功。其次，如果在正方形上应用方法会产生什么结果？您可以参考我在问题中提出的内容。
–开发人员
2012年5月8日，下午2:05
我接受了这个答案，因为它提供了很好的解决方案。即使它不是一个完美的解决方案，但它给了我一些解决方法，希望我能找到一种满足我在问题和评论中提到的要点的解决方案。您也可以考虑对问题[QC]进行评论，这有很多技巧。最后，我应该说，安装了可爱的Scipy软件包后，对python的翻译几乎是简单的。还应将Pablo在质量控制中的评论视为折线（即Bezier曲线）的可能解决方案。大家好运。
–开发人员
2012年5月10日上午11:22
看到您的答案，我后悔没有很好地照顾我的数学！！！
– Vinayan
2012年5月11日在11:51
#2 楼
我知道这是一篇过时的文章，但是它在Google上显示了我正在寻找的内容，因此我认为应该发布解决方案。我不认为这是2D曲线拟合练习，而是3D练习。通过将数据视为3D，我们可以确保曲线永不交叉，并且可以使用其他轮廓的信息来改进对当前轮廓的估计。
以下iPython提取使用了由SciPy。请注意，我绘制的z值并不重要，只要所有轮廓的高度都相等即可。
In [1]: %pylab inline
pylab.rcParams['figure.figsize'] = (10, 10)
Populating the interactive namespace from numpy and matplotlib
In [2]: import scipy.interpolate as si
xs = np.array([0.0, 0.0, 4.5, 4.5,
0.3, 1.5, 2.3, 3.8, 3.7, 2.3,
1.5, 2.2, 2.8, 2.2,
2.1, 2.2, 2.3])
ys = np.array([0.0, 3.0, 3.0, 0.0,
1.1, 2.3, 2.5, 2.3, 1.1, 0.5,
1.1, 2.1, 1.1, 0.8,
1.1, 1.3, 1.1])
zs = np.array([0, 0, 0, 0,
1, 1, 1, 1, 1, 1,
2, 2, 2, 2,
3, 3, 3])
pts = np.array([xs, ys]).transpose()
# set up a grid for us to resample onto
nx, ny = (100, 100)
xrange = np.linspace(np.min(xs[zs!=0])-0.1, np.max(xs[zs!=0])+0.1, nx)
yrange = np.linspace(np.min(ys[zs!=0])-0.1, np.max(ys[zs!=0])+0.1, ny)
xv, yv = np.meshgrid(xrange, yrange)
ptv = np.array([xv, yv]).transpose()
# interpolate over the grid
out = si.griddata(pts, zs, ptv, method='cubic').transpose()
def close(vals):
return np.concatenate((vals, [vals[0]]))
# plot the results
levels = [1, 2, 3]
plt.plot(close(xs[zs==1]), close(ys[zs==1]))
plt.plot(close(xs[zs==2]), close(ys[zs==2]))
plt.plot(close(xs[zs==3]), close(ys[zs==3]))
plt.contour(xrange, yrange, out, levels)
plt.show()
这里的结果看起来不是最好的，但是控制点很少，它们仍然是完全有效的。请注意，如何拉出绿色的拟合线以遵循更宽的蓝色轮廓。
评论
拟合的平滑曲线应尽可能靠近原始多边形/折线。
–开发人员
17年2月28日在18:23
#3 楼
我几乎完全按照您的要求编写了程序包...但是它是在Perl中使用的，并且在十年前：GD :: Polyline。它使用2D三次方Bezier曲线，可以“平滑”任意的Polygon或“ Polyline”（当时叫我的名字，现在通常称为“ LineString”）。该算法分两个步骤：给定在多边形中的点上，在每个点之间添加两个Bezier控制点；然后调用简单的算法对样条曲线进行分段逼近。
第二部分很简单；第一部分是一点艺术。这是见解：将“控制段”视为顶点N：
vN。控制段为三个共线点：[cNa, vN, cNb]。中心点是顶点。该控制段的斜率等于从顶点N-1到顶点N + 1的斜率。该段的左侧部分的长度是从顶点N-1到顶点N的长度的1/3，而该段的右侧部分的长度是从顶点N-1到顶点N + 1的长度的1/3。 br /> 如果原始曲线是四个顶点：
[v1, v2, v3, v4]，则每个顶点现在将获得形式为[c2a, v2, c2b]的控制线段。像这样将它们串在一起：[v1, c1b, c2a, v2, c2b, c3a, v3, c3b, c4a, v4]并一次将它们四个作为贝塞尔曲线的四个点：[v1, c1b, c2a, v2]，然后[v2, c2b, c3a, v3]，依此类推。由于[c2a, v2, c2b]是共线的，因此生成的曲线在每个顶点处都将是平滑的。所以这也满足您参数化曲线“紧密度”的要求：对于，请使用小于1/3的值一条“更紧”的曲线，一条较大的曲线以“更圆滑”的配合。无论哪种情况，生成的曲线总是穿过原始给定的点。
这导致平滑的曲线“外接”了原始多边形。我也有某种方法可以“刻画”一条平滑的曲线...但是我在CPAN代码中看不到它。
无论如何，目前我还没有可用的Python版本，也没有任何数字。但是...如果/当我将其移植到Python时，我一定会在这里发布。
评论
无法评估Perl代码，请添加图形以演示其工作方式。
–开发人员
17年2月28日在18:24
评论
您如何生成“多边形”轮廓？难道它们始终不是直线，因为与DEM边缘相交的轮廓永远不会闭合吗？我已经使用GRASS中的v.generalize函数对轮廓线进行了平滑处理，尽管结果非常不错，但是对于轮廓非常密集的地图可能需要一些时间。
@pistachionut您可能会认为轮廓线是折线。我在第一阶段正在寻找纯代码。如果不可用，请使用适用于Python的light软件包。
也许看看scipy.org/Cookbook/Interpolation，因为这听起来像是要花键
链接中的@Pablo Bezier曲线适用于折线。 Whuber的多边形几乎可以正常工作。因此他们可以共同解决这个问题。非常感谢您免费分享您的知识。