ついにWATLABブログから書籍「いきなりプログラミングPython」が発売しました！
PyVistaを使えばSTLファイルを簡単に表示したり情報を取得したりする事ができるとわかりました。ここでは取得したデータを加工してみます。まずは座標変換として拡大・縮小、平行移動、回転移動をPythonで実装してみます。
こんにちは。wat(@watlablog)です。ここではPythonとPyVistaを駆使してSTLファイルの座標変換をやってみます！
3D座標変換の概要と種類
座標変換の種類
拡大・縮小
ここでは今回実装する座標変換のイメージ図を示します。まずは拡大縮小から。
3D形状の拡大・縮小とは、下図のような変形の事です。3軸方向に倍率をかけてあげれば達成できるため、非常に簡単な変形です。
3軸の倍率、すなわちスケールファクターを\(S\)とすると、\(xyz\)座標値に対して各スケールファクターとのアダマール積をとれば良いです（行列積の方法の方が一般的かもですが）。
※アダマール積とは、行列の順番に則って積を演算するのではなく、行列の対応している位置同士で積をとる方法。Numpyでそのまま「*」を使ってやるのがアダマール積、np.dot()で演算するのが行列積です。
平行移動
平行移動はもっと簡単かも知れません。以下の図のように一様にずらすだけです。
平行移動も移動量\(T\)を\(xyz\)座標値毎にそれぞれ定義し、単純に足すだけで十分です（こちらも行列積の方法もありますが、今回はこの方法で）。
回転移動
回転移動だけはちょっとやっかいですが、イメージ図は下図になります。基本は原点周りで回転させる事を前提とします。
回転移動の場合は以下のように\(R_{x},R_{y},R_{z}\)と回転行列をそれぞれ定義した方がやりやすいです。
定義した回転行列を各点の\(xyz\)座標に対し、行列式で以下のように計算すると、回転させる事が出来ます。但し、この方法は原点周りのみです。任意点周りの場合は予め原点に平行移動させ、回転させてからまた逆の平行移動をさせる方法をとるか、予め回転行列に平行移動分を加味しておくかです。
座標変換には他にもアフィン変換や射影変換といった様々な手法がありますが、今回はこの3つをPythonで実装していきます。
PythonでSTLの座標変換をするコード
PyVistaを使う
STLファイルの取り扱いにはPyVistaを使います。PyVistaについては、以下の記事にまとめましたので参照下さい。
「PyVistaをインストールしてPythonでSTLを扱う備忘録」
拡大・縮小を行うコード
まずは拡大・縮小を行うコードです。
以下のコードはdef magnification()で拡大・縮小の操作を関数化しています。引数はSTLの点情報（mesh.pointsで取得するもの）と、各方向の倍率です。
今回はモデルを表示させずにpv.save_meshio()を使ってSTLを保存する所まで記載しています。今回はy方向のみを2倍にしてみます。
上記説明で行列積を用いなかったのはNumpyのスライスを使って各成分毎に操作したかったからです（もっと効率良い書き方ありそうですが）。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 |
import pyvista as pv import numpy as np # STLの拡大・縮小をするコード def magnification(points, Mx, My, Mz): M = np.array([Mx, My, Mz]) for i in range(3): points[:, i] = points[:, i] * M[i] return points filename = 'sample.stl' # STLファイル名 mesh = pv.read(filename) # STLファイルを読み込む points = mesh.points # 点の座標情報 mesh.points = magnification(points, 1, 2, 1) # 拡大・縮小を行う関数を実行 pv.save_meshio('new.stl', mesh) # 座標変換後のSTLファイルを保存 |
以下が変換前（before）と変換後（after）です。見事STLモデルのy方向のみが大きくなりました。
平行移動を行うコード
以下が平行移動を行うコードです。
こちらも上記で行列積を用いなかった理由はNumpyスライスを使いたかったからという理由だけです。やり方は拡大・縮小の時と同じですね。
今回はx方向にのみ100だけずらしてみましょう。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 |
# STLの平行移動をするコード def translation(points, Tx, Ty, Tz): T = np.array([Tx, Ty, Tz]) for i in range(3): points[:, i] = points[:, i] + T[i] return points filename = 'sample.stl' # STLファイル名 mesh = pv.read(filename) # STLファイルを読み込む points = mesh.points # 点の座標情報 print('移動前＝', mesh.points) pv.save_meshio('new.stl', mesh) # 座標変換後のSTLファイルを保存 print('移動後＝', mesh.points) mesh.points = translation(points, 100, 0, 0) # 平行移動を行う関数を実行 |
以下が結果です。見事xだけ100ずれました（軸に着目して下さい）。しかし、平行移動なのに対比物がないのでよくわかりませんね。
以下が点座標print文の結果です。移動前と移動後を比べると、こちらであれば確かに意図した平行移動がされているという事がわかると思います。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 |
移動前＝ [[ 0. 1. 2.1903157 ] [-1. 1. 1. ] [-1. -1. 1. ] [ 0. -1. 2.1903157 ] [ 0. -1. 0.19031572] [-1. -1. -1. ] [-1. 1. -1. ] [ 0. 1. 0.19031572] [ 1. 1. -1. ] [ 1. -1. -1. ] [ 1. 1. 1. ] [ 1. -1. 1. ]] 移動後＝ [[100. 1. 2.1903157 ] [ 99. 1. 1. ] [ 99. -1. 1. ] [100. -1. 2.1903157 ] [100. -1. 0.19031572] [ 99. -1. -1. ] [ 99. 1. -1. ] [100. 1. 0.19031572] [101. 1. -1. ] [101. -1. -1. ] [101. 1. 1. ] [101. -1. 1. ]] |
回転移動を行うコード
任意点周りの回転を行うコードを以下に示します。原点周りでは無いので、先ほどの平行移動を行う関数と組み合わせています。
各点毎に回転させるためにfor文を使って行列積を行っている所が苦肉の策なのですが、for文を使わないで書くにはどうしたら良いのだろう？？
（「Pythonでfor書いたら負け」という言葉があるようです…。）
今回はx軸のみを90度回転させています。
入力は度ですが、計算時はラジアンに直しています。
関数の引数にcenterがありますが、こちらが任意の回転中心の座標です。
今はmesh.centerでモデルの重心を入れています。そのため、回転移動後でも重心は変わりません。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 |
import pyvista as pv import numpy as np # STLの平行移動をするコード def translation(points, Tx, Ty, Tz): T = np.array([Tx, Ty, Tz]) for i in range(3): points[:, i] = points[:, i] + T[i] return points # STLの回転をするコード def rot(points, center, theta_x, theta_y, theta_z): theta_x = np.radians(theta_x) theta_y = np.radians(theta_y) theta_z = np.radians(theta_z) rx = [[1,0,0], [0,np.cos(theta_x),-np.sin(theta_x)], [0,np.sin(theta_x),np.cos(theta_x)]] ry = [[np.cos(theta_y),0,np.sin(theta_y)], [0,1,0], [-np.sin(theta_y),0,np.cos(theta_y)]] rz = [[np.cos(theta_z),-np.sin(theta_z),0], [np.sin(theta_z),np.cos(theta_z),0], [0,0,1]] points = translation(points, -center[0], -center[1], -center[2]) for i in range(len(points)): points[i] = np.dot(rz,np.dot(ry, np.dot(rx, points[i]))) points = translation(points, center[0], center[1], center[2]) return points filename = 'sample.stl' # STLファイル名 mesh = pv.read(filename) # STLファイルを読み込む points = mesh.points # 点の座標情報 mesh.points = rot(points, mesh.center, 90, 0, 0) # 回転を行う関数を実行 pv.save_meshio('new.stl', mesh) # 座標変換後のSTLファイルを保存 |
以下が結果です。こちらも意図通りにx軸で90度回転させる事が出来ました。
おまけ：結果表示のコード
先ほどまででbefore, afterの図を載せていましたが、PyVistaはmatplotlibのような感覚でsubplotを作って複数モデルを一度に表示させる事が可能です。
以下にコードを載せておきますので、参考までに。
といっても、公式（https://docs.pyvista.org/examples/02-plot/multi-window.html）のページに記載のものをいじっただけです。是非公式もご覧下さい。
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 |
import pyvista as pv filename1 = 'sample.stl' mesh1 = pv.read(filename1) filename2 = 'new.stl' mesh2 = pv.read(filename2) plotter = pv.Plotter(shape=(1, 2)) plotter.subplot(0, 0) plotter.add_text("before", font_size=20) plotter.add_mesh(mesh1) plotter.show_bounds(all_edges=True) plotter.subplot(0, 1) plotter.add_text("after", font_size=20) plotter.add_mesh(mesh2) plotter.show_bounds(all_edges=True) plotter.show() |
まとめ
本記事ではSTLファイルのような3D形状を有するモデルの座標変換の概要を学び、その方法を記録しました。
また、PyVistaを使う事で簡単にSTLファイルの座標変換をするPythonコードが書けました。
STLファイルの座標変換がプログラムで自動で行えるようになるという事は、単品モデルを自動でアセンブリ化する際に利用できそうです。
本ページのコードはMITライセンスに準拠させていただきます（自己責任の範囲でご自由にお使い下さい）。
座標変換も非常に簡単に出来ました！
Twitterでも関連情報をつぶやいているので、wat(@watlablog)のフォローお待ちしています！