在任意X坐标处确定B样条的值

Question

我已经在数据点xp和yp上创建了一个开放、夹紧、立方B样条。

样条是由跨越xp域的向量u参数化的。

我的目标是确定在xp域中给定“x”坐标处的B样条的“y”坐标。

当我在计算tck之后将值“4”传递给splev时，与参数4对应的x和y坐标的值都会返回，这是生成参数曲线时预期的行为。

我能够使用牛顿法确定给定“x”坐标的参数u的值；然而，这是间接的，并且需要比我的最终应用程序允许的更多的计算时间。

有人能建议一种更直接的方法来确定给定“x”的B样条上的“y”坐标吗？

import numpy as np
import matplotlib.pyplot as plt
from scipy import interpolate

xp = [0., 0.71428571, 1.42857143, 2.14285714, 2.85714286, 3.57142857, 4.28571429, 5.]
yp = [0., -0.86217009, -2.4457478, -2.19354839, -2.32844575, -0.48680352, -0.41055718, -3.]

length = len(xp)
t = np.linspace(0., xp[-1], length - 2, endpoint=True)
t = np.append([0, 0, 0], t)
t = np.append(t, [xp[-1], xp[-1], xp[-1]])

tck = [t, [xp, yp], 3]
u = np.linspace(0, 5., 1000, endpoint=True)
out = interpolate.splev(u, tck)

x_value_in_xp_domain = 4.
y_value_out = interpolate.splev(x_value_in_xp_domain, tck)

plt.plot(xp, yp, linestyle='--', marker='o', color='purple')
plt.plot(out[0], out[1], color = 'teal')
plt.plot(x_value_in_xp_domain, y_value_out[1], marker='o', color = 'orangered')
plt.plot(y_value_out[0], y_value_out[1], marker='o', color = 'black')
plt.axvline(x=x_value_in_xp_domain, color = 'orangered')
plt.show()

下面的图像显示了指导多边形和由上述代码生成的B样条曲线。在x=4处的橙色点对应于我希望直接确定B样条曲线上y值的点。黑色点是将4作为参数传递时B样条曲线的值。

提供一些有用的参考资料：

使用numpy/scipy实现快速b样条算法

https://github.com/kawache/Python-B-spline-examples

https://pages.mtu.edu/~shene/COURSES/cs3621/NOTES/spline/B-spline/bspline-curve.html

http://web.mit.edu/hyperbook/Patrikalakis-Maekawa-Cho/node17.html

- xemx

---

一般来说，样条曲线在给定X值时可能有多个Y值，这会使直接计算变得复杂。是否有更好的表示曲线的方法？ - Martin Stone

@MartinStone 很正确。在这个应用程序中，对于给定的X值，样条曲线保证是单值的。关于应用程序的更多背景，一组样条曲线将被用来表示解决方案，然后传递给优化器。在此解决方案空间中，端点是固定的，因此开放、夹紧、立方、B 样条似乎是表示任意函数的好方法，同时避免了朗格现象。虽然没有反对提出其他方法来近似任意函数，但仍然对原始问题的解决方案感兴趣。 - xemx

Answer 1

让我建议一个框架变更。如你所知，在平面上创建样条曲线有两种方式：作为函数和作为参数曲线。你选择了后者，但在问题中我没有看到任何强制要求这样做的理由；而选择前者将使得解决你的问题几乎是轻而易举的。

我们从一些输入数据开始。

xdata  = [0.,  0.60626102,  1.36904762,  2.14285714,  2.85714286,  3.63095238,  4.39373898,  5.]
ydata  = [0., -0.82509685, -2.00782014, -2.25806452, -1.99902249, -0.77468231, -1.19975743, -3.]

请注意，与您的xp和yp不同，在我的圈子里被称为控制点，我收集了一些点xdata，ydata，我们将进行插值。我通过在Greville abscissae处对您的函数进行采样来实现这一点，以获得类似于您的输入的东西。然而，确切类型的输入取决于您的应用程序；您在评论中提到了函数逼近，所以我想这可能是一种方法。

现在，我们可以使用UnivariateSpline创建插值数据的样条函数：

spline = interpolate.InterpolatedUnivariateSpline(xdata, ydata)

有一些替代方案；我建议看一下插值教程和其中的参考资料。

大致上就是这样了。您现在可以在许多点上取样值来绘制样条曲线：

xfunc  = np.linspace(0, 5., 1000, endpoint=True)
yfunc  = spline(xfunc)
plt.plot(xfunc, yfunc,  color = 'green')

然而，我们之所以要经过所有的这一切，是为了回答你的问题：如何找到对应于给定 x 坐标的 y 坐标？ 这很简单：

x_value_in_xp_domain = 4.
y_value_out = spline(x_value_in_xp_domain)

这里有一张图片展示了我的输入数据（蓝色；请注意，它们与您的输入不完全相同），插值样条函数（绿色），x_value_in_xp_domain（橙色曲线）和y_value_out（橙色点）。

最后，这里是完整的源代码。

import numpy as np
import matplotlib.pyplot as plt
from scipy import interpolate

xdata  = [0.,  0.60626102,  1.36904762,  2.14285714,  2.85714286,  3.63095238,  4.39373898,  5.]
ydata  = [0., -0.82509685, -2.00782014, -2.25806452, -1.99902249, -0.77468231, -1.19975743, -3.]

spline = interpolate.InterpolatedUnivariateSpline(xdata, ydata)
xfunc  = np.linspace(0, 5., 1000, endpoint=True)
yfunc  = spline(xfunc)

x_value_in_xp_domain = 4.
y_value_out = spline(x_value_in_xp_domain)

plt.plot(xfunc, yfunc,  color = 'green')
plt.plot(xdata, ydata, 'o', color = 'blue')
plt.plot(x_value_in_xp_domain, y_value_out, marker='o', color = 'orangered')
plt.axvline(x=x_value_in_xp_domain, color = 'orangered')
plt.show()

- Dominik Mokriš

Answer 2

< p > PPoly 类的 roots 方法可能是您问题的解决方案。这比牛顿法更快，并且如果有多个解决方案，它将为您提供所有解决方案。

sx = interpolate.BSpline(t, xp, 3)
sy = interpolate.BSpline(t, yp, 3)

x0 = 4
u0 = interpolate.PPoly.from_spline((sx.t, sx.c - x0, 3)).roots()
sy(u0)

plt.plot(xp, yp, linestyle='--', marker='o', color='purple')
plt.plot(out[0], out[1], color = 'teal')
plt.plot(sx(u0), sy(u0), 'o')
plt.show()

- UNAI SAN MIGUEL