Как лучше всего создать заполненный объем Аттрактора Лоренца в python для объемного рендеринга

Question

Я пытаюсь создать трехмерный массив, чтобы затем выполнить объемный рендеринг (в другом программном обеспечении или пакетах объемного рендеринга) странного аттрактора, такого как Аттрактор Лоренца. Достаточно легко построить аттрактор из точек данных и предоставить значение для назначения цвета и визуализации, например, в matplotlib.

Однако я бы хотел массив заполненных объемов. Я пробовал такие методы интерполяции, как griddata, но они не давали желаемого результата. Я представляю себе что-то вроде:

Это со страницы википедии.

Вот то, что я пробовал, но если вы откроете результат в простой программе просмотра, он будет выглядеть не очень хорошо. Я думаю, вместо этого, возможно, выполняю интерполяцию только между точками, которые составляют массив x, y, z ... Я немного потерялся после того, как поиграл с этим в течение нескольких часов. Я думаю, что мне нужно взять точки и выполнить какую-то интерполяцию или заполнение массива, здесь я вызываю interp_im. Затем это можно просмотреть в объемном рендеринге. Любая помощь приветствуется!

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import odeint
from scipy.interpolate import griddata
from scipy.interpolate import LinearNDInterpolator
from skimage.external import tifffile

rho = 28.0
sigma = 10.0
beta = 8.0 / 3.0

def f(state, t):
    x, y, z = state  # Unpack the state vector
    return sigma * (y - x), x * (rho - z) - y, x * y - beta * z  # Derivatives

state0 = [1.0, 1.0, 1.0]
t = np.arange(0.0, 40.0, 0.01) #t = np.arange(0.0, 40.0, 0.01)

states = odeint(f, state0, t)

# shift x,y,z positions to int for regular image volume

x = states[:, 0]
y = states[:, 1]
z = states[:, 2]
x_min = x.min()
y_min = y.min()
z_min = z.min()
states_int = states + [abs(x_min),abs(y_min),abs(z_min)] + 1
states_int = states_int * 10
states_int = states_int.astype(int)

# values will be in order of tracing for color
values = []
for i,j in enumerate(states_int):
    values.append(i*10)

values = np.asarray(values)

fig = plt.figure()
ax = fig.gca(projection='3d')
sc = ax.scatter(states_int[:, 0], states_int[:, 1], states_int[:, 2],c=values)
plt.colorbar(sc)
plt.draw()
plt.show()

#print(x.shape, y.shape, z.shape, values.shape)

#Interpolate for volume rendering
x_ = np.linspace(0,999,500)
y_ = np.linspace(0,999,500)
z_ = np.linspace(0,999,500)
xx,yy,zz = np.meshgrid(x_,y_,z_, sparse = True)
#
# X = states_int.tolist()
#
interp_im = griddata(states_int, values, (xx,yy,zz), method='linear')
interp_im = interp_im.astype(np.uint16)

np.save('interp_im.npy', interp_im)

tifffile.imsave('LorenzAttractor.tif', interp_im)

Answer 1

Ваши данные в томе, они просто пикселированы. Если вы размываете объем, например, с помощью гаусса, вы получаете что-то гораздо более удобное. Например:

from scipy import ndimage

vol = np.zeros((512, 512, 512), dtype=states_int.dtype)
# add data to vol
vol[tuple(np.split(states_int, vol.ndim, axis=1))] = values[:, np.newaxis]
# apply gaussian filter, sigma=5 in this case
vol = ndimage.gaussian_filter(vol, 5)

Затем я бы использовал что-то вроде napari для просмотра данных в 3D:

import napari

with napari.gui_qt():
    napari.view_image(v)

Чтобы сделать громкость более плавной, вы можете уменьшить размер шага интегрирования.