Mlab3Dと2Dの例¶
特定のシーンについて，3D全体座標から2D表示座標(ピクセル座標)への投影を計算するスクリプト．
イメージプレーンのオブジェクトの2次元ピクセル位置は，一連の線形変換によって3次元ワールド座標に関連付けられます．特定の変換は射影変換として知られるグループに属する．この集合は純粋な射影性，affine 変換，遠近変換，及びeuclid変換を含みます．mlab(他の多くの可視化ソフトウェアは)の場合は，パースペクティブとeuclidの場合のみを扱います．射影空間の概要は http://en.wikipedia.org/wiki/Projective_space にあり，射影幾何学の完全な取り扱いはRichard Hartleyの本 "Multiple View Geometry in Computer Vision" にあります．
この例で知っておくべき重要なことは，3つの空間内の点は，4x4行列の一連の遠近法とeuclid変換の乗算によって2つの空間内の点に関連付けられるということです．4x4行列は，点を表すために長さ4のベクトルを使用することを述語とします．この表現は同次座標として知られており，最初はそのように見えますが，関係するすべての数学を真に単純化します．要するに，同次座標はあなたの友人であり，あなたはそれらについてここで読むべきです http://en.wikipedia.org/wiki/Homogeneous_coordinates
通常のピンホールカメラモデル(理想的な現実世界モデル)において,3次元ワールド点は,透視変換とユークリッド変換の組合せである '必須の' 行列と呼ばれる行列により2次元画像点に関係づけられます．パースペクティブ変換はカメラの組み込み関数(焦点距離，撮像素子オフセットなど．)によって定義され，ユークリッド変換はカメラの位置と方向によって定義されます．コンピュータグラフィックスでは，問題はそれほど単純ではありません．これは，コンピュータグラフィックスには，クリッピング平面，オフセット投影中心，任意の歪みなど，実世界では不可能なことを実行できるという利点があるためです．したがって，若干異なるモデルが使用されます．
これに続くのは，OpenGLのカメラ/ビューモデル，つまりVTKです．他のパッケージがこのモデルに従うことを保証することはできません．
3Dワールド座標を適用して2Dピクセル座標にマッピングする変換には，4種類あります．モデル変換，ビュー変換，パース変換，ビューポートまたは表示変換です．
OpenGLでは，最初の2つの変換が連結され，modelview変換(VTKでは単にビュー変換と呼ばれます)が生成されます．modelview変換は，任意のスケーリングと歪みをモデルに適用し(指定されている場合)，方向が負のZ軸を見下ろすのと同じになるようにモデルを変換します．たとえば，カメラを移動して負のZ軸を見下ろすようにしてから，ワールド内のすべてのものを移動して，カメラを移動する前と同じように見えるようにしたとします．作成された座標は，OpenGLでは "eye" 座標と呼ばれます(VTKに名前があるかは知りません)．
パース変換は，カメラのパースを目の座標に適用します．この変換によって，フォアグラウンドのオブジェクトがバックグラウンドの同等のオブジェクトより大きく表示されます．ピンホールカメラモデルでは,この変換はカメラの焦点距離と3空間におけるその位置によって一意的に決定される．Vtk/OpenGLでは，フラスタムによって決まります．フラスタムとは，頂点が切り取られたピラミッドのことです．ピラミッドの最上部(一点)はカメラ位置であり，ピラミッドの底面は，ファークリッピング距離と呼ばれる距離で主カメラ光線に垂直に定義される平面(ファークリッピングプレーン)であり，フラスタムの最上部(切断された部分)は，ファークリッピング平面と同様の定義でニアクリッピング平面である．フラスタムの側面は，カメラのアスペクト比(幅/高さ)と視野によって決まります．フラスタム内にないポイントは画面にマッピングされません(表示領域の外側にある場合と同じです．)．パース変換には，等質座標で表される (-1,1)(-1,1)(-1,1) の範囲で定義される立方体に収まるように，フラスタム内のすべてをスケールする効果があります．最後のフレーズが重要で，最初の3つの座標は，通常，最後の座標(これが紛らわしい場合は，同次座標のwikipediaを参照してください．)で除算するまでは単一の範囲内にありません．得られた座標は，(適切に)正規化ビュー座標と呼ばれる．
最後の変換(ビューポート変換)は，正規化されたビュー座標から座標を表示します．この時点で，「座標を直接表示するだけでなく，正規化されたビュー座標が必要になるのはなぜですか? 」と自問しているかもしれませんが，答えは，特定のウィンドウに複数のビューを埋め込みたい場合があるため，各ビューをウィンドウ内の適切な位置にサイズ変更する方法が異なるということです．正規化されたビュー座標は，いわば良い共通基盤を提供します．いずれにしても，ビューポートの変換は，単純に，規格化されたビュー座標のXおよびY座標を適切なピクセル座標に尺度変更および変換する．この例ではZ値を使用していません．他にもいろいろな用途に使われています．
これで全部です，かなり簡単ですよね?そうですね．概要は次のとおりです．
- 指定された3Dワールド座標系のセット
modelview変換(VTKのビュー変換)を適用してeye coordinatesを取得します．
透視変換を適用して，正規化されたビュー座標を取得します．
viewport変換を適用して表示座標を取得
VTKには，最初の2つの演算を組み合わせた4x4行列を取得するための優れたメソッドが用意されています．私の知る限り，VTKではビューポート変換を表す4x4行列を取得するメソッドをエクスポートしていないため，ここでは(心配はありませんが，ご覧のように難しくはありません．)を作成します．
準備が終わりましたので，始めましょう．
Pythonソースコード: mlab_3D_to_2D.py
# Author: S. Chris Colbert <sccolbert@gmail.com>
# Copyright (c) 2009, S. Chris Colbert
# License: BSD Style
# this import is here because we need to ensure that matplotlib uses the
# wx backend and having regular code outside the main block is PyTaboo.
# It needs to be imported first, so that matplotlib can impose the
# version of Wx it requires.
import matplotlib
matplotlib.use('WXAgg')
import pylab as pl
import numpy as np
from mayavi import mlab
from mayavi.core.ui.mayavi_scene import MayaviScene
def get_world_to_view_matrix(mlab_scene):
"""returns the 4x4 matrix that is a concatenation of the modelview transform and
perspective transform. Takes as input an mlab scene object."""
if not isinstance(mlab_scene, MayaviScene):
raise TypeError('argument must be an instance of MayaviScene')
# The VTK method needs the aspect ratio and near and far clipping planes
# in order to return the proper transform. So we query the current scene
# object to get the parameters we need.
scene_size = tuple(mlab_scene.get_size())
clip_range = mlab_scene.camera.clipping_range
aspect_ratio = float(scene_size[0])/float(scene_size[1])
# this actually just gets a vtk matrix object, we can't really do anything with it yet
vtk_comb_trans_mat = mlab_scene.camera.get_composite_projection_transform_matrix(
aspect_ratio, clip_range[0], clip_range[1])
# get the vtk mat as a numpy array
np_comb_trans_mat = vtk_comb_trans_mat.to_array()
return np_comb_trans_mat
def get_view_to_display_matrix(mlab_scene):
""" this function returns a 4x4 matrix that will convert normalized
view coordinates to display coordinates. It's assumed that the view should
take up the entire window and that the origin of the window is in the
upper left corner"""
if not (isinstance(mlab_scene, MayaviScene)):
raise TypeError('argument must be an instance of MayaviScene')
# this gets the client size of the window
x, y = tuple(mlab_scene.get_size())
# normalized view coordinates have the origin in the middle of the space
# so we need to scale by width and height of the display window and shift
# by half width and half height. The matrix accomplishes that.
view_to_disp_mat = np.array([[x/2.0, 0., 0., x/2.0],
[ 0., -y/2.0, 0., y/2.0],
[ 0., 0., 1., 0.],
[ 0., 0., 0., 1.]])
return view_to_disp_mat
def apply_transform_to_points(points, trans_mat):
"""a function that applies a 4x4 transformation matrix to an of
homogeneous points. The array of points should have shape Nx4"""
if not trans_mat.shape == (4, 4):
raise ValueError('transform matrix must be 4x4')
if not points.shape[1] == 4:
raise ValueError('point array must have shape Nx4')
return np.dot(trans_mat, points.T).T
if __name__ == '__main__':
f = mlab.figure()
N = 4
# create a few points in 3-space
X = np.random.random_integers(-3, 3, N)
Y = np.random.random_integers(-3, 3, N)
Z = np.random.random_integers(-3, 3, N)
# plot the points with mlab
pts = mlab.points3d(X, Y, Z)
# now were going to create a single N x 4 array of our points
# adding a fourth column of ones expresses the world points in
# homogenous coordinates
W = np.ones(X.shape)
hmgns_world_coords = np.column_stack((X, Y, Z, W))
# applying the first transform will give us 'unnormalized' view
# coordinates we also have to get the transform matrix for the
# current scene view
comb_trans_mat = get_world_to_view_matrix(f.scene)
view_coords = \
apply_transform_to_points(hmgns_world_coords, comb_trans_mat)
# to get normalized view coordinates, we divide through by the fourth
# element
norm_view_coords = view_coords / (view_coords[:, 3].reshape(-1, 1))
# the last step is to transform from normalized view coordinates to
# display coordinates.
view_to_disp_mat = get_view_to_display_matrix(f.scene)
disp_coords = apply_transform_to_points(norm_view_coords, view_to_disp_mat)
# at this point disp_coords is an Nx4 array of homogenous coordinates
# where X and Y are the pixel coordinates of the X and Y 3D world
# coordinates, so lets take a screenshot of mlab view and open it
# with matplotlib so we can check the accuracy
img = mlab.screenshot()
pl.imshow(img)
for i in range(N):
print('Point %d: (x, y) ' % i, disp_coords[:, 0:2][i])
pl.plot([disp_coords[:, 0][i]], [disp_coords[:, 1][i]], 'ro')
pl.show()
# you should check that the printed coordinates correspond to the
# proper points on the screen
mlab.show()
#EOF
