如何生成等间距插值数值

13

我有一组(x,y)值,这些值不是等距的。 此处 是本问题中使用的归档。

我能够在这些值之间进行插值,但我得到的插值点不是等间距的。我所做的是:

x_data = [0.613,0.615,0.615,...]
y_data = [5.919,5.349,5.413,...]

# Interpolate values for x and y.
t = np.linspace(0, 1, len(x_data))
t2 = np.linspace(0, 1, 100)
# One-dimensional linear interpolation.
x2 = np.interp(t2, t, x_data)
y2 = np.interp(t2, t, y_data)

# Plot x,y data.
plt.scatter(x_data, y_data, marker='o', color='k', s=40, lw=0.)

# Plot interpolated points.
plt.scatter(x2, y2, marker='o', color='r', s=10, lw=0.5)

这导致结果如下:
如图所示,在原始点分布更密的图形部分,红点彼此更接近。
我需要一种方法来生成插值点,使其在x、y方向上 等距,步长给定为0.1。
正如 askewchan 正确指出的那样,当我说"等距在x,y时,我的意思是曲线中两个相邻的插值点之间应该有相同的距离(欧几里得直线距离)。
我尝试了unubtu的答案,对于平滑曲线效果很好,但对于不太平滑的曲线似乎会出现问题:
发生这种情况是因为代码以欧几里得方式计算点距离,而不是直接在曲线上进行计算,我需要距离曲线上的点之间的距离相同。有什么方法可以解决这个问题吗?

1
equispaced in x, y” 是指曲线在切线方向上等距离分布吗? - askewchan
@askewchan 抱歉,我的表达含义不太明确。我会重构问题,让它更加清晰明了。 - Gabriel
我不确定我理解你的意思@askewchan,你能详细说明一下你的评论吗? - Gabriel
2
你需要将插值建立为一个名为 y(x) 的函数,并在一些非线性间隔的数组 x 上调用它,但步长(dx)应该是所需步长(ds)除以根号下 (1 + 斜率平方)。这是因为曲线段的弧长长度为 ds = sqrt(1 + slope^2)*dx (使用勾股定理 ds^2 = dx^2 + dy^2)。 - askewchan
1
"...应该相互保持距离...", 你如何测量距离?直线还是沿曲线的距离? - elyase
显示剩余2条评论
5个回答

18

将你的xy数据转换为参数化曲线,即计算点之间的所有距离,并通过累积求和生成曲线上的坐标。然后独立地相对于新坐标插值x和y坐标。

import numpy as np
from matplotlib import pyplot as plt

data = '''0.615   5.349
    0.615   5.413
    0.617   6.674
    0.617   6.616
    0.63    7.418
    0.642   7.809
    0.648   8.04
    0.673   8.789
    0.695   9.45
    0.712   9.825
    0.734   10.265
    0.748   10.516
    0.764   10.782
    0.775   10.979
    0.783   11.1
    0.808   11.479
    0.849   11.951
    0.899   12.295
    0.951   12.537
    0.972   12.675
    1.038   12.937
    1.098   13.173
    1.162   13.464
    1.228   13.789
    1.294   14.126
    1.363   14.518
    1.441   14.969
    1.545   15.538
    1.64    16.071
    1.765   16.7
    1.904   17.484
    2.027   18.36
    2.123   19.235
    2.149   19.655
    2.172   20.096
    2.198   20.528
    2.221   20.945
    2.265   21.352
    2.312   21.76
    2.365   22.228
    2.401   22.836
    2.477   23.804'''

data = np.array([line.split() for line in data.split('\n')],dtype=float)

x,y = data.T
xd = np.diff(x)
yd = np.diff(y)
dist = np.sqrt(xd**2+yd**2)
u = np.cumsum(dist)
u = np.hstack([[0],u])

t = np.linspace(0,u.max(),10)
xn = np.interp(t, u, x)
yn = np.interp(t, u, y)

f = plt.figure()
ax = f.add_subplot(111)
ax.set_aspect('equal')
ax.plot(x,y,'o', alpha=0.3)
ax.plot(xn,yn,'ro', markersize=8)
ax.set_xlim(0,5)

plot generated by code


在我的情况下,它在这一行给了我 ValueError: object too deep for desired array - xn = numpy.interp(t, u, x) - Agniva De Sarker
我刚刚运行了代码,它正常工作。numpy的版本是1.10.1,scipy的版本是0.15.1。 - Christian K.
1
谢谢这个,这正是我也在寻找的。您有什么想法可以将其适应于更高维度的数据(也许使用scipy.interpolate.interpn)?我不想循环每个维度并独立执行np.interp - Maghoumi
1
这条评论发表后仅10分钟,我就找到了解决方案 :). 已在此处发布为答案。 - Maghoumi

10

首先考虑一个简单的情况。假设你的数据看起来像下面的蓝线。

enter image description here

如果你想选择等距离的点,这些点之间的距离为r,那么会有一些关键值r,其中在(1,2)处的尖点是第一个等距离点。

如果你想要大于这个关键距离的点,则第一个等距离点将从(1,2)跳到非常不同的位置——由绿色弧与蓝线相交所示。这种变化并不是逐渐的。

这个玩具案例表明,参数r的微小变化可能对解决方案产生根本性、不连续的影响。

它还表明,在确定(i+1)-th等距离点的位置之前,你必须知道第i个等距离点的位置。

因此,似乎需要一个迭代解决方案:

import numpy as np
import matplotlib.pyplot as plt
import math

x, y = np.genfromtxt('data', unpack=True, skip_header=1)
# find lots of points on the piecewise linear curve defined by x and y
M = 1000
t = np.linspace(0, len(x), M)
x = np.interp(t, np.arange(len(x)), x)
y = np.interp(t, np.arange(len(y)), y)
tol = 1.5
i, idx = 0, [0]
while i < len(x):
    total_dist = 0
    for j in range(i+1, len(x)):
        total_dist += math.sqrt((x[j]-x[j-1])**2 + (y[j]-y[j-1])**2)
        if total_dist > tol:
            idx.append(j)
            break
    i = j+1

xn = x[idx]
yn = y[idx]
fig, ax = plt.subplots()
ax.plot(x, y, '-')
ax.scatter(xn, yn, s=50)
ax.set_aspect('equal')
plt.show()

这里输入图片描述

注意:我将宽高比设置为“equal”以更清楚地显示点是等距离的。


这是一个很好的答案,抱歉我回复得这么晚(被另一个问题卡住了)。这段代码适用于相对平滑的曲线,但似乎对不太规则的曲线无法奏效。请参见编辑后的问题。 - Gabriel
1
@Gabriel:我已经编辑了我的答案。基本上是一样的——你只需要在跳转点之间跟踪距离的累加和。 - unutbu
非常感谢 @unubtu,现在这个答案完美了! - Gabriel

2

在 @Christian K. 的回答基础上,这里介绍如何使用 scipy.interpolate.interpn 对高维数据进行重采样。假设我们想要将数据重采样为10个等间距点:

import numpy as np
import scipy
# Assuming that 'data' is rows x dims (where dims is the dimensionality)
diffs = data[1:, :] - data[:-1, :]
dist = np.linalg.norm(diffs, axis=1)
u = np.cumsum(dist)
u = np.hstack([[0], u])
t = np.linspace(0, u[-1], 10)
resampled = scipy.interpolate.interpn((u,), pts, t)

2
以下脚本将使用等步长x_max - x_min / len(x) = 0.04438插值点。
import numpy as np
from scipy.interpolate import interp1d
import matplotlib.pyplot as plt

data = np.loadtxt('data.txt')
x = data[:,0]
y = data[:,1]

f = interp1d(x, y)
x_new = np.linspace(np.min(x), np.max(x), x.shape[0])
y_new = f(x_new)

plt.plot(x,y,'o', x_new, y_new, '*r')
plt.show()

enter image description here


这个解决方案在x轴上产生等间距点,但不是在y轴上,因此曲线上不是等距的。无论如何,还是谢谢你的答案。 - Gabriel

0

在曲线上生成等距离点是可能的。但是需要更明确地定义您想要什么样的实际答案。很抱歉,我为此任务编写的代码是用MATLAB编写的,但我可以描述一般思路。有三种可能性。

首先,点是否真正与邻居在简单欧几里得距离方面等距离?要实现这一点,需要找到曲线上任意点与固定半径圆的交点。然后沿着曲线步进。

接下来,如果您希望距离是沿着曲线本身的距离,如果曲线是分段线性的,则问题再次容易解决。只需沿着曲线步进,因为测量线段上的距离很容易。

最后,如果您希望曲线是一个三次样条,这也不是非常困难,但需要更多的工作。这里的诀窍是:

  • 计算曲线上从点到点的分段线性弧长。将其称为t。
  • 生成一对三次样条,x(t),y(t)。
  • 将x和y作为t的函数进行微分。由于这些是三次线段,因此很容易。导数函数将是分段二次的。
  • 使用ode求解器沿着曲线移动,积分微分弧长函数。在MATLAB中,ODE45效果很好。

因此,要进行积分:

sqrt((x')^2 + (y')^2)

在MATLAB中,ODE45可以被设置为识别函数穿过某些指定点的位置。

如果你的MATLAB技能足够强,你可以查看interparc中的代码以获得更多解释。这是相当好注释的代码。


网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接