三维数据在Python插值/子采样不VTK插值、数据、Python、VTK

2023-09-08 00:44:25 作者:夜巷

我想要做的是相当简单的,但我还没有找到一个简单的方法迄今:

What I want to do is rather simple but I havent found a straightforward approach thus far:

我有一个三维直线网格与浮点值(因此3坐标轴-1D numpy的arrays-网格单元和3D numpy的阵列,为每个单元中心的值相应形状的中心),我想插值(或者你可以称之为子样本)这整个阵列子采样阵列(5如大小的因素)与线性插值。 所有我见过这么远涉及2D,然后一维插值或VTK的技巧而宁愿不使用(便携),该方法

I have a 3D rectilinear grid with float values (therefore 3 coordinate axes -1D numpy arrays- for the centers of the grid cells and a 3D numpy array with the corresponding shape with a value for each cell center) and I want to interpolate (or you may call it subsample) this entire array to a subsampled array (e.g. size factor of 5) with linear interpolation. All the approaches I've seen this far involve 2D and then 1D interpolation or VTK tricks which Id rather not use (portability).

可能有人建议将是服用5×5×细胞同时的三维阵列中,平均和返回的阵列5倍小的在每个方向上的等效的方法?

Could someone suggest an approach that would be the equivalent of taking 5x5x5 cells at the same time in the 3D array, averaging and returning an array 5times smaller in each direction?

感谢您事先的任何建议

编辑: 这里的数据是什么样子,D是一款3D阵列重新presenting细胞的三维网格。每个单元都有一个标量浮点值(在我的情况下,pressure)和X,Y和Z是包含每一个细胞的细胞的空间坐标(见图形3维数组,以及如何' X'阵列的样子)

Here's what the data looks like, 'd' is a 3D array representing a 3D grid of cells. Each cell has a scalar float value (pressure in my case) and 'x','y' and 'z' are three 1D arrays containing the spatial coordinates of the cells of every cell (see the shapes and how the 'x' array looks like)

In [42]: x.shape
Out[42]: (181L,)

In [43]: y.shape
Out[43]: (181L,)

In [44]: z.shape
Out[44]: (421L,)

In [45]: d.shape
Out[45]: (181L, 181L, 421L)

In [46]: x
Out[46]: 
array([-0.410607  , -0.3927568 , -0.37780656, -0.36527296, -0.35475321,
       -0.34591168, -0.33846866, -0.33219107, -0.32688467, -0.3223876 ,
        ...
        0.34591168,  0.35475321,  0.36527296,  0.37780656,  0.3927568 ,
        0.410607  ])

我想要做的就是创建一个三维阵列,可以说的90x90x210形状(约缩小了2倍)由第一子采样来自轴坐标上的阵列,以上尺寸再插在3D数据数组。林不知道是否'插'是正确的术语,但。向下取样?平均? 下面是数据的二维切片:

What I want to do is create a 3D array with lets say a shape of 90x90x210 (roughly downsize by a factor of 2) by first subsampling the coordinates from the axes on arrays with the above dimensions and then 'interpolating' the 3D data to that array. Im not sure whether 'interpolating' is the right term though. Downsampling? Averaging? Here's an 2D slice of the data:

推荐答案

下面是3D插值对不规则网格使用的 scipy.interpolate.griddata 。

Here is an example of 3D interpolation on an irregular grid using scipy.interpolate.griddata.

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


def func(x, y, z):
    return x ** 2 + y ** 2 + z ** 2

# Nx, Ny, Nz = 181, 181, 421
Nx, Ny, Nz = 18, 18, 42

subsample = 2
Mx, My, Mz = Nx // subsample, Ny // subsample, Nz // subsample

# Define irregularly spaced arrays
x = np.random.random(Nx)
y = np.random.random(Ny)
z = np.random.random(Nz)

# Compute the matrix D of shape (Nx, Ny, Nz).
# D could be experimental data, but here I'll define it using func
# D[i,j,k] is associated with location (x[i], y[j], z[k])
X_irregular, Y_irregular, Z_irregular = (
    x[:, None, None], y[None, :, None], z[None, None, :])
D = func(X_irregular, Y_irregular, Z_irregular)

# Create a uniformly spaced grid
xi = np.linspace(x.min(), x.max(), Mx)
yi = np.linspace(y.min(), y.max(), My)
zi = np.linspace(y.min(), y.max(), Mz)
X_uniform, Y_uniform, Z_uniform = (
    xi[:, None, None], yi[None, :, None], zi[None, None, :])

# To use griddata, I need 1D-arrays for x, y, z of length 
# len(D.ravel()) = Nx*Ny*Nz.
# To do this, I broadcast up my *_irregular arrays to each be 
# of shape (Nx, Ny, Nz)
# and then use ravel() to make them 1D-arrays
X_irregular, Y_irregular, Z_irregular = np.broadcast_arrays(
    X_irregular, Y_irregular, Z_irregular)
D_interpolated = interpolate.griddata(
    (X_irregular.ravel(), Y_irregular.ravel(), Z_irregular.ravel()),
    D.ravel(),
    (X_uniform, Y_uniform, Z_uniform),
    method='linear')

print(D_interpolated.shape)
# (90, 90, 210)

# Make plots
fig, ax = plt.subplots(2)

# Choose a z value in the uniform z-grid
# Let's take the middle value
zindex = Mz // 2
z_crosssection = zi[zindex]

# Plot a cross-section of the raw irregularly spaced data
X_irr, Y_irr = np.meshgrid(sorted(x), sorted(y))
# find the value in the irregular z-grid closest to z_crosssection
z_near_cross = z[(np.abs(z - z_crosssection)).argmin()]
ax[0].contourf(X_irr, Y_irr, func(X_irr, Y_irr, z_near_cross))
ax[0].scatter(X_irr, Y_irr, c='white', s=20)   
ax[0].set_title('Cross-section of irregular data')
ax[0].set_xlim(x.min(), x.max())
ax[0].set_ylim(y.min(), y.max())

# Plot a cross-section of the Interpolated uniformly spaced data
X_unif, Y_unif = np.meshgrid(xi, yi)
ax[1].contourf(X_unif, Y_unif, D_interpolated[:, :, zindex])
ax[1].scatter(X_unif, Y_unif, c='white', s=20)
ax[1].set_title('Cross-section of downsampled and interpolated data')
ax[1].set_xlim(x.min(), x.max())
ax[1].set_ylim(y.min(), y.max())

plt.show()