scipy Python：如何在网格上插值Angular ？

hmmo2u0o 于 2023-04-30 发布在 Python

关注(0)|答案(2)|浏览(118)

我有一个网格，里面有一些给定的数据。该数据由其Angular （从0到π）给出。在这个网格中，我有另一个较小的网格。
这可能看起来像这样：

现在我想在这个网格上插入Angular 。
我使用scipy.interpolate.griddata进行了尝试，结果很好。但是当Angular 从几乎0变到几乎π时就有问题了（因为中间是π/2。..）
这是结果，很容易看出哪里出了问题。

我该怎么处理这个问题呢？Thank you！：）
下面是代码复制：

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

ax = plt.subplot()
ax.set_aspect(1)

# Simulate some given data.
x, y = np.meshgrid(np.linspace(-10, 10, 20), np.linspace(-10, 10, 20))
data = np.arctan(y / 10) % np.pi
u = np.cos(data)
v = np.sin(data)

ax.quiver(x, y, u, v, headlength=0.01, headaxislength=0, pivot='middle', units='xy')

# Create a smaller grid within.
x1, y1 = np.meshgrid(np.linspace(-1, 5, 15), np.linspace(-6, 2, 20))
# ax.plot(x1, y1, '.', color='red', markersize=2)

# Interpolate data on grid.
interpolation = griddata((x.flatten(), y.flatten()), data.flatten(), (x1.flatten(), y1.flatten()))
u1 = np.cos(interpolation)
v1 = np.sin(interpolation)
ax.quiver(x1, y1, u1, v1, headlength=0.01, headaxislength=0, pivot='middle', units='xy',
          color='red', scale=3, width=0.03)

plt.show()

编辑：
感谢@bubble，有一种方法可以在插值之前调整给定的Angular ，以便结果符合要求。因此：
1.定义校正函数：

def RectifyData(data):
    for j in range(len(data)):
        step = data[j] - data[j - 1]
        if abs(step) > np.pi / 2:
            data[j] += np.pi * (2 * (step < 0) - 1)
    return data

1.插值如下：

interpolation = griddata((x.flatten(), y.flatten()),
                         RectifyData(data.flatten()),
                         (x1.flatten(), y1.flatten()))
u1 = np.cos(interpolation)
v1 = np.sin(interpolation)

scipy

来源：https://stackoverflow.com/questions/55858802/python-how-to-interpolate-angles-on-a-grid

2条答案

按热度按时间

yh2wf1be1#

我尝试直接插值cos(angle)和sin(angle)值，但这仍然会导致中断，导致错误的行方向。其主要思想在于减少中断，e。例如，[2.99,3.01, 0.05,0.06]应该转换为如下所示：[2.99, 3.01, pi+0.05, pi+0.06]。这是正确应用2D插值算法所需的。在下面的post中出现了几乎相同的问题。

def get_rectified_angles(u, v):
    angles = np.arcsin(v)
    inds = u < 0
    angles[inds] *= -1
# Direct approach of removing discontinues 
#     for j in range(len(angles[1:])):  
#         if abs(angles[j] - angles[j - 1]) > np.pi / 2:
#             sel = [abs(angles[j] + np.pi - angles[j - 1]), abs(angles[j] - np.pi - angles[j-1])]
#             if np.argmin(sel) == 0:
#                 angles[j] += np.pi
#             else:
#                 angles[j] -= np.pi
    return angles

ax.quiver(x, y, u, v, headlength=0.01, headaxislength=0, pivot='middle', units='xy')

# # Create a smaller grid within.
x1, y1 = np.meshgrid(np.linspace(-1, 5, 15), np.linspace(-6, 2, 20))

angles = get_rectified_angles(u.flatten(), v.flatten())

interpolation = griddata((x.flatten(), y.flatten()), angles, (x1.flatten(), y1.flatten()))
u1 = np.cos(interpolation)
v1 = np.sin(interpolation)
ax.quiver(x1, y1, u1, v1, headlength=0.01, headaxislength=0, pivot='middle', units='xy',
          color='red', scale=3, width=0.03)

也许，numpy.unwrap函数可以用来修复中断。在1d数据的情况下，numpy.interp具有关键字period来处理周期性数据。

赞(0）回复(0）举报 2023-04-30

osh3o9ms2#

导入库

# cmath is mathematical functions for complex numbers
import cmath 
import numpy as np
from scipy.interpolate import griddata

让我们使用此示例数据集

# example dataset
x = [1,2, 1, 2]
y = [1, 1, 2, 2]
xi = [1,   1.5,  2,   1.5]
yi = [1.5, 1,    1.5, 2]
angles = [30, 350, 270, 130]

4个蓝色点上的Angular 将在4个红色点上插值，如图所示：

插值

# initialize an empty list named cmplx to store values 
cmplx=[]

# represent angles with unit vectors using complex numbers 
# i.e. radius=1, angle=a (angle needs to be in radians).
for a in angles:
    cmplx.append(cmath.rect( 1, np.deg2rad(a) ))

# interpolate values to grid (xi, yi)
cmplx_nums = griddata((x, y), cmplx, (xi, yi))

# phase of the complex number is the angle in radians. 
# convert to degrees and use modulus operator to complement to 360 i.e % 360 
angles_2 =[]
for cmplx in cmplx_nums:
    angle_rad = cmath.phase(cmplx)
    angles_2.append(np.rad2deg(angle_rad) % 360)

打印结果

for x in angles_2:
print('{:.0f}'.format(x))

330
10
60
200

赞(0）回复(0）举报 2023-04-30

我来回答

scipy Python：如何在网格上插值Angular ？

2条答案

导入库

让我们使用此示例数据集

插值

打印结果

相关问题

热门标签

最新问答