Shade3D 公式

26 滑らかな掃引( 1/2 ) [ Shade Labo ]


#1

関連記事:
  • 25 Custom 掃引

26 - 1  概要


通常の掃引では、掃引中心線の anchor point 上に掃引断面が配置され、それらをできるだけスムースな線で結んだ自由曲面として作られます。

しかし、中心線が anchor point 間で激しく曲がっている場合、得られる形状には大きな歪みが生じます。

滑らかな掃引とは、中心線の曲がりに応じて断面配置の数を増やし、歪みの少ない滑らかな形状の掃引体を得ようとするものです。

           [ fig 26 - 1 ]


この手順は次の2つのステップに分けられます。
  1. 曲がりに応じて中心線を分割

  2. ポイントが追加された中心線をベースにして掃引体を作成


後半の掃引処理については **25 Custom 掃引** で作成した **custom_sweep( )** を用いますので、ここでは前半のポイント追加をどのように処理するかについて紹介することになります。

また、実際に Shade 上で中心線を分割するのではなく、 読み込まれた中心線を script 内部で分割する処理を custom_sweep( ) に組み込み、関数 c_sweep( ) として作成します。


曲線の曲がり具合に応じてそれを分割するための key word は 曲率変曲点 の二つです。


#2

26 - 2  曲率


線の曲がり具合を表すのが **曲率 κ** で、2次元曲線では左右の曲がり方に応じて正負の符号が付きますが、3次元曲線ではκ >= 0 で表されます。
  • 曲率が大きい      曲線の曲がり具合が急激
  • 曲率が小さい      曲線の曲がり具合が緩やか
  • 曲率 = 0              曲がっていない = 直線

          [ fig 26 - 2 ]


三次元曲線の bezier line P( t ) の parameter t における曲率 κ は次式で与えられます。

          [ 式 26 - 1 ]


#3

26 - 3  全曲率


曲線のトータルでの曲がり量を **全曲率** ( >= 0 )といいます。

bezier 曲線が 平面上にある場合、その曲がりが C curve を描くなら inhandle と outhandle のなす角 θ が 全曲率になり、S curve を描くなら θ1 + θ2 となります.。

          [ fig 26 - 3 ]


全曲率は曲率を積分することで求められます。

          [ 式 26 - 2 ]


また全曲線は次のようにも表現できます。

ある一点を始点として曲線の接線単位ベクトルをプロットしていくとき、プロットした単位ベクトルの先端が描く軌跡の長さ

          [ fig 26 - 4 ]


三次元曲線の場合は、 **単位半径の球面上の軌跡の長さ** になり、全曲率が可視化されます。

          [ fig 26 - 5 ]


#4

26 - 4  変曲点


**20 - 1 S curve** で、二次元の bezier line の曲がり方について観察しました。

      [ fig 20 - 1 ]


E ~ H では曲線は S curve を描きますが、この S curve の曲がる方向が変化する点が **変曲点** で、その位置の求め方を考えます。

変曲点は平面曲線に対して定義され 三次元曲線には存在しません が、ここでは後述する 屈曲点 のからみから、三次元も含めた一般の bezier line の変曲点を求めるというストーリーで話を進めていきます。




変曲点では 曲率 = 0 となりますから、その点での bezier parameter t を求めるには、[ 式 26 - 1 ] より、


          [ 式 26 - 3 ]


          [ 式 26 - 4 ]


[ 式 26 - 4 ] をそのまま展開すると極めて煩雑な式になるので、それを避けるため **V1**, **V2** を求める matrix 演算の部分を次のように変形します。

     

     [ 式 26 - 5 ]


[ Df ] は [ **P** ] を階差表現に変換する matrix で、[ **P** ] と [ **Q** ] との間には次のような関係があります。
          [ 式 26 - 6 ]
一方、

          [ 式 26 - 7 ]


よって [ 式 26 - 4 ] の V1, V2 はより簡潔に表記され、


          [ 式 26 - 8 ]


V1 x V2 を展開すれば、


          [ 式 26 - 9 ]


| **V1** x **V2** | = 0 より、**V1** x **V2** の各方向成分が全て 0 ですから、
          [ 式 26 - 10 ]

[ 式 26 - 10 ] の t に解が存在するには、**A**, **B**, **C** の各ベクトルが平行であることが必要条件になります。

一方、ベクトル A, B, C は ベクトル q1, q2, q3 間の相互積で与えられますから、q1, q2, q3 は同一平面上になくてはなりません。

この q1, q2, q3 は [ 式 26 - 6 ] が示すように、p0, p1, p2, p3 の階差で与えられますから、結局、p0, p1, p2, p3 が同一平面上にあることが要求されます。


つまり、bezier line が平面上に存在している ことが変曲点が存在しうる必要条件となり、平面曲線でなければ変曲点は存在しないという当然の帰結に行き着きます。

平面曲線では、[ 式 26 - 10 ] で与えられる二次方程式の解の内、0 < t < 1 なる範囲にある解 が変曲点の位置となります。




bezier line がある平面内に存在しているなら、[ 式 26 - 10 ] のベクトル A, B, C の X, Y, Z いずれかの成分に対して二次方程式を解けば、他の成分も必ず同じ解を満足します。

ただし、ベクトル A, B, C の X, Y, Z 成分から適当に選んでしまうと、bezier line が X-Y, Y-Z, Z-X 平面のいずれかに平行な場合に解が求められないケースが生じてしまいます。

例えば X-Z 平面に平行なら、A, B, C は共に q1, q2, q3 の内の2つのベクトルの外積ですから、それらの X, Z 方向成分の積で与えられる Y 成分を用いなければ方程式は解けません。


これを一般化すれば、次の要領で解を求めることになります。
  • [ Q ] の bounding box size を求め、その X, Y, Z size の内、最小の size を与える成分 W ( = X or Y or Z ) を求める

  • [ 式 26 - 10 ] の二次方程式を ベクトル A, B, C それぞれの W 成分を用いて解く


#5

26 - 5  屈曲点


大きくS字形に曲がった平面状の bezier line に対して滑らかな掃引を作る際には、**曲がり方が変化するポイント** として変曲点を検出することが求められますが、handle が僅かに平面上からずれただけで変曲点がないことになってしまうのは、滑らかさを維持する上で好ましいことではありません。

          [ fig 26 - 6 ]


そこで三次元曲線でも **曲がり方が大きく変化するポイント** を求めることにし、それを **屈曲点** と呼ぶことにします。

これはオーソライズな用語ではなく単に便宜的に呼び表すもので、推奨できるような用語ではありません



屈曲点を次のように定義します。

bezier line が平面曲線であれば、ここで定義される屈曲点は変曲点に一致します。


[ 定義 26 - 1 ]

4つの bezier control points [ P0 P1 P2 P3 ] に対して、方向ベクトル Vn を次のように定義する

          W1 = Va x Vb = ( P1 - P0 ) x ( p2 - p3 )

          W2 = Vc x Vd = ( P2 - P1 ) x ( P3 - P0 )

          W3 = Ve x Vf = ( P2 - P0 ) x ( P1 - P3 )

W1, W2, W3 の内、絶対値の最大なるものを Vn とする

bezier line を面法線 Vn なる面上に投影した平面曲線に変曲点があれば、その bezier parameter で与えられるポイントを 屈曲点 とする


          [ fig 26 - 7 ]


#6

26 - 6  屈曲点を求める script


以上より、屈曲点の bezier parameter を返す関数 **bezier_inflection_point( )** を作成し、[ script 26 - 3 ] ではこれを用いて、Shade 上で選択した線形状の **最初のブロック** について屈曲点を求め、そこにマーカーを出力します。

[ script 26 - 3 ]

import labo
from vec3 import *
from quaternion import *
	
		

#  bezier_inflection_point( bz as vec3 list) as float list
#
#		bz : bezier control point 座標を格納した vec3 list
#		bezier line bz の屈曲点位置の bezier parameter を返す。
#		屈曲点がない場合は、空 list が返される

def bezier_inflection_point(bz) :
	w1 = (bz[1] - bz[0])*(bz[2] - bz[3])
	w2 = (bz[2] - bz[1])*(bz[3] - bz[0])
	w3 = (bz[2] - bz[0])*(bz[1] - bz[3])
	r1 = w1.norm()
	r2 = w2.norm()
	r3 = w3.norm()
	L = [[r1, w1], [r2, w2], [r3, w3]]
	L.sort()
	vn = L[2][1]
	if vn.norm() == 0 :				#  bz は直線であり、変曲点は存在しない
		return []

	#  bz2 : 原点を通り面法線 vn なる面に bz を投影	
	bz2 = []
	for i in range(4) :
		d = vn.dot(bz[i])
		bz2.append(bz[i] - d*vn)
	
	#  Q : bz2 の階差表現
	bz2 = matrix(False, bz2)
	Md = labo.bezier_d()
	Q = bz2*Md
	
	#  こちらの方法でも可
#	bz2 = matrix(True, bz2)
#	Md = labo.bezier_d()
#	Md.transpose()
#	Q = Md*bz2			

	A = Q[2]*Q[3]
	B = Q[1]*Q[3]
	C = Q[1]*Q[2]
	
	bb = labo.vec3_bounding_box_size(Q)
	bb.pop()
	bb2 = sorted(bb)
	for i in range(3) :
		if bb2[0] == bb[i] :
			a = A[i]
			b = B[i]
			c = C[i]
			t = []
			if a != 0 :
				d = b**2 - 4*a*c
				if d > 0 :
					tt = (-b - d**0.5)/(2*a)
					if tt > 0 and tt < 1 :
						t.append(tt)
					tt = (-b + d**0.5)/(2*a)
					if tt >= 0 and tt <= 1 :
						t.append(tt)
				elif d == 0 :
					tt = -b/(2*a)
					if tt > 0 and tt < 1 :
						t.append(tt)
			elif b != 0 :
				tt = -c/b
				if tt > 0 and tt < 1 :
					t.append(tt)
			
			if len(t) == 2 :
				t.sort()	
			return t
				
					

		

scene = xshade.scene()

bz = labo.get_bezier(xshade, 0)
t = bezier_inflection_point(bz)
n = len(t)
if n > 0 :
	obL = []
	obL.append(scene.active_shape())
	
	bb = labo.vec3_bounding_box_size(bz)
	size = bb[3]/32
	for i in range(n) :
		p = labo.bezier_line_position(bz, t[i])
		scene.create_sphere(str(t[i]), p, size)
		obL.append(scene.active_shape())
		
	scene.active_shapes = obL
else :
	print '屈曲点はありません'