5

(64,64,64) などの形状の 3D データを含む配列があります。このデータセットを使用して、点と法線 (結晶学の hkl 平面に似ています) で指定された平面をどのようにプロットしますか? MayaVi でデータを介して平面を回転させることでできることと同様です。

結果のプロットには、ほとんどの場合、非正方形の平面が含まれます。それらはmatplotlib(ある種の非長方形パッチ)で行うことができますか?

編集:私はこれを自分でほとんど解決しました(以下を参照)が、非長方形のパッチをmatplotlibでどのようにプロットできるかまだ疑問に思っています...?

編集:以下の議論のために、私は質問を言い直しました。

4

5 に答える 5

2

これは面白いです。私が今日回答した同様の質問です。行く方法は次のとおりです。補間。scipy.interpolate から griddata を使用できます。

グリッドデータ

このページには非常に優れた例があり、関数のシグネチャは実際のデータに非常に近いものです。

データを補間したい平面上のポイントをどうにかして定義する必要があります。私はこれを見てみましょう、数年前の私の線形代数のレッスン

于 2013-01-03T12:58:13.597 に答える
2

この問題の最後から 2 番目の解決策があります。Plot a plane based on a normal vector and a point in Matlab または matplotlibへの 2 番目の回答を使用して部分的に解決しました。

# coding: utf-8
import numpy as np
from matplotlib.pyplot import imshow,show

A=np.empty((64,64,64)) #This is the data array
def f(x,y):
    return np.sin(x/(2*np.pi))+np.cos(y/(2*np.pi))
xx,yy= np.meshgrid(range(64), range(64))
for x in range(64):
    A[:,:,x]=f(xx,yy)*np.cos(x/np.pi)

N=np.zeros((64,64)) 
"""This is the plane we cut from A. 
It should be larger than 64, due to diagonal planes being larger. 
Will be fixed."""

normal=np.array([-1,-1,1]) #Define cut plane here. Normal vector components restricted to integers
point=np.array([0,0,0])
d = -np.sum(point*normal)

def plane(x,y): # Get plane's z values
    return (-normal[0]*x-normal[1]*y-d)/normal[2]

def getZZ(x,y): #Get z for all values x,y. If z>64 it's out of range
    for i in x:
        for j in y:
            if plane(i,j)<64:
                N[i,j]=A[i,j,plane(i,j)]

getZZ(range(64),range(64))
imshow(N, interpolation="Nearest")
show()

プロットは z 値を持つポイントに制限されず、64 * 64 より大きい平面は考慮されず、平面は (0,0,0) で定義する必要があるため、これは究極の解決策ではありません。

于 2013-01-03T15:33:19.423 に答える
1

要件を減らすために、簡単な例を用意しました

import numpy as np
import pylab as plt

data = np.arange((64**3))
data.resize((64,64,64))

def get_slice(volume, orientation, index):
    orientation2slicefunc = {
        "x" : lambda ar:ar[index,:,:], 
        "y" : lambda ar:ar[:,index,:],  
        "z" : lambda ar:ar[:,:,index]
    }
    return orientation2slicefunc[orientation](volume)

plt.subplot(221)
plt.imshow(get_slice(data, "x", 10), vmin=0, vmax=64**3)

plt.subplot(222)
plt.imshow(get_slice(data, "x", 39), vmin=0, vmax=64**3)

plt.subplot(223)
plt.imshow(get_slice(data, "y", 15), vmin=0, vmax=64**3)
plt.subplot(224)
plt.imshow(get_slice(data, "z", 25), vmin=0, vmax=64**3)  

plt.show()  

これにより、次のプロットが作成されます。

4つのスライス例

主なトリックは、辞書マッピングのオリエンテーションをラムダメソッドにマッピングすることです。これにより、煩わしいif-then-else-blockを作成する必要がなくなります。もちろん、向きに別の名前、たとえば番号を付けることもできます。

多分これはあなたを助けます。

Thorsten

PS:「IndexOutOfRange」については気にしませんでした。このコンテキストでは完全に理解できるので、この例外をポップアウトさせてもかまいません。

于 2013-01-03T16:48:57.960 に答える