7

我有一个形状数组(201,201),我想通过数据绘制一些横截面,但我无法访问相关点。例如说我想绘制由产生的图中的线给出的横截面,

from pylab import *
Z = randn(201,201)
x = linspace(-1,1,201)
X,Y = meshgrid(x,x)
pcolormesh(X,Y,Z)
plot(x,x*.5)

我想在不同的方向绘制这些图,但如果有帮助的话,它们总是会穿过原点......

4

1 回答 1

9

基本上,您想沿一条线(或任意路径)插入一个 2D 网格。

首先,您应该决定是要插入网格还是只进行最近邻采样。如果你想做后者,你可以使用索引。

如果您想插值,请查看scipy.ndimage.map_coordinates. 一开始有点难以理解,但它非常适合这个。(这比使用假设数据点随机分布的插值例程要高效得多。)

我将举两个例子。这些改编自我对另一个问题的回答。 但是,在这些示例中,所有内容都绘制在“像素”(即行、列)坐标中。

在您的情况下,您在与“像素”坐标不同的坐标系中工作,因此您需要将“世界”(即 x,y)坐标转换为“像素”坐标以进行插值。

首先,这是一个使用三次插值的示例map_coordinates

import numpy as np
import scipy.ndimage
import matplotlib.pyplot as plt

# Generate some data...
x, y = np.mgrid[-5:5:0.1, -5:5:0.1]
z = np.sqrt(x**2 + y**2) + np.sin(x**2 + y**2)

# Coordinates of the line we'd like to sample along
line = [(-3, -1), (4, 3)]

# Convert the line to pixel/index coordinates
x_world, y_world = np.array(zip(*line))
col = z.shape[1] * (x_world - x.min()) / x.ptp()
row = z.shape[0] * (y_world - y.min()) / y.ptp()

# Interpolate the line at "num" points...
num = 1000
row, col = [np.linspace(item[0], item[1], num) for item in [row, col]]

# Extract the values along the line, using cubic interpolation
zi = scipy.ndimage.map_coordinates(z, np.vstack((row, col)))

# Plot...
fig, axes = plt.subplots(nrows=2)
axes[0].pcolormesh(x, y, z)
axes[0].plot(x_world, y_world, 'ro-')
axes[0].axis('image')

axes[1].plot(zi)

plt.show()

在此处输入图像描述

或者,我们可以使用最近邻插值。一种方法是在上面的示例中传递order=0map_coordinates。相反,我将使用索引来展示另一种方法。如果我们只是换行

# Extract the values along the line, using cubic interpolation
zi = scipy.ndimage.map_coordinates(z, np.vstack((row, col)))

至:

# Extract the values along the line, using nearest-neighbor interpolation
zi = z[row.astype(int), col.astype(int)]

我们会得到:

在此处输入图像描述

于 2013-09-20T16:26:26.170 回答