Shade online フォーラム
ログイン
ユーザ名:

パスワード:

IDとパスワードを記憶

パスワード紛失
スレッド表示 | 古いものから 投稿するには登録が必要です 前のスレッド | 次のスレッド | 下へ
投稿者 スレッド
投稿数: 189
投稿日時: 2011-11-24 21:12
Re: 線形状(ベジェ曲線)の分割について教えてください
   
< Bezier 曲面の区分曲面 >


Bezier 曲面の U 方向分割については、Bezier 曲線の分割と同じ方法で求められます。

     



次のプロセスで V 方向での分割が得られます。

   1)[ Bp ] の転置マトリックスをとり、

   2)区分マトリックス [ Df ] あるいは [ Dr ] を求め、

   3)それを乗じて得られたマトリックスの転置マトリックスをとる

     




以上を Python で書き表せば、

def matrix_cross_20(A, B) :	#  4行, 4列 のベクトルマトリックス A と、 4行, 4列 のマトリックス B  との積
	C = [
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(4) :
				for m in range(3) :
					C[i][j][m] = C[i][j][m] + A[i][k][m]*B[k][j]
	return C




def matrix_cross_20t(A, B) :	#  4行, 4列 のベクトルマトリックス A の転置マトリックス と、 4行, 4列 のマトリックス B  との積 の転置マトリックス
	C = [
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		]

	for i in range(4) :
		for j in range(4) :
			for k in range(4) :
				for m in range(3) :
					C[j][i][m] = C[j][i][m] + A[k][i][m]*B[k][j]
	return C


	
def get_Df_matrix(t) :		#  parameter 0〜 t 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	1,	1 - t,	(1 - t)**2,	(1 - t)**3	),
		(	0,	t,	2*t*(1 - t),	3*t*(1 - t)**2	),
		(	0,	0,	t**2,		3*t**2*(1 - t)	),
		(	0,	0,	0,		t**3		)
		)
		
		
def get_Dr_matrix(t) :		#  parameter t〜 1 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	(1 - t)**3,	0,		0,	0	),
		(	3*t*(1 - t)**2,	(1 - t)**2,	0,	0	),
		(	3*t**2*(1 - t),	2*t*(1 - t),	1 - t,	0	),
		(	t**3,		t**2,		t,	1	)
		)
		

def subdivide_bezier_patch_U_1(Bp, t) :		#  bezier patch Bp を、Uパラメータ 0 〜 t で分割した patch を返す
	Df = get_Df_matrix(t)
	return matrix_cross_20(Bp, Df)


	
def subdivide_bezier_patch_U_2(Bp, t) :		#  bezier patch Bp を、Uパラメータ t 〜1 で分割した patch を返す
	Dr = get_Dr_matrix(t)
	return matrix_cross_20(Bp, Dr)	



def subdivide_bezier_patch_V_1(Bp, s) :		#  bezier patch Bp を、Vパラメータ 0 〜 s で分割した patch を返す
	Df = get_Df_matrix(s)
	return matrix_cross_20t(Bp, Df)


	
def subdivide_bezier_patch_V_2(Bp, s) :		#  bezier patch Bp を、Vパラメータ s 〜1 で分割した patch を返す
	Dr = get_Dr_matrix(s)
	return matrix_cross_20t(Bp, Dr)	
	

		
def subdivide_bezier_patch(Bp, t1, t2, s1, s2) :		#  bezier patch bP を、区間 t1〜t2, s1〜s2 で分割した bezier patch を返す
	if t1 == 0 :
		Bp1 = subdivide_bezier_patch_U_1(Bp, t2)			#  bezier patch Bp を、Uパラメータ 0〜 t2 で分割した patch
	elif t2 == 1 :
		Bp1 = subdivide_bezier_patch_U_2(Bp, t1)			#  bezier patch Bp を、Uパラメータ t1〜 1 で分割した patch
	else :
		Bp1 = subdivide_bezier_patch_U_2(Bp, t1)			#  bezier patch Bp を、Uパラメータ t1〜 1 で分割した patch
		Bp1 = subdivide_bezier_patch_U_1(Bp1, (t2 - t1)/(1 - t1))	#  bezier patch Bp を、Uパラメータ t1〜 t2 で分割した patch
		
	if s1 == 0 :
		return subdivide_bezier_patch_V_1(Bp1, s2)			#  bezier patch Bp1 を、Vパラメータ 0〜 s2 で分割した patch
	elif s2 == 1 :
		return subdivide_bezier_patch_V_2(Bp1, s1)			#  bezier patch Bp1 を、Vパラメータ s1〜 1 で分割した patch
	else :
		Bp2 = subdivide_bezier_patch_V_2(Bp1, s1)			#  bezier patch Bp1 を、Vパラメータ s1〜 1 で分割した patch
		Bp2 = subdivide_bezier_patch_V_1(Bp2, (s2 - s1)/(1 - s1))	#  bezier patch Bp1 を、Vパラメータ s1〜 s2 で分割した patch
		return Bp2



以下のスクリプトではサンプル自由曲面を指定するパラメータ値の位置で分割します。


自由曲面の分割に対して、insert_control_point( ) を用いたり、マニュアルで行った場合、Shade の Bezier Patch 組み立てのルールから、


   再分割を繰り返すと分割線すらオリジナルの曲面をトレースしなくなる


となってしまいますが、上式を用いたスクリプトでは分割後の Bezier Patch の制御点座標を本来の定義のまま保ちますので、再分割を繰り返しても、


   その分割線は常にオリジナルの曲面をトレースし続ける


ことが可能になります。


          左 : 分割前    中 : insert_control_point( ) 使用  右 : マトリックス演算使用

     


#  サンプル自由曲面を指定する Bezier パラメータ値で分割する。


####    入力     #################################################
p1= ((0, 1000, 0), None, (0, 1000, 550), None, None)
p2 = ((0, 0, 1000), (0, 550, 1000), None, None, (550, 0, 1000))
p3 = ((0, 1000, 0), None, (550, 1000, 0), None, None)
p4 = ((1000, 0, 0), (1000, 550, -0), None, (1000, 0, 550), None)

bS = ((p1, p2), (p3, p4))				#  サンプル自由曲面			

paramL = [[0.225, 0.5, 0.775], [0.225, 0.5, 0.775]]	#  分割位置を指定する Bezier パラメータ			
############################################################


#----------  vector 演算  ------------------------------------------------------------------------------------

def sub_vec(u, v):	#  ベクトル減算
	return (u[0] - v[0], u[1] - v[1], u[2] - v[2])
#-------------------------------------------------------------------------------------------------------------


#----------  matrix 演算  ------------------------------------------------------------------------------------

def matrix_cross_20(A, B) :	#  4行, 4列 のベクトルマトリックス A と、 4行, 4列 のマトリックス B  との積
	C = [
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(4) :
				for m in range(3) :
					C[i][j][m] = C[i][j][m] + A[i][k][m]*B[k][j]
	return C




def matrix_cross_20t(A, B) :	#  4行, 4列 のベクトルマトリックス A の転置マトリックス と、 4行, 4列 のマトリックス B  との積 の転置マトリックス
	C = [
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		[[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]],
		]

	for i in range(4) :
		for j in range(4) :
			for k in range(4) :
				for m in range(3) :
					C[j][i][m] = C[j][i][m] + A[k][i][m]*B[k][j]
	return C
#-------------------------------------------------------------------------------------------------------------


#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------
	
def get_Df_matrix(t) :		#  parameter 0〜 t 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	1,	1 - t,	(1 - t)**2,	(1 - t)**3	),
		(	0,	t,	2*t*(1 - t),	3*t*(1 - t)**2	),
		(	0,	0,	t**2,		3*t**2*(1 - t)	),
		(	0,	0,	0,		t**3		)
		)
		
		
def get_Dr_matrix(t) :		#  parameter t〜 1 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	(1 - t)**3,	0,		0,	0	),
		(	3*t*(1 - t)**2,	(1 - t)**2,	0,	0	),
		(	3*t**2*(1 - t),	2*t*(1 - t),	1 - t,	0	),
		(	t**3,		t**2,		t,	1	)
		)
		

def subdivide_bezier_patch_U_1(Bp, t) :		#  bezier patch Bp を、Uパラメータ 0 〜 t で分割した patch を返す
	Df = get_Df_matrix(t)
	return matrix_cross_20(Bp, Df)


	
def subdivide_bezier_patch_U_2(Bp, t) :		#  bezier patch Bp を、Uパラメータ t 〜1 で分割した patch を返す
	Dr = get_Dr_matrix(t)
	return matrix_cross_20(Bp, Dr)	



def subdivide_bezier_patch_V_1(Bp, s) :		#  bezier patch Bp を、Vパラメータ 0 〜 s で分割した patch を返す
	Df = get_Df_matrix(s)
	return matrix_cross_20t(Bp, Df)


	
def subdivide_bezier_patch_V_2(Bp, s) :		#  bezier patch Bp を、Vパラメータ s 〜1 で分割した patch を返す
	Dr = get_Dr_matrix(s)
	return matrix_cross_20t(Bp, Dr)	
	

		
def subdivide_bezier_patch(Bp, t1, t2, s1, s2) :		#  bezier patch bP を、区間 t1〜t2, s1〜s2 で分割した bezier patch を返す
	if t1 == 0 :
		Bp1 = subdivide_bezier_patch_U_1(Bp, t2)			#  bezier patch Bp を、Uパラメータ 0〜 t2 で分割した patch
	elif t2 == 1 :
		Bp1 = subdivide_bezier_patch_U_2(Bp, t1)			#  bezier patch Bp を、Uパラメータ t1〜 1 で分割した patch
	else :
		Bp1 = subdivide_bezier_patch_U_2(Bp, t1)			#  bezier patch Bp を、Uパラメータ t1〜 1 で分割した patch
		Bp1 = subdivide_bezier_patch_U_1(Bp1, (t2 - t1)/(1 - t1))	#  bezier patch Bp を、Uパラメータ t1〜 t2 で分割した patch
		
	if s1 == 0 :
		return subdivide_bezier_patch_V_1(Bp1, s2)			#  bezier patch Bp1 を、Vパラメータ 0〜 s2 で分割した patch
	elif s2 == 1 :
		return subdivide_bezier_patch_V_2(Bp1, s1)			#  bezier patch Bp1 を、Vパラメータ s1〜 1 で分割した patch
	else :
		Bp2 = subdivide_bezier_patch_V_2(Bp1, s1)			#  bezier patch Bp1 を、Vパラメータ s1〜 1 で分割した patch
		Bp2 = subdivide_bezier_patch_V_1(Bp2, (s2 - s1)/(1 - s1))	#  bezier patch Bp1 を、Vパラメータ s1〜 s2 で分割した patch
		return Bp2
#-------------------------------------------------------------------------------------------------------------


#----------  Bezier Patch  ---------------------------------------------------------------------------------

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])


def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)


def get_standard_length(bP) :	#  bP の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if bP[i][j] != None :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
	
	
	
def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])
	
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	
	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	
#-------------------------------------------------------------------------------------------------------------



#----------  出力  ----------------------------------------------------------------------------------------	

def write_surface(bS) :		#  自由曲面出力
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	for bL in bS :
		scene.begin_line(None, False)
		for p in bL :
			scene.append_point(p[0], p[1], p[2], p[3], p[4])
		scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	
	
	return ob 	#  作成した自由曲面を返す
	
	
	
def write_surface_of_bezier_patch(bP) :	 #  Bezier Patch bP で表された自由曲面を出力
	p1= (bP[0][0], None, bP[0][1], None, bP[1][0])
	p2 = (bP[0][3], bP[0][2], None, None, bP[1][3])
	p3 = (bP[3][0], None, bP[3][1], bP[2][0], None)
	p4 = (bP[3][3], bP[3][2], None, bP[2][3], None)

	bS = ((p1, p2), (p3, p4))
	
	return write_surface(bS)
	

	
def divide_surface(ob, paramL) :	#  自由曲面 ob を paramL 位置で分割

	ob.copy_object(((1,0,0,0), (0,1,0,0), (0,0,1,0), (0,0,0,1)))
	ob2 = ob.bro

	obL = []
	
	for i in range(2) :
		ob3 = ob2.son.bro
		paramL[i].reverse()
		for j in range(ob2.number_of_sons) :
			prev = 1.
			for p in paramL[i] :
				ob3.insert_control_point(p/prev)
				prev = p
			ob3 = ob3.bro
		paramL[i].reverse()
		ob2.switch()

	return ob2
#-------------------------------------------------------------------------------------------------------------	

	
def main(bS, paramL) :			#  paramL は sort されていることを前提とする
	scene.inhibit_update()

	scene.create_part()
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob = write_surface(bS)			#  サンプル自由曲面を出力
	bP = get_bezier_patch(ob)		#  サンプル自由曲面の Bezier Patch
	ob2 = divide_surface(ob, paramL)	#  サンプル自由曲面の複製を paramL 位置で分割
	ob2.name = 'insert_control_point( ) 使用'
	
	ob2.activate()
	scene.create_part()
	ob3 = ob2.bro
	ob3.activate()
	
	for i in range(2) :			#  paramL の再構成
		if paramL[i][0] != 0 :		
			paramL[i].insert(0, 0)
		if paramL[i][len(paramL[i]) - 1] != 1 :
			paramL[i].append(1)
	
	for i in range(len(paramL[0]) - 1) :
		for j in range(len(paramL[1]) - 1) :
			bP1 = subdivide_bezier_patch(bP, paramL[0][i], paramL[0][i + 1], paramL[1][j], paramL[1][j + 1])	#  bP を分割した Bezier Patchを得る
			write_surface_of_bezier_patch	(bP1)

		
	ob0.activate()
	ob0.reset_transformation()
		
	scene.allow_update()
	scene.force_update()
		
	
scene = xshade.scene()
main(bS, paramL)
投稿数: 189
投稿日時: 2011-11-24 21:11
Re: 線形状(ベジェ曲線)の分割について教えてください
   
< Bezier 曲線の区分曲線 >


Bezier 曲線の分割については、このスレッドの最初の方で De Castelijau のアルゴリズム について触れましたが、階差表記を用いると Bezier パラメータ 0 〜 t の区間で区分される Bezier 曲線は次のようにとても簡単に表されます。

     



そこでこれを用いて Bezier 曲線の分割を マトリックス計算を用いた求め方で表してみます。

   1)階差表記から区分曲線を求める

     


   2)通常表記の制御点座標を階差表記に変換する

     


   3)階差表記の制御点座標を通常表記に変換する

     


   ここで、通常表記を階差表記に変換するマトリックスと 階差表記を通常表記に変換するマトリックスとは
   互いに逆マトリックスの関係にあります。

     


   以上から、通常表記で表された制御点座標 [ P ] のパラメータ区間 0 〜 t における
   区分曲線の制御点座標 [ Q ]( 通常表記 )は次のように表されます。

     

   ここでマトリックス計算では結合法則が成立しますから

     [ A ] × [ b ] × [ C ] = ( [ A ] × [ b ] ) × [ C ] = [ A ] × ( [ b ] × [ C ] )

   よって、

     


   上式の( )内を展開すれば

     




一方、次のプロセスで Bezier パラメータ t 〜 1 で区分される曲線が求められます。

   1)[ P ] の並び順を反転し、

   2)s = 1 - t として区分マトリックス [ Df ] 求め、

   3)それを乗じて得られたマトリックスの並び順を再び反転する

   
   ここで、並び順反転のためのマトリックスを [ R ] とすれば、

     


   よって、通常表記で表された制御点座標 [ P ] のパラメータ区間 t 〜 1 における
   区分曲線の制御点座標 [ Q ]( 通常表記 )は次のように表されます。

     


          [ R ] × [ A ] × [ R ] は、[ A ] の行と列の順番を反転させることになりますから、

     

          s を 1 - t に置き換えれば、

     




以上を Python で書き表せば、



def matrix_cross_10(A, B) :	#  1行, 4列 のベクトルマトリックス A と、 4行, 4列 のマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j][k]*B[j][i]
	return C


	
def get_Df_matrix(t) :		#  parameter 0〜 t 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	1,	1 - t,	(1 - t)**2,	(1 - t)**3	),
		(	0,	t,	2*t*(1 - t),	3*t*(1 - t)**2	),
		(	0,	0,	t**2,		3*t**2*(1 - t)	),
		(	0,	0,	0,		t**3		)
		)
		
		
def get_Dr_matrix(t) :		#  parameter t〜 1 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	(1 - t)**3,	0,		0,	0	),
		(	3*t*(1 - t)**2,	(1 - t)**2,	0,	0	),
		(	3*t**2*(1 - t),	2*t*(1 - t),	1 - t,	0	),
		(	t**3,		t**2,		t,	1	)
		)
		

def subdivide_bezier_line_1(p, t) :		#  bezier line p を、区間 0〜 t で分割した線形状を返す
	Df = get_Df_matrix(t)
	return matrix_cross_10(p, Df)


	
def subdivide_bezier_line_2(p, t) :		#  bezier line p を、区間 t〜 1 で分割した線形状を返す
	Dr = get_Dr_matrix(t)
	return matrix_cross_10(p, Dr)	



def subdivide_bezier_line(p, t1, t2) :	#  bezier line p を、区間 t1〜 t2 で分割した線形状を返す
	if t1 == 0 :
		return subdivide_bezier_line_1(p, t2)			#  bezier line p を、区間 0〜 t2 で分割した線形状を返す
	elif t2 == 1 :
		return subdivide_bezier_line_2(p, t1)			#  bezier line p を、区間 t1〜 1 で分割した線形状を返す
	else :
		bL = subdivide_bezier_line_2(p, t1)			#  bezier line p を、区間 t1〜 1 で分割した線形状を得る
		bL = subdivide_bezier_line_1(bL, (t2 - t1)/(1 - t1))	#  さらに 区間 t1〜 t2 で分割した線形状を得る
		return bL





以下のスクリプトではサンプル線形状を指定するパラメータ値の位置で分割します。

#  サンプル線形状を指定する Bezier パラメータ値で分割する。


####    入力     #################################################
p = (				#  サンプル線形状の Bezier 制御点座標
	(-2000, 1000, 0),		
	(0, 3000, 0),
	(1700, 1000, 0),
	(2000, -1000, 0)
	)			

paramL = [0.3, 0.7]		#  分割位置を指定する Bezier パラメータ			
############################################################



#----------  matrix 演算  ------------------------------------------------------------------------------------

def matrix_cross_10(A, B) :	#  1行, 4列 のベクトルマトリックス A と、 4行, 4列 のマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j][k]*B[j][i]
	return C
#-------------------------------------------------------------------------------------------------------------


#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------
	
def get_Df_matrix(t) :		#  parameter 0〜 t 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	1,	1 - t,	(1 - t)**2,	(1 - t)**3	),
		(	0,	t,	2*t*(1 - t),	3*t*(1 - t)**2	),
		(	0,	0,	t**2,		3*t**2*(1 - t)	),
		(	0,	0,	0,		t**3		)
		)
		
		
def get_Dr_matrix(t) :		#  parameter t〜 1 区間の bezier 制御点座標を与えるための マトリックス
	return (
		(	(1 - t)**3,	0,		0,	0	),
		(	3*t*(1 - t)**2,	(1 - t)**2,	0,	0	),
		(	3*t**2*(1 - t),	2*t*(1 - t),	1 - t,	0	),
		(	t**3,		t**2,		t,	1	)
		)
		

def subdivide_bezier_line_1(p, t) :		#  bezier line p を、区間 0〜 t で分割した線形状を返す
	Df = get_Df_matrix(t)
	return matrix_cross_10(p, Df)


	
def subdivide_bezier_line_2(p, t) :		#  bezier line p を、区間 t〜 1 で分割した線形状を返す
	Dr = get_Dr_matrix(t)
	return matrix_cross_10(p, Dr)	



def subdivide_bezier_line(p, t1, t2) :	#  bezier line p を、区間 t1〜 t2 で分割した線形状を返す
	if t1 == 0 :
		return subdivide_bezier_line_1(p, t2)			#  bezier line p を、区間 0〜 t2 で分割した線形状を返す
	elif t2 == 1 :
		return subdivide_bezier_line_2(p, t1)			#  bezier line p を、区間 t1〜 1 で分割した線形状を返す
	else :
		bL = subdivide_bezier_line_2(p, t1)			#  bezier line p を、区間 t1〜 1 で分割した線形状を得る
		bL = subdivide_bezier_line_1(bL, (t2 - t1)/(1 - t1))	#  さらに 区間 t1〜 t2 で分割した線形状を得る
		return bL	
#-------------------------------------------------------------------------------------------------------------


#----------  出力  ----------------------------------------------------------------------------------------	

def write_line(p) :	 #  線形状出力
			
	scene.begin_creating()
	scene.begin_line(None, False)
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
	scene.end_line()	
	scene.end_creating()
	
	ob = scene.active_shape()
	
	return ob
	

	
def divide_line(ob, paramL) :	#  線形状 ob を paramL 位置で分割

	ob.copy_object(((1,0,0,0), (0,1,0,0), (0,0,1,0), (0,0,0,1)))
	ob2 = ob.bro
	
	paramL.reverse()	#  リスト順を反転
	prev = 1.
	for p in paramL :
		ob2.insert_control_point(p/prev)
		prev = p
	
	paramL.reverse()	#  リスト順を復帰
	
	return ob2
#-------------------------------------------------------------------------------------------------------------	

	
def main(p, paramL) :			#  paramL は sort されていることを前提とする
	scene.inhibit_update()

	scene.create_part()
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob = write_line(p)		#  サンプル線形状を描く
	ob2 = divide_line(ob, paramL)		#  サンプル線形状の複製を paramL 位置で分割
	
	ob2.activate()
	scene.create_part()
	
	
	if paramL[0] != 0 :		#  paramL の再構成
		paramL.insert(0, 0)
	if paramL[len(paramL) - 1] != 1	:
		paramL.append(1)
		
		
	for i in range(len(paramL) - 1) :
		bL = subdivide_bezier_line(p, paramL[i], paramL[i + 1])		#  bezier line p を、区間 para1〜 para2 で分割した線形状を得る
		ob = write_line(bL)						#  分割線形状を描く
		
	ob0.activate()
	ob0.reset_transformation()
	
	if not scene.is_modify_mode :
		scene.enter_modify_mode()
		
	scene.allow_update()
	scene.force_update()
		
	
scene = xshade.scene()
main(p, paramL)
投稿数: 189
投稿日時: 2011-11-12 22:52
Re: 線形状(ベジェ曲線)の分割について教えてください
  
< Bezier 曲面の U/V 接線 - 4/4 >


以上で Shade の Bezier Patch の U/V 接線ベクトルの求め方がブラッシュアップされたので、最後に自由曲面の U/V接線ベクトルと面法線ベクトルを出力するスクリプトを書いておきます。


#  U/V 接線と面法線ベクトルを求める



####    入力     #################################################
#  bS :  サンプル自由曲面


p1= ((0, 1000, 0), None, None, None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (600, 0, 1000))
p3 = ((0, 1000, 0), None, None, None, None)
p4 = ((1000, 0, 0), (1000, 500,0), None, (1000, 0, 600), None)

bS = ((p1, p2), (p3, p4))


#  paramL :  U/V 接線と法線のパラメータ値 (t, s)
paramL = []
for i in range(6) :
	for j in range(6) :
		paramL.append((i/5., j/5.))

# vlen : 法線ベクトルの長さ
vlen = 300
############################################################



#----------  vector  演算  -----------------------------------------------------------------------------

def add_vec(u, v):	#  ベクトル加算
	return (u[0] + v[0], u[1] + v[1], u[2] + v[2])
	
	
def sub_vec(u, v):	#  ベクトル減算
	return (u[0] - v[0], u[1] - v[1], u[2] - v[2])
	

def times_vec(v, a):	# ベクトル乗算
	return (v[0]*a, v[1]*a, v[2]*a)
	
	
def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	
	

def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)


def cross_vec(v1, v2) :		#  ベクトルの外積
	return[v1[1]*v2[2] - v1[2]*v2[1], v1[2]*v2[0] - v1[0]*v2[2], v1[0]*v2[1] - v1[1]*v2[0]]
	
#-------------------------------------------------------------------------------------------------------------


#----------  出力  ---------------------------------------------------------------------------------	

def write_surface(bS) :		#  自由曲面出力
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	
	for bL in bS :
		scene.begin_line(None, False)
		for p in bL :
			scene.append_point(p[0], p[1], p[2], p[3], p[4])
		scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	
	return ob 	#  作成した自由曲面を返す
		


	
def write_line(p, name) :	#  線形状出力	p : 制御点座標
	
	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
		
	
	
	
def write_tangent(p, v, vlen, name) :		#  接線ベクトル出力

	q = [0, 0, 0]
	for i in range(3) :
		q[i] = p[i] + vlen*v[i]

	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p, None, None, None, None)
	scene.append_point(q, None, None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す

#-------------------------------------------------------------------------------------------------------------



#----------  Bezier Patch  ---------------------------------------------------------------------------------

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])


def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)


def get_standard_length(bP) :	#  bP の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if bP[i][j] != None :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
	
	
	
def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])
	
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	
	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	
#-------------------------------------------------------------------------------------------------------------
	
	
	

#----------  matrix 演算  ------------------------------------------------------------------------------------
	
def matrix_cross_1(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[j][i][k]
	return C
	
	

def matrix_cross_1t(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B の転置マトリックス  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[i][j][k]	
	return C



def matrix_cross_2(A, B) :	#  4列 のベクトルトリックス A と、 4列 のマトリックス B の転置マトリックス  との積

	C = [0, 0, 0]
	
	for i in range(4) :
		for k in range(3) :
			C[k] = C[k] + A[i][k]*B[i]	
				
	return C
	


def matrix_cross_0(A, B, C) :	#  4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B と、 4列のマトリックス C の転置マトリックス  との積
	Vs = matrix_cross_1(A, B)
	return matrix_cross_2(Vs, C)	
#-------------------------------------------------------------------------------------------------------------




#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------

def get_T_matrix(t) :
	return ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	
	
def get_T_difference_matrix(t) :
	return (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	

def get_point_from_bezier_patch(Bp, t, s) :		#  指定する U, V パラメータ t, s における座標値を返す
	Tu = get_T_matrix(t)
	Tv = get_T_matrix(s)
	return matrix_cross_0(Tv, Bp, Tu)

	
def get_U_line_from_bezier_patch(Bp, s) :	#  Bezier Patch bP の Vパラメータ s における U線形状の制御点座標を返す
	Tv = get_T_matrix(s)				
	return matrix_cross_1(Tv, Bp)
	

def get_V_line_from_bezier_patch(Bp, t) :	#  Bezier Patch bP の Uパラメータ t における V線形状の制御点座標を返す
	Tu = get_T_matrix(t)				
	return matrix_cross_1t(Tu, Bp)
			

def get_std_L(p) :				#  Bezier line  の bounding box 基準の 基準長を返す
	a = [p[0][0], p[0][1], p[0][2]]
	b = [p[0][0], p[0][1], p[0][2]]
	for i in range(1, 4) :
		for k in range(3) :
			a[k] = max(a[k], p[i][k])
			b[k] = min(b[k], p[i][k])
	return a[0] + a[1] + a[2] - b[0] - b[1] - b[2]
	


def check_tangent(v, p, t) :			#  v = (0, 0, 0) となって計算上接線が求められない場合の処理				
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = sub_vec(p[2], p[0])
		else :
			v = sub_vec(p[3], p[1])
		v= normalize_vec(v)
			
	if v == (0, 0, 0) :
		v = sub_vec(p[3], p[0])
		v= normalize_vec(v)		
	return v	

	

def get_tangent(p, t, std_L) :	
	T_ = get_T_difference_matrix(t)
	v = matrix_cross_2(p, T_)
	
	if abs_vec(v)/std_L>= 1.5e-9 :
		v = normalize_vec(v)
	else :
		tt = round(t)			#  parameter t を 0 または 1 として再計算 ( ここに入ってくる場合 t は 0 または 1 に非常に近い )		
		if tt != t :
			T_ = get_T_difference_matrix(tt)
			v = matrix_cross_2(p, T_)
		v = normalize_vec(v)
		v = check_tangent(v, p, tt)
	return v



def get_U_tangent(Bp, t, s) :					#  Bezier Patch bP の U, Vパラメータ t, s における U方向接線単位ベクトルを返す

	P1 = (Bp[0])						#  Bezier Patch を構成する最初の U方向線形状の制御点座標
	L1 = get_std_L(P1)					#  P1 の基準長
	P2 = (Bp[3])						#  Bezier Patch を構成する二番目の U方向線形状の制御点座標
	L2 = get_std_L(P2)					#  P2 の基準長
		
	if (L1 > 0 and L2 > 0) or (L1 == 0 and L2 == 0) :
		Pu = get_U_line_from_bezier_patch(Bp, s)	#  Bezier Patch Bp の Vパラメータ s における U線形状の制御点座標
		L = get_std_L(Pu)				#  Pu の基準長
		if L == 0 :
			return (0, 0, 0)				#  接線ベクトルが求まらない
	elif L1 == 0 :						#  P1 のみが極
		Pu = P2
		L = L2
	else :							#  P2 のみが極
		Pu = P1
		L = L1
	v = get_tangent(Pu, t, L)	

	return v


	
def get_V_tangent(Bp, t, s) :					#  Bezier Patch bP の U, Vパラメータ t, s における V方向接線単位ベクトルを返す

	P1 = (Bp[0][0], Bp[1][0], Bp[2][0], Bp[3][0])		#  Bezier Patch を構成する最初の V方向線形状の制御点座標
	L1 = get_std_L(P1)					#  P1 の基準長
	P2 = (Bp[0][3], Bp[1][3], Bp[2][3], Bp[3][3])		#  Bezier Patch を構成する二番目の V方向線形状の制御点座標
	L2 = get_std_L(P2)					#  P2 の基準長
		
	if (L1 > 0 and L2 > 0) or (L1 == 0 and L2 == 0) :
		Pv = get_V_line_from_bezier_patch(Bp, t)	#  Bezier Patch Bp の Uパラメータ t における V線形状の制御点座標
		L = get_std_L(Pv)				#  Pv の基準長
		if L == 0 :
			return (0, 0, 0)				#  接線ベクトルが求まらない	
	elif L1 == 0 :						#  P1 のみが極
		Pv = P2
		L = L2
	else :							#  P2 のみが極
		Pv = P1
		L = L1
	
	v = get_tangent(Pv, s, L)	

	return v
	
	
	
def get_normal(Bp,t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における法線単位ベクトルを返す
	Vu = get_U_tangent(Bp, t, s)			#  Bezier Patch Bp の U, Vパラメータ t, s における U方向接線単位ベクトル
	Vv = get_V_tangent(Bp, t, s)			#  Bezier Patch Bp の U, Vパラメータ t, s における V方向接線単位ベクトル
	v = cross_vec(Vu, Vv)
	return normalize_vec(v)
#-------------------------------------------------------------------------------------------------------------	



	
	
def main(bS, paramL, vlen) :

	scene.inhibit_update()
	
	scene.create_part()			#  全体を格納するパート
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob = write_surface(bS)			#  サンプル自由曲面を出力	
	Bp = get_bezier_patch(ob)		#  サンプル自由曲面の Bezier Patch を求める
	
	scene.create_part('U/V接線と法線')
	ob1 = scene.active_shape()
	ob1.disclosed = False
	
	for (t, s) in paramL :
		scene.create_part(str((t, s)))
		ob2 = scene.active_shape()
		ob2.disclosed = False
					
		p = get_point_from_bezier_patch(Bp, t, s)	#  指定する U, V パラメータ値における座標を取得
		
		Vu = get_U_tangent(Bp, t, s)			#  U, Vパラメータ t, s における U方向接線単位ベクトルを取得
		ob3 = write_tangent(p, Vu, vlen, 'U')		#  接線ベクトルを出力
	
		Vv = get_V_tangent(Bp, t, s)			#  U, Vパラメータ t, s における V方向接線単位ベクトルを取得		
		ob4 = write_tangent(p, Vv, vlen, 'V')		#  接線ベクトルを出力
	
		Vn= get_normal(Bp, t, s)			#  U, Vパラメータ t, s における法線単位ベクトルを取得		
		ob5 = write_tangent(p, Vn, vlen, 'N')		#  法線ベクトルを出力

		ob2.activate()
		
	ob1.activate()
	
	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	
	
	
scene = xshade.scene()
main(bS, paramL, vlen)
投稿数: 189
投稿日時: 2011-11-12 22:51
Re: 線形状(ベジェ曲線)の分割について教えてください
 
< Bezier 曲面の U/V 接線 - 3/4 >


球の極における経線方向( V方向)の接線ベクトルを求めるような場合、基準となる線形状に長さがないので、そのまま計算しても求められません。


球の極における経線方向( V方向)の接線ベクトルの求め方について、考えてみます。


   


     L1 = ( P0, P1, P2, P3 ) Bezier 曲線制御点座標

     L2 = ( Q0, Q1, Q2, Q3 ) Bezier 曲線制御点座標

     L1 が U 方向であるとし、U方向の parameter : t, V方向の parameter : s = 0 で与えられる

     Bezier Patch 上の点の V 方向接線の向きは、

     曲線 L2, L1 上の parameter t で与えられる点を A, B とすれば、

       tangent(t) = A - B
             = ((1 - t)^3*Q0 + 3*t*(1 - t)^2*Q1 + 3*t^2*(1 - t)*Q2 + t^3*Q3)
             - ((1 - t)^3*P0 + 3*t*(1 - t)^2*P1 + 3*t^2*(1 - t)*P2 + t^3*P3)
             = (1 - t)^3*(Q0 - P0) + 3*t*(1 - t)^2*(Q1 - P1) + 3*t^2*(1 - t)*(Q2 - P2) + t^3*(Q3 - P3)

     ここで、

       Q0 = P0
       Q0 - P0 = 0

     さらに Shade の Bezier Patch の作り方から、

       Q2 - P2 = Q3 - P3
       Q1 - P1 (Q2 - P2)/2

     つまり、s = 0 で与えられる点の V 方向接線の向きは t の値には関係なく、Q3 - P3 で与えられることが解ります。




次のスクリプトでは、上記の考え方を組み入れて接線ベクトルを求めるプログラムを書き直してみました。


#  V接線ベクトルを求める



####    入力     #################################################
List = []

p1= ((0, 1000, 0), None, None, None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (600, 0, 1000))
p3 = ((0, 1000, 0), None, None, None, None)
p4 = ((1000, 0, 0), (1000, 500,0), None, (1000, 0, 600), None)

paramL2 = []
paramL2.append((0, 0.7))
paramL2.append((0.5, 0.7))
paramL2.append((1, 0.7))

List.append((((p1, p2), (p3, p4)), paramL2))


p1= ((0, 1000, 0), None, None, None, (0, 1000, 0.000001))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (600, 0, 1000))
p3 = ((0, 1000, 0), None, None, (0.000001, 1000, 0), None)
p4 = ((1000, 0, 0), (1000, 500,0), None, (1000, 0, 600), None)

paramL2 = []
paramL2.append((0, 0))
paramL2.append((0, 0.5))
paramL2.append((0, 1))

List.append((((p1, p2), (p3, p4)), paramL2))


p1= ((-1000, 1000, 0), None, None, None, None)
p2 = ((1000, -1000, 0), None, None, None, None)
p3 = ((1000, 1000, 0), None, None, None, None)
p4 = ((-1000, -1000, 0), None, None, None, None)

paramL2 = []
paramL2.append((0.5, 0))
paramL2.append((0.5, 0.5))
paramL2.append((0.5, 1))

List.append((((p1, p2), (p3, p4)), paramL2))


# vlen : 接線ベクトルの長さ
vlen = 500
############################################################



#----------  vector  演算  -----------------------------------------------------------------------------

def add_vec(u, v):	#  ベクトル加算
	return (u[0] + v[0], u[1] + v[1], u[2] + v[2])
	
	
def sub_vec(u, v):	#  ベクトル減算
	return (u[0] - v[0], u[1] - v[1], u[2] - v[2])
	

def times_vec(v, a):	# ベクトル乗算
	return (v[0]*a, v[1]*a, v[2]*a)
	
	
def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	
	

def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)

#-------------------------------------------------------------------------------------------------------------


#----------  出力  ---------------------------------------------------------------------------------	

def write_surface(bS) :		#  自由曲面出力
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	
	for bL in bS :
		scene.begin_line(None, False)
		for p in bL :
			scene.append_point(p[0], p[1], p[2], p[3], p[4])
		scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	
	return ob 	#  作成した自由曲面を返す
		


	
def write_line(p, name) :	#  線形状出力	p : 制御点座標
	
	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
		
	
	
	
def write_tangent(p, v, vlen, name) :		#  接線ベクトル出力

	q = [0, 0, 0]
	for i in range(3) :
		q[i] = p[i] + vlen*v[i]

	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p, None, None, None, None)
	scene.append_point(q, None, None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す

#-------------------------------------------------------------------------------------------------------------



#----------  Bezier Patch  ---------------------------------------------------------------------------------

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])


def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)


def get_standard_length(bP) :	#  bP の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if bP[i][j] != None :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
	
	
	
def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])
	
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	
	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	
#-------------------------------------------------------------------------------------------------------------
	
	
	

#----------  matrix 演算  ------------------------------------------------------------------------------------
	
def matrix_cross_1(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[j][i][k]
	return C
	
	

def matrix_cross_1t(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B の転置マトリックス  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[i][j][k]	
	return C



def matrix_cross_2(A, B) :	#  4列 のベクトルトリックス A と、 4列 のマトリックス B の転置マトリックス  との積

	C = [0, 0, 0]
	
	for i in range(4) :
		for k in range(3) :
			C[k] = C[k] + A[i][k]*B[i]	
				
	return C
	


def matrix_cross_0(A, B, C) :	#  4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B と、 4列のマトリックス C の転置マトリックス  との積
	Vs = matrix_cross_1(A, B)
	return matrix_cross_2(Vs, C)	
#-------------------------------------------------------------------------------------------------------------




#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------

def get_T_matrix(t) :
	return ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	
	
def get_T_difference_matrix(t) :
	return (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	

def get_point_from_bezier_patch(Bp, t, s) :		#  指定する U, V パラメータ t, s における座標値を返す
	Tu = get_T_matrix(t)
	Tv = get_T_matrix(s)
	return matrix_cross_0(Tv, Bp, Tu)

	
def get_U_line_from_bezier_patch(Bp, s) :	#  Bezier Patch bP の Vパラメータ s における U線形状の制御点座標を返す
	Tv = get_T_matrix(s)				
	return matrix_cross_1(Tv, Bp)
	

def get_V_line_from_bezier_patch(Bp, t) :	#  Bezier Patch bP の Uパラメータ t における V線形状の制御点座標を返す
	Tu = get_T_matrix(t)				
	return matrix_cross_1t(Tu, Bp)
				
	
def get_std_L(p) :				#  Bezier line  の bounding box 基準の 基準長を返す
	a = [p[0][0], p[0][1], p[0][2]]
	b = [p[0][0], p[0][1], p[0][2]]
	for i in range(1, 4) :
		for k in range(3) :
			a[k] = max(a[k], p[i][k])
			b[k] = min(b[k], p[i][k])
	return a[0] + a[1] + a[2] - b[0] - b[1] - b[2]
	


def check_tangent(v, p, t) :			#  v = (0, 0, 0) となって計算上接線が求められない場合の処理				
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = sub_vec(p[2], p[0])
		else :
			v = sub_vec(p[3], p[1])
		v= normalize_vec(v)
			
	if v == (0, 0, 0) :
		v = sub_vec(p[3], p[0])
		v= normalize_vec(v)		
	return v	

	

def get_tangent(p, t, std_L) :	
	T_ = get_T_difference_matrix(t)
	v = matrix_cross_2(p, T_)
	
	if abs_vec(v)/std_L>= 1.5e-9 :
		v = normalize_vec(v)
	else :
		tt = round(t)			#  parameter t を 0 または 1 として再計算 ( ここに入ってくる場合 t は 0 または 1 に非常に近い )		
		if tt != t :
			T_ = get_T_difference_matrix(tt)
			v = matrix_cross_2(p, T_)
		v = normalize_vec(v)
		v = check_tangent(v, p, tt)
	return v



def get_U_tangent(Bp, t, s) :					#  Bezier Patch bP の U, Vパラメータ t, s における U方向接線単位ベクトルを返す

	P1 = (Bp[0])						#  Bezier Patch を構成する最初の U方向線形状の制御点座標
	L1 = get_std_L(P1)					#  P1 の基準長
	P2 = (Bp[3])						#  Bezier Patch を構成する二番目の U方向線形状の制御点座標
	L2 = get_std_L(P2)					#  P2 の基準長
		
	if (L1 > 0 and L2 > 0) or (L1 == 0 and L2 == 0) :
		Pu = get_U_line_from_bezier_patch(Bp, s)	#  Bezier Patch Bp の Vパラメータ s における U線形状の制御点座標
		L = get_std_L(Pu)				#  Pu の基準長
		if L == 0 :
			return (0, 0, 0)				#  接線ベクトルが求まらない
	elif L1 == 0 :						#  P1 のみが極
		Pu = P2
		L = L2
	else :							#  P2 のみが極
		Pu = P1
		L = L1
	v = get_tangent(Pu, t, L)	

	return v


	
def get_V_tangent(Bp, t, s) :					#  Bezier Patch bP の U, Vパラメータ t, s における V方向接線単位ベクトルを返す

	P1 = (Bp[0][0], Bp[1][0], Bp[2][0], Bp[3][0])		#  Bezier Patch を構成する最初の V方向線形状の制御点座標
	L1 = get_std_L(P1)					#  P1 の基準長
	P2 = (Bp[0][3], Bp[1][3], Bp[2][3], Bp[3][3])		#  Bezier Patch を構成する二番目の V方向線形状の制御点座標
	L2 = get_std_L(P2)					#  P2 の基準長
		
	if (L1 > 0 and L2 > 0) or (L1 == 0 and L2 == 0) :
		Pv = get_V_line_from_bezier_patch(Bp, t)	#  Bezier Patch Bp の Uパラメータ t における V線形状の制御点座標
		L = get_std_L(Pv)				#  Pv の基準長
		if L == 0 :
			#---------------------------------------------------
			print '     接線ベクトルが求まらない'		#  report
			#---------------------------------------------------
			return (0, 0, 0)				#  接線ベクトルが求まらない	
	elif L1 == 0 :				#  P1 のみが極
		Pv = P2
		L = L2
		#---------------------------------------------------
		print '     L1 = 0',				#  report
		#---------------------------------------------------
	else :					#  P2 のみが極
		Pv = P1
		L = L1
		#---------------------------------------------------
		print '     L2 = 0',				#  report
		#---------------------------------------------------
			
	v = get_tangent(Pv, s, L)	
	#---------------------------------------------------
	print '     接線ベクトル = ', v			#  report
	#---------------------------------------------------

	return v
#-------------------------------------------------------------------------------------------------------------	


	
	
def main(List, vlen) :

	scene.inhibit_update()
	
	scene.create_part()			#  全体を格納するパート
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	for ((bS, paramL2)) in List :
		scene.create_part('No. ' + str(List.index(((bS, paramL2)))))
		ob1 = scene.active_shape()
	
		ob = write_surface(bS)			#  サンプル自由曲面を出力	
		Bp = get_bezier_patch(ob)		#  サンプル自由曲面の Bezier Patch を求める
	
		#---------------------------------------------------
		print ''
		print '- - - - - - - - - - - - - - - - - -'
		print ''
		print 'No. ' + str(List.index(((bS, paramL2))))
		print ''
		print 'V方向接線'			# report
		#---------------------------------------------------
		scene.create_part('V方向接線')
		ob2 = scene.active_shape()
		ob2.disclosed = False
	
		for (t, s) in paramL2 :
			scene.create_part('t = ' + str(t))
			ob3 = scene.active_shape()
					
			p = get_point_from_bezier_patch(Bp, t, s)	#  指定する U, V パラメータ値における座標を取得
		
			Pv = get_V_line_from_bezier_patch(Bp, t)	#  指定する Vパラメータ値における V線形状の制御点座標
			write_line(Pv, 't = ' +str(t))
		
			#---------------------------------------------------
			print 't = ', t			# report
			#---------------------------------------------------
			
			Vv = get_V_tangent(Bp, t, s)			#  U, Vパラメータ値 における V接線ベクトルを取得
			ob4 = write_tangent(p, Vv, vlen, 'V')		#  接線ベクトルを出力
	
			ob3.activate()
		ob1.activate()
		
	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	
	
	
scene = xshade.scene()
main(List, vlen)
投稿数: 189
投稿日時: 2011-11-12 22:51
Re: 線形状(ベジェ曲線)の分割について教えてください
  
< Bezier 曲面の U/V 接線 - 2/4 >


前回 Bezier 曲面の U/V 方向の接線ベクトルの求め方について述べましたが、U/V handle が出ていない場合に、数値の丸め誤差に起因する精度上の問題が顕れてきます。


これに対処するため、つぎのように考えてみます。

   接線ベクトルを求めるベースとなる線形状の bounding box 基準の基準長 std_L を求め、

   計算された単位ベクトル化する前の接線ベクトルの長さ abs(v) との比を考え、

   誤差が float の有効桁数を超えて影響し始める abs(v)/std_L の限界値 epsilon を見つけ、

   abs(v)/std_L < epsilon ならば、Bezier Parameter t が 0 に近ければ t = 0 として再計算、1 に近ければ t = 1 として再計算する
   ( t = 0 or 1 では丸め誤差の影響は最小化される )

   そして t = 0 or 1 で得られた近似値が、想定される真の値と float の有効桁数内で遜色ないものであるかを確認する



次のスクリプトでは上の考え方に基づいてサンプル自由曲面の V方向パラメータ s = 0.7 における U 方向接線ベクトルを求めています。

( U方向パラメータ t = 0.1, 1e-9, 7e-10, 5.4e-10, 3.5e-10, 1e-10, 1e-13, 1e-16, 0 )

   


このスクリプトでは、

   U方向線形状の の bounding box 基準の基準長 std_L を求め、

   計算された単位ベクトル化する前の接線ベクトルの長さ abs(v) との比をとり、

   abs(v)/std_L < 1.5e-9 ならば、Bezier Parameter t が 0 に近ければ t = 0 として再計算、1 に近ければ t = 1 として再計算する



スクリプトは次のレポートを返します。

   


レポート結果から、abs(v)/std_L = 1.5e-9 あたりで丸め誤差が最小化される近似値に切り替えても、float の有効桁数内であまり問題ない値になっていることが解ります。

サンプル自由曲面のハンドル長さや向きを変えても同じ傾向を示します。


#  U接線ベクトルを求める



####    入力     #################################################
#  bS :  サンプル自由曲面


p1= ((0, 1000, 0), None, None, None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (600, 0, 1000))
p3 = ((0, 1000, 0), None, None, None, None)
p4 = ((1000, 0, 0), (1000, 500,0), None, (1000, 0, 600), None)

bS = ((p1, p2), (p3, p4))


#  paramL1 :  U接線 パラメータ値 (t, s)
paramL1 = []
paramL1.append((1e-1, 0.7))
paramL1.append((1e-9, 0.7))
paramL1.append((7e-10, 0.7))
paramL1.append((5.4e-10, 0.7))
paramL1.append((3.5e-10, 0.7))
paramL1.append((1e-10, 0.7))
paramL1.append((1e-13, 0.7))
paramL1.append((1e-16, 0.7))
paramL1.append((0, 0.7))

#  paramL2 :  V接線 パラメータ値 (t, s)
paramL2 = []


# vlen : 接線ベクトルの長さ
vlen = 500
############################################################



#----------  vector  演算  -----------------------------------------------------------------------------

def add_vec(u, v):	#  ベクトル加算
	return (u[0] + v[0], u[1] + v[1], u[2] + v[2])
	
	
def sub_vec(u, v):	#  ベクトル減算
	return (u[0] - v[0], u[1] - v[1], u[2] - v[2])
	

def times_vec(v, a):	# ベクトル乗算
	return (v[0]*a, v[1]*a, v[2]*a)
	
	
def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	
	

def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)

#-------------------------------------------------------------------------------------------------------------


#----------  出力  ---------------------------------------------------------------------------------	

def write_surface(bS) :		#  自由曲面出力
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	
	for bL in bS :
		scene.begin_line(None, False)
		for p in bL :
			scene.append_point(p[0], p[1], p[2], p[3], p[4])
		scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	
	return ob 	#  作成した自由曲面を返す
		


	
def write_line(p, name) :	#  線形状出力	p : 制御点座標
	
	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
		
	
	
	
def write_tangent(p, v, vlen, name) :		#  接線ベクトル出力

	q = [0, 0, 0]
	for i in range(3) :
		q[i] = p[i] + vlen*v[i]

	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p, None, None, None, None)
	scene.append_point(q, None, None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す

#-------------------------------------------------------------------------------------------------------------



#----------  Bezier Patch  ---------------------------------------------------------------------------------

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])


def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)


def get_standard_length(bP) :	#  bP の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if bP[i][j] != None :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
	
	
	
def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])
	
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	
	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	
#-------------------------------------------------------------------------------------------------------------
	
	
	

#----------  matrix 演算  ------------------------------------------------------------------------------------
	
def matrix_cross_1(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[j][i][k]
	return C
	
	

def matrix_cross_1t(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B の転置マトリックス  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[i][j][k]	
	return C



def matrix_cross_2(A, B) :	#  4列 のベクトルトリックス A と、 4列 のマトリックス B の転置マトリックス  との積

	C = [0, 0, 0]
	
	for i in range(4) :
		for k in range(3) :
			C[k] = C[k] + A[i][k]*B[i]	
				
	return C
	


def matrix_cross_0(A, B, C) :	#  4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B と、 4列のマトリックス C の転置マトリックス  との積
	Vs = matrix_cross_1(A, B)
	return matrix_cross_2(Vs, C)	
#-------------------------------------------------------------------------------------------------------------




#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------

def get_T_matrix(t) :
	return ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	
	
def get_T_difference_matrix(t) :
	return (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	

def get_point_from_bezier_patch(Bp, t, s) :		#  指定する U, V パラメータ t, s における座標値を返す
	Tu = get_T_matrix(t)
	Tv = get_T_matrix(s)
	return matrix_cross_0(Tv, Bp, Tu)

	
def get_U_line_from_bezier_patch(Bp, s) :	#  Bezier Patch bP の Vパラメータ s における U線形状の制御点座標を返す
	Tv = get_T_matrix(s)				
	return matrix_cross_1(Tv, Bp)
	

def get_V_line_from_bezier_patch(Bp, t) :	#  Bezier Patch bP の Uパラメータ t における V線形状の制御点座標を返す
	Tu = get_T_matrix(t)				
	return matrix_cross_1t(Tu, Bp)
			

def get_std_L(p) :				#  Bezier line  の bounding box 基準の 基準長を返す
	a = [p[0][0], p[0][1], p[0][2]]
	b = [p[0][0], p[0][1], p[0][2]]
	for i in range(1, 4) :
		for k in range(3) :
			a[k] = max(a[k], p[i][k])
			b[k] = min(b[k], p[i][k])
	return a[0] + a[1] + a[2] - b[0] - b[1] - b[2]
	
	

def check_tangent(v, p, t) :			#  v = (0, 0, 0) となって計算上接線が求められない場合の処理				
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = sub_vec(p[2], p[0])
		else :
			v = sub_vec(p[3], p[1])
		v= normalize_vec(v)
			
	if v == (0, 0, 0) :
		v = sub_vec(p[3], p[0])
		v= normalize_vec(v)		
	return v	

	

def get_tangent(p, t) :	
	std_L = get_std_L(p)			#  線形状 p の基準長
	if std_L == 0 :
		return (0, 0, 0)			#  接線ベクトルが求まらない
	else :
		T_ = get_T_difference_matrix(t)
		v = matrix_cross_2(p, T_)
	
		#---------------------------------------------------
		print '     abs(v)/std_L = ', abs_vec(v)/std_L		#  report 
		#---------------------------------------------------
	
		if abs_vec(v)/std_L>= 1.5e-9 :
			v = normalize_vec(v)
			#---------------------------------------------------
			print '     接線ベクトル = ', v			#  report
			#---------------------------------------------------
		
		else :
			#--------------------------------------------------- 
			print '     補正前の接線ベクトル = ', normalize_vec(v)	#  report
			#---------------------------------------------------
		
			tt = round(t)			#  parameter t を 0 または 1 として再計算 ( ここに入ってくる場合 t は 0 または 1 に非常に近い )		
			if tt != t :
				T_ = get_T_difference_matrix(tt)
				v = matrix_cross_2(p, T_)
			v = normalize_vec(v)
			v = check_tangent(v, p, tt)
			#---------------------------------------------------
			print '     補正後の接線ベクトル = ', v		#  report
			#---------------------------------------------------
		
		return v


def get_U_tangent(Bp, t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における U方向接線単位ベクトルを返す
	Pu = get_U_line_from_bezier_patch(Bp, s)	#  Bezier Patch Bp の Vパラメータ s における U線形状の制御点座標
	v = get_tangent(Pu, t)
	return v

	
def get_V_tangent(Bp, t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における V方向接線単位ベクトルを返す
	Pv = get_V_line_from_bezier_patch(Bp, t)	#  Bezier Patch Bp の Uパラメータ t における V線形状の制御点座標
	v = get_tangent(Pv, s)	
	return v
#-------------------------------------------------------------------------------------------------------------	



	
	
def main(bS, paramL1, paramL2, vlen) :

	scene.inhibit_update()
	
	scene.create_part()			#  全体を格納するパート
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob = write_surface(bS)			#  サンプル自由曲面を出力	
	Bp = get_bezier_patch(ob)		#  サンプル自由曲面の Bezier Patch を求める
	
	#---------------------------------------------------
	print ''
	print ''
	print 'U方向接線'			# report
	#---------------------------------------------------
	scene.create_part('U方向接線')
	ob1 = scene.active_shape()
	ob1.disclosed = False
	
	for (t, s) in paramL1 :
		scene.create_part('t = ' + str(t))
		ob2 = scene.active_shape()
					
		p = get_point_from_bezier_patch(Bp, t, s)	#  指定する U, V パラメータ値における座標を取得
		
		Pu = get_U_line_from_bezier_patch(Bp, s)	#  指定する Vパラメータ値における U線形状の制御点座標
		write_line(Pu, 's = ' +str(s))
		
		#---------------------------------------------------
		print 't = ', t			# report
		#---------------------------------------------------
			
		Vu = get_U_tangent(Bp, t, s)			#  U, Vパラメータ値 における U接線ベクトルを取得		
		ob3 = write_tangent(p, Vu, vlen, 'U')		#  接線ベクトルを出力	

		ob2.activate()
		
		
	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	
	
	
scene = xshade.scene()
main(bS, paramL1, paramL2, vlen)
投稿数: 189
投稿日時: 2011-11-07 13:48
Re: 線形状(ベジェ曲線)の分割について教えてください
< スクリプト訂正 >

一つ前の書き込みで書いておいた script の中で、Bezier Patch を求める部分で間違って has_handle を使用したものを書いてしまいました。

失礼しました。

改めて書き直しておきます。


#  サンプル自由曲面の 指定する U パラメータ値 t,  V パラメータ値&#160;s における U, V 曲線を出力し、
#  U, V 方向の接線ベクトルを出力する



####    入力     #################################################
#  bS :  サンプル自由曲面

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (300, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (900, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -200), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 400), None)

bS = ((p1, p2), (p3, p4))


#  param :  U,V パラメータ値 (t, s)
param = (0.3, 0.7)


# vlen : 接線ベクトルの長さ
vlen = 500
############################################################



#----------  vector  演算  -----------------------------------------------------------------------------

def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	
	

def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)
		
		
def cross_vec(v1, v2) :		#  ベクトルの外積
	return[v1[1]*v2[2] - v1[2]*v2[1], v1[2]*v2[0] - v1[0]*v2[2], v1[0]*v2[1] - v1[1]*v2[0]]
	


#-------------------------------------------------------------------------------------------------------------


#----------  出力  ---------------------------------------------------------------------------------	

def write_surface(bS) :		#  自由曲面出力
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	
	for bL in bS :
		scene.begin_line(None, False)
		for p in bL :
			scene.append_point(p[0], p[1], p[2], p[3], p[4])
		scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	
	
	return ob 	#  作成した自由曲面を返す
		


	
def write_line(p, name) :	#  線形状出力	p : 制御点座標
	
	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
		

	
def write_tangent(p, v, vlen, name) :		#  接線ベクトル出力

	q = [0, 0, 0]
	for i in range(3) :
		q[i] = p[i] + vlen*v[i]

	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p, None, None, None, None)
	scene.append_point(q, None, None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
#-------------------------------------------------------------------------------------------------------------



#----------  Bezier Patch  ---------------------------------------------------------------------------------

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])


def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)


def get_standard_length(bP) :	#  bP の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if not ((i == 1 and j == 1) or (i == 1 and j == 2) or (i == 2 and j == 1) or (i == 2 and j == 2)) :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
	
	
	
def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP	
#-------------------------------------------------------------------------------------------------------------



#----------  matrix 演算  ------------------------------------------------------------------------------------
	
def matrix_cross_1(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[j][i][k]
	return C
	
	

def matrix_cross_1t(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B の転置マトリックス  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[i][j][k]	
	return C



def matrix_cross_2(A, B) :	#  4列 のベクトルトリックス A と、 4列 のマトリックス B の転置マトリックス  との積

	C = [0, 0, 0]
	
	for i in range(4) :
		for k in range(3) :
			C[k] = C[k] + A[i][k]*B[i]	
				
	return C
	


def matrix_cross_0(A, B, C) :	#  4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B と、 4列のマトリックス C の転置マトリックス  との積
	Vs = matrix_cross_1(A, B)
	return matrix_cross_2(Vs, C)	
#-------------------------------------------------------------------------------------------------------------




#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------

def get_T_matrix(t) :
	return ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	
	
def get_T_difference_matrix(t) :
	return (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	

def get_point_from_bezier_patch(Bp, t, s) :		#  指定する U, V パラメータ t, s における座標値を返す
	Tu = get_T_matrix(t)
	Tv = get_T_matrix(s)
	return matrix_cross_0(Tv, Bp, Tu)
	
	
def get_U_line_from_bezier_patch(Bp, s) :	#  Bezier Patch bP の Vパラメータ s における U線形状の制御点座標を返す
	Tv = get_T_matrix(s)				
	return matrix_cross_1(Tv, Bp)
		
	
def get_V_line_from_bezier_patch(Bp, t) :	#  Bezier Patch bP の Uパラメータ t における V線形状の制御点座標を返す
	Tu = get_T_matrix(t)				
	return matrix_cross_1t(Tu, Bp)
			

def check_tangent(v, p, t) :			#  v = (0, 0, 0) となって計算上接線が求められない場合の処理
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = sub_vec(p[2], p[0])
		else :
			v = sub_vec(p[3], p[1])
		v= normalize_vec(v)
			
	if v == (0, 0, 0) :
		v = sub_vec(p[3], p[0])
		v= normalize_vec(v)		
	return v	


def get_tangent(p, t) :
	T_ = get_T_difference_matrix(t)
	v = matrix_cross_2(p, T_)
	v = normalize_vec(v)
	return check_tangent(v, p, t)
	
	
def get_U_tangent(Bp, t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における U方向接線単位ベクトルを返す
	Pu = get_U_line_from_bezier_patch(Bp, s)	#  Bezier Patch Bp の Vパラメータ s における U線形状の制御点座標
	return get_tangent(Pu, t)
	
	
def get_V_tangent(Bp, t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における V方向接線単位ベクトルを返す
	Pv = get_V_line_from_bezier_patch(Bp, t)	#  Bezier Patch Bp の Uパラメータ t における V線形状の制御点座標
	return get_tangent(Pv, s)
	
	
def get_normal(Bp,t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における法線単位ベクトルを返す
	Vu = get_U_tangent(Bp, t, s)			#  Bezier Patch Bp の U, Vパラメータ t, s における U方向接線単位ベクトル
	Vv = get_V_tangent(Bp, t, s)			#  Bezier Patch Bp の U, Vパラメータ t, s における V方向接線単位ベクトル
	v = cross_vec(Vu, Vv)
	return normalize_vec(v)

#-------------------------------------------------------------------------------------------------------------	


	
	
def main(bS, (t, s), vlen) :

	scene.inhibit_update()
	
	scene.create_part()			#  全体を格納するパート
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob = write_surface(bS)			#  サンプル自由曲面を出力	
	Bp = get_bezier_patch(ob)		#  サンプル自由曲面の Bezier Patch を求める
	
	p = get_point_from_bezier_patch(Bp, t, s)	#  指定する U, V パラメータ値における座標を取得
	
	Pu = get_U_line_from_bezier_patch(Bp, s)	#  指定する Vパラメータ値における U線形状の制御点座標を取得
	ob1 = write_line(Pu, "s = " + str(s))		#  U 線形状を出力
	
	Pv = get_V_line_from_bezier_patch(Bp, t,)	#  指定する Uパラメータ値における V線形状の制御点座標を取得
	ob2 = write_line(Pv, "t = " + str(t))		#  V 線形状を出力
	
	Vu = get_U_tangent(Bp, t, s)			#  U, Vパラメータ t, s における U方向接線単位ベクトルを取得
	ob3 = write_tangent(p, Vu, vlen, 'u - tangent')	#  接線ベクトルを出力
	
	Vv = get_V_tangent(Bp, t, s)			#  U, Vパラメータ t, s における V方向接線単位ベクトルを取得		
	ob4 = write_tangent(p, Vv, vlen, 'v - tangent')	#  接線ベクトルを出力
	
	Vn= get_normal(Bp, t, s)			#  U, Vパラメータ t, s における法線単位ベクトルを取得		
	ob5 = write_tangent(p, Vn, vlen, 'normal')	#  法線ベクトルを出力
	
	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main(bS, param, vlen)
投稿数: 189
投稿日時: 2011-11-06 18:25
Re: 線形状(ベジェ曲線)の分割について教えてください
    
< Bezier 曲面の接線, 法線 >

Bezier 曲線の指定する U, V パラメータ t, s における ポイント座標や接線ベクトルを求める場合、マトリックス計算を用いるのが便利です。



[ マトリックスの積 ]

1 行, 4列 のマトリックス A と、4行, 4列のベクトルマトリックス B の積から 1行, 4列のベクトルマトリックス C が得られる。

   
def matrix_cross_1(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[j][i][k]
	return C
	
	

def matrix_cross_1t(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B の転置マトリックス  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[i][j][k]	
	return C



1 行, 4列 のベクトルマトリックス A と、4行, 1列のマトリックス B の積から ベクトル C が得られる。

   
def matrix_cross_2(A, B) :	#  4列 のベクトルトリックス A と、 4列 のマトリックス B の転置マトリックス  との積

	C = [0, 0, 0]
	
	for i in range(4) :
		for k in range(3) :
			C[k] = C[k] + A[i][k]*B[i]	
				
	return C



1 行, 4列 のマトリックス A と、4行, 4列のベクトルマトリックス B と、4行, 1列のマトリックス C の積から ベクトル D が得られる。

   
def matrix_cross_0(A, B, C) :	#  4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B と、 4列のマトリックス C の転置マトリックス  との積
	Vs = matrix_cross_1(A, B)
	return matrix_cross_2(Vs, C)






Matrix 計算を用いれば、Bezier 関連の諸量は次のように表されます。



[ Bezier 曲線上のポイント座標 ]

制御点座標 [ p0, p1, p2, p3 ] で与えられる Bezier 曲線上のパラメータ t でのポイント座標 Q

  matrix_cross_2( ) で計算

   



[ Bezier 曲線の接線ベクトル ]

制御点座標 [ p0, p1, p2, p3 ] で与えられる Bezier 曲線上のパラメータ t での接線ベクトル V

  matrix_cross_2( ) で計算 

   



[ Bezier 曲面の U 方向線形状 ]

Bezier Patch [ Bp ] で与えられる Bezier 曲面の V パラメータ s で指定される U 方向 Bezier line の制御点座標 Pu = [ p0, p1, p2, p3 ]

  matrix_cross_1( ) で計算

   



[ Bezier 曲面の V 方向線形状 ]

Bezier Patch [ Bp ] で与えられる Bezier 曲面の U パラメータ t で指定される V 方向 Bezier line の制御点座標 Pv = [ p0, p1, p2, p3 ]

  matrix_cross_1_t( ) で計算

   



[ Bezier 曲面上の座標値 ]

Bezier Patch [ Bp ] で与えられる Bezier 曲面の U, V パラメータ t, s で指定されるポイントの座標値 Q

  matrix_cross_0( ) で計算

   



[ Bezier 曲面の U 方向接線ベクトル ]

Bezier Patch [ Bp ] で与えられる Bezier 曲面の U,V パラメータ t, s で指定される U 方向接線ベクトル Vu

  matrix_cross_0( ) で計算

   



[ Bezier 曲面の V 方向接線ベクトル ]

Bezier Patch [ Bp ] で与えられる Bezier 曲面の U,V パラメータ t, s で指定される V 方向接線ベクトル Vv

  matrix_cross_0( ) で計算

   



[ Bezier 曲面の 面法線ベクトル ]

Bezier Patch [ Bp ] で与えられる Bezier 曲面の U,V パラメータ t, s で指定される面法線 Vn
  
  面法線 = U方向接線ベクトル と V方向接線ベクトルの外積

   



以下の script ではサンプル自由曲面の U, V パラメータ値における U, V 曲線と接線, 面法線 を出力します。

接線ベクトルは制御点座標の相対位置の差で決まるため、数値計算の丸め誤差の影響が大きく現れてしまうケースがあります。
また、球の極における V方向( 経線方向 )の接線が存在しないので、極における面法線が求められないなどの興味深い話がありますが、それらに対する対処については次回で調べます。


#  サンプル自由曲面の 指定する U パラメータ値 t,  V パラメータ値&#160;s における U, V 曲線を出力し、
#  U, V 方向の接線ベクトルを出力する



####    入力     #################################################
#  bS :  サンプル自由曲面

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (300, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (900, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -200), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 400), None)

bS = ((p1, p2), (p3, p4))


#  param :  U,V パラメータ値 (t, s)
param = (0.3, 0.7)


# vlen : 接線ベクトルの長さ
vlen = 500
############################################################



#----------  vector  演算  -----------------------------------------------------------------------------

def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	
	

def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)
		
		
def cross_vec(v1, v2) :		#  ベクトルの外積
	return[v1[1]*v2[2] - v1[2]*v2[1], v1[2]*v2[0] - v1[0]*v2[2], v1[0]*v2[1] - v1[1]*v2[0]]
	


#-------------------------------------------------------------------------------------------------------------


#----------  出力  ---------------------------------------------------------------------------------	

def write_surface(bS) :		#  自由曲面出力
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	
	for bL in bS :
		scene.begin_line(None, False)
		for p in bL :
			scene.append_point(p[0], p[1], p[2], p[3], p[4])
		scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	
	
	return ob 	#  作成した自由曲面を返す
		


	
def write_line(p, name) :	#  線形状出力	p : 制御点座標
	
	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
		

	
def write_tangent(p, v, vlen, name) :		#  接線ベクトル出力

	q = [0, 0, 0]
	for i in range(3) :
		q[i] = p[i] + vlen*v[i]

	scene.begin_creating()
	scene.begin_line(None, False)
	
	scene.append_point(p, None, None, None, None)
	scene.append_point(q, None, None, None, None)
		
	scene.end_line()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.name = name
	
	return ob 	#  作成した線形状を返す
#-------------------------------------------------------------------------------------------------------------



#----------  Bezier Patch  ---------------------------------------------------------------------------------

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])



def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(ob2.control_point(0).has_out_handle, ob2.control_point(1).has_in_handle),
		(ob2.control_point(2).has_out_handle, ob2.control_point(3).has_in_handle))
	has_vH = (
		(ob2.control_point(0).has_lateral_out_handle, ob2.control_point(1).has_lateral_out_handle),
		(ob2.control_point(2).has_lateral_in_handle, ob2.control_point(3).has_lateral_in_handle))
	
	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP	
#-------------------------------------------------------------------------------------------------------------



#----------  matrix 演算  ------------------------------------------------------------------------------------
	
def matrix_cross_1(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]

	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[j][i][k]
	return C
	
	

def matrix_cross_1t(A, B) :	#  1行, 4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B の転置マトリックス  との積

	C = [[0, 0, 0], [0, 0, 0], [0, 0, 0], [0, 0, 0]]
	
	for i in range(4) :
		for j in range(4) :
			for k in range(3) :
				C[i][k] = C[i][k] + A[j]*B[i][j][k]	
	return C



def matrix_cross_2(A, B) :	#  4列 のベクトルトリックス A と、 4列 のマトリックス B の転置マトリックス  との積

	C = [0, 0, 0]
	
	for i in range(4) :
		for k in range(3) :
			C[k] = C[k] + A[i][k]*B[i]	
				
	return C
	


def matrix_cross_0(A, B, C) :	#  4列 のマトリックス A と、 4行, 4列 のベクトルマトリックス B と、 4列のマトリックス C の転置マトリックス  との積
	Vs = matrix_cross_1(A, B)
	return matrix_cross_2(Vs, C)	
#-------------------------------------------------------------------------------------------------------------




#----------  Bezier 関連の諸量  ----------------------------------------------------------------------------

def get_T_matrix(t) :
	return ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	
	
def get_T_difference_matrix(t) :
	return (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	

def get_point_from_bezier_patch(Bp, t, s) :		#  指定する U, V パラメータ t, s における座標値を返す
	Tu = get_T_matrix(t)
	Tv = get_T_matrix(s)
	return matrix_cross_0(Tv, Bp, Tu)
	
	
def get_U_line_from_bezier_patch(Bp, s) :	#  Bezier Patch bP の Vパラメータ s における U線形状の制御点座標を返す
	Tv = get_T_matrix(s)				
	return matrix_cross_1(Tv, Bp)
		
	
def get_V_line_from_bezier_patch(Bp, t) :	#  Bezier Patch bP の Uパラメータ t における V線形状の制御点座標を返す
	Tu = get_T_matrix(t)				
	return matrix_cross_1t(Tu, Bp)
			

def check_tangent(v, p, t) :			#  v = (0, 0, 0) となって計算上接線が求められない場合の処理
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = sub_vec(p[2], p[0])
		else :
			v = sub_vec(p[3], p[1])
		v= normalize_vec(v)
			
	if v == (0, 0, 0) :
		v = sub_vec(p[3], p[0])
		v= normalize_vec(v)		
	return v	


def get_tangent(p, t) :
	T_ = get_T_difference_matrix(t)
	v = matrix_cross_2(p, T_)
	v = normalize_vec(v)
	return check_tangent(v, p, t)
	
	
def get_U_tangent(Bp, t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における U方向接線単位ベクトルを返す
	Pu = get_U_line_from_bezier_patch(Bp, s)	#  Bezier Patch Bp の Vパラメータ s における U線形状の制御点座標
	return get_tangent(Pu, t)
	
	
def get_V_tangent(Bp, t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における V方向接線単位ベクトルを返す
	Pv = get_V_line_from_bezier_patch(Bp, t)	#  Bezier Patch Bp の Uパラメータ t における V線形状の制御点座標
	return get_tangent(Pv, s)
	
	
def get_normal(Bp,t, s) :				#  Bezier Patch Bp の U, Vパラメータ t, s における法線単位ベクトルを返す
	Vu = get_U_tangent(Bp, t, s)			#  Bezier Patch Bp の U, Vパラメータ t, s における U方向接線単位ベクトル
	Vv = get_V_tangent(Bp, t, s)			#  Bezier Patch Bp の U, Vパラメータ t, s における V方向接線単位ベクトル
	v = cross_vec(Vu, Vv)
	return normalize_vec(v)

#-------------------------------------------------------------------------------------------------------------	


	
	
def main(bS, (t, s), vlen) :

	scene.inhibit_update()
	
	scene.create_part()			#  全体を格納するパート
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob = write_surface(bS)			#  サンプル自由曲面を出力	
	Bp = get_bezier_patch(ob)		#  サンプル自由曲面の Bezier Patch を求める
	
	p = get_point_from_bezier_patch(Bp, t, s)	#  指定する U, V パラメータ値における座標を取得
	
	Pu = get_U_line_from_bezier_patch(Bp, s)	#  指定する Vパラメータ値における U線形状の制御点座標を取得
	ob1 = write_line(Pu, "s = " + str(s))		#  U 線形状を出力
	
	Pv = get_V_line_from_bezier_patch(Bp, t,)	#  指定する Uパラメータ値における V線形状の制御点座標を取得
	ob2 = write_line(Pv, "t = " + str(t))		#  V 線形状を出力
	
	Vu = get_U_tangent(Bp, t, s)			#  U, Vパラメータ t, s における U方向接線単位ベクトルを取得
	ob3 = write_tangent(p, Vu, vlen, 'u - tangent')	#  接線ベクトルを出力
	
	Vv = get_V_tangent(Bp, t, s)			#  U, Vパラメータ t, s における V方向接線単位ベクトルを取得		
	ob4 = write_tangent(p, Vv, vlen, 'v - tangent')	#  接線ベクトルを出力
	
	Vn= get_normal(Bp, t, s)			#  U, Vパラメータ t, s における法線単位ベクトルを取得		
	ob5 = write_tangent(p, Vn, vlen, 'normal')	#  法線ベクトルを出力
	
	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main(bS, param, vlen)
投稿数: 189
投稿日時: 2011-10-29 23:59
Re: 線形状(ベジェ曲線)の分割について教えてください
 
< 制御点座標の階差表記の利点 >


Bezier 制御点座標の階差表記については、参考書でもネット上でも見かけたことがなく、いささかナニコレなのですが、この階差表記の利点をもう一つ示しておきます。


下記のスクリプトでは、制御点座標を通常表記と階差表記の両ケースで接線ベクトルを求めます。

サンプル Bezier 曲線は outhandle のない曲線であり、これに対して微小のパラメータ値における接線ベクトルを求めています。




下に示したとおり、両者の違いは明らかです。

通常表記では微小のパラメータ値に対しては計算中に発生する誤差によって、接線の方向があらぬ向きになってしまいますが、階差表記ではこれをほぼ回避できています。





レポート結果から次のことが解ります。

  1)単位ベクトル化する前の P'(t) = dF/dt の長さが微小になる場合、
    通常表記では丸め誤差によって接線の方向があらぬ向きになってしまうことがある。

  2)階差表記ではこの誤差を非常に小さく抑えておくことができ、格段の注意を払う必要がない。

  3)abs(P'(t)) < 1e-6 であれば、正しい接線の方向は t = 0 として求めた接線の方向と極めて近く、
    実用上問題の無い近似解として用いることができる。( 通常表記ではこれを利用すればよい )


#  サンプル線形状の指定する位置の接線方向を Bezier 制御点座標の 通常表記と階差表記の両ケースで求める。

#####  入力  #####
v_len = 2000			#  接線を表示する線分の長さ

p = []				#  サンプル線形状の Bezier 制御点座標
p.append((0, 1000, 0))
p.append((0, 1000, 0))
p.append((1000, 600, 0))
p.append((1000, 0, 0))

paramL = []			#  接線を求める位置を指定する Bezier パラメータ
paramL .append(1e-8)
paramL .append(1e-9)	
paramL .append(1e-10)	
paramL .append(1e-11)	
paramL .append(1e-14)	
paramL .append(1e-15)	
paramL .append(1e-16)	
paramL .append(0)			
#################


def add_vec(u, v):	#  ベクトル加算
	return (u[0] + v[0], u[1] + v[1], u[2] + v[2])
	
	
def sub_vec(u, v):	#  ベクトル減算
	return (u[0] - v[0], u[1] - v[1], u[2] - v[2])
	

def times_vec(v, a):	# ベクトル乗算
	return (v[0]*a, v[1]*a, v[2]*a)
	
	
def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	


def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)
		


def write_line(p) :	 #  サンプル線形状を描く
			
	scene.begin_creating()
	scene.begin_line(None, False)
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
	scene.end_line()	
	scene.end_creating()
	
	ob = scene.active_shape()
	
	return ob

	
	
def differences_expression(p) :		#  Bezier 制御点座標を階差表示形式に変換
	d11 = sub_vec(p[1],  p[0])
	d12 = sub_vec(p[2],  p[1])
	d13 = sub_vec(p[3],  p[2])
	d21 = sub_vec(d12, d11)
	d22 = sub_vec(d13, d12)
	d31 = sub_vec(d22, d21)

	return (p[0], d11, d21, d31)
	


def matrix_calc_2(m, mv) :	#  1行4列のマトリックスと、4行1列のベクトルマトリックスの積
	v = (0, 0, 0)
	for i in range(4) :
		v = add_vec(v, times_vec(mv[i], m[i]))
	return v
	
	
	
def get_point_1(p, t) :		#  パラメータ t の位置の座標値を返す ( p は通常の座標値 )
	Tu = ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	return matrix_calc_2(Tu, p)
	
	
def get_point_2(p, t) :		#  パラメータ t の位置の座標値を返す ( p は制御点座標の階差表示 )
	Tu = (1, 3*t, 3*t**2, t**3)
	return matrix_calc_2(Tu, p)
	


def get_tangent_1(p, t) :		#  パラメータ t の位置の接線単位ベクトルを返す ( p は通常の座標値 )
	Tu = (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	v = matrix_calc_2(Tu, p)
	
	print 't = ', t, '   ベクトル長さ = ', abs_vec(v)	#  通常表記で計算されたベクトル長さをレポート
	
	v = normalize_vec(v)
	
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = sub_vec(p[2], p[0])
		else :
			v = sub_vec(p[3], p[1])
		v= normalize_vec(v)
	
	if v == (0, 0, 0) :
		v = sub_vec(p[3], p[0])
		v= normalize_vec(v)
		
	return v
	

	
def get_tangent_2(p, t) :		#  パラメータ t の位置の接線単位ベクトルを返す ( p は制御点座標の階差表示 )
	Tu = (0, 1, 2*t, t**2)
	v = matrix_calc_2(Tu, p)
	v = normalize_vec(v)
	
	if v == (0, 0, 0) :
		if t < 0.5 :
			v = p[2]
		else :
			v = add_vec(p[2], p[3])
			v = times_vec(v, -1)
		v= normalize_vec(v)
	
	if v == (0, 0, 0) :
		v = p[3]
		v= normalize_vec(v)
		
	return v

	

def write_tangent(p0, p1, name) :
	scene.begin_creating()
	scene.begin_line(name, False)
	scene.append_point(p0, None, None, None, None)
	scene.append_point(p1, None, None, None, None)
	scene.end_line()
	scene.end_creating()
	


	
def main(p, paramL, v_len) :
	scene.inhibit_update()

	scene.create_part()
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	print ''
	print '通常表記における単位ベクトル化する前のベクトル長さ'
	
	ob = write_line(p)		#  サンプル線形状を描く
	
	for t in paramL :
		scene.create_part('t = ' + str(t))
		ob1 = scene.active_shape()
		
		p0 = get_point_1(p, t)		#  パラメータ t の位置の座標値
		v = get_tangent_1(p, t)	#  パラメータ t の位置の接線単位ベクトル
		v = times_vec(v, v_len)
		p1 = add_vec(p0, v)
		write_tangent(p0, p1, '通常表記')
		
		pp = differences_expression(p)	#  階差表示
		p0 = get_point_2(pp, t)	#  パラメータ t の位置の座標値
		v = get_tangent_2(pp, t)	#  パラメータ t の位置の接線単位ベクトル
		v = times_vec(v, v_len)
		p1 = add_vec(p0, v)
		write_tangent(p0, p1, '階差表記')
	
		ob1.activate()
		
	ob0.activate()
	ob0.reset_transformation()
	
	scene.allow_update()
	scene.force_update()
		
	
scene = xshade.scene()
main(p, paramL, v_len)
投稿数: 189
投稿日時: 2011-10-21 01:12
Re: 線形状(ベジェ曲線)の分割について教えてください

< 接線ベクトル >


曲線における接線ベクトルとは、幾何数学では暗黙の了解事項として、

    弧長( 曲線の長さ )をパラメータとした弧長関数 F(s)で表された曲線の
    弧長 s に関する導関数 dF/ds で表される単位ベクトル

のことをいいます。( dF/ds はその定義から必ず単位ベクトルになります )



「 弧長パラメータ 」などと聞き慣れない言葉が出てきますが、曲線の長さをパラメータとしたもので、下図中央のような性質を持ちます。

  


Bezier 曲線を弧長パラメータで表すこともできますが、解析的な数式では表されないので実用的ではありません。


曲線を F、弧長を s、Bezier パラメータを t とすれば、

    接線ベクトル = dF/ds

s を t に変数変換すると

    dF/ds = (dF/dt) × (dt/ds)


ここで dF/dt が先に接線の向きを求めるときに使用した P '(t) と同じことを意味しますから、

    P '(t) = dF/dt に
    dt/ds を乗じたものが
    接線ベクトル dF/ds に等しい

ちなみに dF/dt, dF/ds はベクトル関数、dt/ds はスカラー関数です。


先に outhandle 長さが 0 の時に、P '(t) = (0, 0, 0) となって接線の方向が求まらないと書きましたが、これを弧長 s を使って説明すると、

    P '(t) = dF/dt
       = (dF/ds) / (dt/ds)
       = (dF/ds) * (ds/dt)
       = ( 接線ベクトル ) x (ds/dt) において ds/dt = 0 となるため、
       = (0, 0, 0)

接線ベクトルは存在するが、それに乗ぜられる ds/dt が 0 となるため P '(t) = (0, 0, 0) となってしまうという意味になります。





3D の専門書の多くでは、解析的な Bezier パラメータ t を用いた P '(t) = dF/dt にしか興味は無く、弧長関数などには無関心ですので、「 接線ベクトル 」という言葉を dF/ds ではなく、dF/dt の意味で使用しているケースが見られます。

「 正しい、間違っている 」という話ではなく、「 分野が違えば使用する言葉の意味も違ってくる 」ということで、たまに混用されて誤解の素となっているみたいです。





ついでに、弧長 s は何も怪しげなモノではないという証拠を一つ記しておきます。

    P '(t) = dF/dt に
    dt/ds を乗じたものが
    接線ベクトル dF/ds に等しい

    接線ベクトルは単位ベクトル

ということから、dt/ds は ベクトル P '(t) の長さの逆数になっています。

ですから、ベクトル P '(t) の長さを次のように表せます。

    abs(P '(t)) = 1/(dt/ds)  ( dt/ds > 0 )
         = ds/dt

これを t で積分すると、弧長 s つまり曲線の長さが得られることになります。

    ∫ abs(P '(t)) dt = ∫ (ds/dt) dt
           = ∫ ds
           = s

これはこのスレッドの最初の方で Bezier 曲線の長さを数値積分で求めるときに用いた式にほかなりません。

    s = ∫ abs(P '(t)) dt
投稿数: 189
投稿日時: 2011-10-21 01:11
Re: 線形状(ベジェ曲線)の分割について教えてください
 
< Bezier 曲線の接線 >


Bezier 曲線の4つの制御点座標 P0, P1, P2, P3 から、パラメータ t における座標値と導関数は次のようになります


  パラメータ t における座標 P(t) は、

    P(t) = ( (1 - t)^3*P0 + 3*t*(1 - t)^2*P1 + 3*t^2*(1 - t)*P2 + t^3*P3 )

  パラメータ t における導関数 P '(t) は、

    P '(t) = ( -3*(1 - t)^2*P0 + 3*(1 - t)*(1 - 3*t)*P1 + 3*t*(2 - 3*t)*P2 + 3*t^2*P3 )


P0〜P3 の代わりに階差表示 d0〜d3 を用いれば、

    d0 = P0
    d1 = P1 - P0
    d2 = (P2 - P1) - (P1 - P0)
    d3 = (P3 - P2) - (P2 - P1) - (P2 - P1) + (P1 - P0)

  パラメータ t における座標 P(t) は、

    P(t) = ( d0 + 3*t*d1 + 3*t^2*d2 + t^3*d3 )

  パラメータ t における導関数 P '(t) は、

    P '(t) = 3*(d1 + 2*t*d2 + t^2*d3)



この導関数 P '(t) が接線の方向を表します。


ここで階差表示の式から t = 0 つまり anchor point P0 における接線の方向を見ると、

    P '(0) = 3*d1
       = 3*(P1 - P0)

となり、outhandle の方向が接線の方向になるという、よく知られた結果通りになります。

  


ところで outhandle の長さが 0 の場合には d1 = (0, 0, 0) なので、接線の方向は上の式からは P '(0) = (0, 0, 0) となってしまいます。

これは「 接線の方向が求められない 」のであって、「 接線が定義できない 」ことではありません。

何故求められないのかについては後で話すことにして、このケースでは次のように考えてみます。

    P '(t) = 3*(d1 + 2*t*d2 + t^2*d3) において、d1 = (0, 0, 0) ならば、

    P '(t) = 3*(2*t*d2 + t^2*d3) となり、

    t を限りなく 0 に近づけると、t^2 は t に比べて急速に 0 に近づくので無視することができ、

    t → 0 では P '(t) = 6*t*d2 となり、

    接線の方向は
      d2 = (P2 - P1) - (P1 - P0)
        = P2 - P1
        = P2 - P0
    で与えられる


つまり、outhandle の長さが 0 の場合には P2 - P0 で P0 における接線の方向が求められます。

  


さらに outhandle 長さ 0 に加えて、P3 の inhandle P2 が P0 と一致してしまう場合には、d2 = (0, 0, 0) となるので、

    P '(t) = 3*t^2*d3 となり、

    接線の方向は
      d3 = (P3 - P2) - (P2 - P1) - (P2 - P1) + (P1 - P0)
        = (P3 - P2) - (P2 - P1)
        = P3 - P2
        = P3 - P0
    で与えられる

これは言うまでもないことです。( このような説明に、階差表示はとても解りやすくで便利です )

  


下記のスクリプトでは、サンプル線形状の指定する位置での接線の方向を表示します。

#  サンプル線形状の指定する位置の接線方向を示す
#  サンプル線形状の指定する位置の接線方向を示す

#####  入力  #####
v_len = 2000			#  接線を表示する線分の長さ
Ls = []				#  サンプル線形状の Bezier 制御点座標

L = []
p = []
p.append((0, 0, 0))
p.append((1000, 1000, 0))
p.append((2000, 1000, 0))
p.append((3000, 0, 0))
L.append(p)
L.append(0.3)			#  接線を求める位置を指定する Bezier パラメータ
Ls.append(L)

L = []
p = []
p.append((0, 0, 0))
p.append((1000, 1000, 0))
p.append((2000, 1000, 0))
p.append((3000, 0, 0))
L.append(p)
L.append(0)			#  接線を求める位置を指定する Bezier パラメータ
Ls.append(L)

L = []
p = []
p.append((0, 0, 0))
p.append((0, 0, 0))
p.append((2000, 1000, 0))
p.append((3000, 0, 0))
L.append(p)
L.append(0)			#  接線を求める位置を指定する Bezier パラメータ
Ls.append(L)

L = []
p = []
p.append((0, 0, 0))
p.append((0, 0, 0))
p.append((0, 0, 0))
p.append((3000, 0, 0))
L.append(p)
L.append(0)			#  接線を求める位置を指定する Bezier パラメータ
Ls.append(L)

L = []
p = []
p.append((0, 0, 0))
p.append((1000, 1000, 0))
p.append((3000, 0, 0))
p.append((3000, 0, 0))
L.append(p)
L.append(1)			#  接線を求める位置を指定する Bezier パラメータ
Ls.append(L)

L = []
p = []
p.append((0, 0, 0))
p.append((3000, 0, 0))
p.append((3000, 0, 0))
p.append((3000, 0, 0))
L.append(p)
L.append(1)			#  接線を求める位置を指定する Bezier パラメータ
Ls.append(L)
#################


def add_vec(u, v):	#  ベクトル加算
	return (u[0] + v[0], u[1] + v[1], u[2] + v[2])
	
	
def sub_vec(u, v):	#  ベクトル減算
	return (u[0] - v[0], u[1] - v[1], u[2] - v[2])
	

def times_vec(v, a):	# ベクトル乗算
	return (v[0]*a, v[1]*a, v[2]*a)
	
	
def abs2_vec(v) :		#  ベクトル長さの二乗
	return v[0]**2 + v[1]**2 + v[2]**2
	
	
def abs_vec(v) :		#  ベクトル長さ
	import math
	a = abs2_vec(v)
	return math.sqrt(a)
	


def normalize_vec(v):	# 単位ベクトル化  長さが 0 ならば (0, 0, 0) を返す
	a = abs_vec(v)
	if a != 0 :
		return (v[0]/a, v[1]/a, v[2]/a)	
	else :
		return (0, 0, 0)
		


def write_line(p) :	 #  サンプル線形状を描く
		
	scene.create_part()
	ob = scene.active_shape()
	ob.disclosed = False
		
	scene.begin_creating()
	scene.begin_line(None, False)
	scene.append_point(p[0], None, p[1], None, None)
	scene.append_point(p[3], p[2], None, None, None)
	scene.end_line()	
	scene.end_creating()
	
	ob1 = ob.son.bro
	ob.activate()
	
	return ob1

	
	
def differences_expression(p) :		#  Bezier 制御点座標を階差表示形式に変換
	d11 = sub_vec(p[1],  p[0])
	d12 = sub_vec(p[2],  p[1])
	d13 = sub_vec(p[3],  p[2])
	d21 = sub_vec(d12, d11)
	d22 = sub_vec(d13, d12)
	d31 = sub_vec(d22, d21)

	return (p[0], d11, d21, d31)
	


def matrix_calc_2(m, mv) :	#  1行4列のマトリックスと、4行1列のベクトルマトリックスの積
	v = (0, 0, 0)
	for i in range(4) :
		v = add_vec(v, times_vec(mv[i], m[i]))
	return v
	
	
	
def get_point_1(p, t) :		#  パラメータ t の位置の座標値を返す ( p は通常の座標値 )
	Tu = ((1 - t)**3, 3*t*(1 - t)**2, 3*t**2*(1 - t), t**3)
	return matrix_calc_2(Tu, p)
	
	
def get_point_2(p, t) :		#  パラメータ t の位置の座標値を返す ( p は制御点座標の階差表示 )
	Tu = (1, 3*t, 3*t**2, t**3)
	return matrix_calc_2(Tu, p)
	


def get_tangent_1(p, t) :		#  パラメータ t の位置の接線単位ベクトルを返す ( p は通常の座標値 )
	Tu = (-3*(1 - t)**2, 3*(1 - t)*(1 - 3*t), 3*t*(2 - 3*t), 3*t**2)
	tangent_v = matrix_calc_2(Tu, p)
	tangent_v= normalize_vec(tangent_v)
	
	if tangent_v == (0, 0, 0) :
		if t < 0.5 :
			tangent_v = sub_vec(p[2], p[0])
		else :
			tangent_v = sub_vec(p[3], p[1])
		tangent_v= normalize_vec(tangent_v)
	
	if tangent_v == (0, 0, 0) :
		tangent_v = sub_vec(p[3], p[0])
		tangent_v= normalize_vec(tangent_v)
		
	return tangent_v
	

	
def get_tangent_2(p, t) :		#  パラメータ t の位置の接線単位ベクトルを返す ( p は制御点座標の階差表示 )
	Tu = (0, 1, 2*t, t**2)
	tangent_v = matrix_calc_2(Tu, p)
	tangent_v= normalize_vec(tangent_v)
	
	if tangent_v == (0, 0, 0) :
		if t < 0.5 :
			tangent_v = p[2]
		else :
			tangent_v = add_vec(p[2], p[3])
			tangent_v = times_vec(tangent_v, -1)
		tangent_v= normalize_vec(tangent_v)
	
	if tangent_v == (0, 0, 0) :
		tangent_v = p[3]
		tangent_v= normalize_vec(tangent_v)
		
	return tangent_v

	
	
def write_tangent(ob, p, t, v_len) : 	#  線形状 ob の Bezier パラメータ t の位置の接線を描く
	v0 = get_point_1(p, t)		#  パラメータ t の位置の座標値
	v = get_tangent_1(p, t)		#  パラメータ t の位置の接線単位ベクトル
	
#	p = differences_expression(p)	#  階差表示
#	v0 = get_point_2(p, t)		#  パラメータ t の位置の座標値
#	v = get_tangent_2(p, t)		#  パラメータ t の位置の接線単位ベクトル
	
	v1 = times_vec(v, v_len)
	v1 = add_vec(v1, v0)

	scene.begin_creating_at(ob)
	scene.begin_line('t = ' + str(t), False)
	scene.append_point(v0, None, None, None, None)
	scene.append_point(v1, None, None, None, None)
	scene.end_line()	
	scene.end_creating()


	
def main(Ls, v_len) :
	scene.inhibit_update()

	scene.create_part()
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	for L in Ls :
		p = L[0]
		t = L[1]
		ob = write_line(p)		#  サンプル線形状を描く
		write_tangent(ob, p, t, v_len)	#  サンプル線形状の Bezier パラメータ t の位置の接線を描く
		if (t > 0) and (t < 1) :
			ob.insert_control_point(t)	#  パラメータ t の位置にポイントを追加
		
	ob0.activate()
	ob0.reset_transformation()
	
	scene.allow_update()
	scene.force_update()
		
	
scene = xshade.scene()
main(Ls, v_len)
投稿数: 189
投稿日時: 2011-10-07 14:50
Re: 線形状(ベジェ曲線)の分割について教えてください
< 精度 >


Bezier Patch を求める際の 0 判定に bounding box size の 1/10000 が出てきました。

そして Bezier Line の handle linked 判定で r/L <= 0.00025 も出てきました。

これに対して in handle, out handle を囲む boding box を考えれば、r/L <= 0.00025 は、bounding box の 1/10000 のオーダーと捉えることもできます。

bounding box size の 1/10000 のオーダーとは、 1/10000 〜 1/1000 という幅を持った表現の仕方と思って下さい。

この 「 1/10000 のオーダー 」というは、かっては Shade の色々な所に見られた数値です。

ポリゴンメッシュにもありましたし、昔の仕様では数値入力による縮小では、最小スケールが 1/10000 に制限されたりもしていました。( これは undo の際に逆マトリックスを求める必要からリミッターを与えていたらしいです )



全くの推測なのですが、この 1/10000 という数字は Shade のオリジナルプログラマーである時枝さんによって与えられた Shade 全体に対する精度上の一つの目安とする数値ではないかと。

特に Bezier 周りはこのルールでかっきりと組み立てられていますから、今でも安心して使えます。


Shade の開発体制が強化されてからは、プログラマーの方が増えたせいもあるのでしょうか、この「 1/10000 」ルールもだんだんと崩れて行っている気がします。

特に、ポリゴンメッシュ関連、

以前から、「 あれれ‥‥? 」と気になっていた所があったのですが、今回改めて確認してみたら、ちょっと、これは、‥‥‥‥‥‥‥‥‥‥


ここから先は、何を調べたかだけを記しておきます。
( Shade 12.0.3 にて調査 )



1)非平面ポリゴン表示

下記のスクリプトでは、非平面度が少しずつ大きくなるポリゴンメッシュを作ります。

前半は上位パートにマトリックスがかかっていないケース
後半は掛かっているケースです。


スクリプトを実行した後、出力されたポリゴンメッシュを選択し、編集モードに入って全ての面を選択し、

  メニュー > メッシュ >エラー診断選択 > 非平面

を選んで非平面と判定されたポリゴンを確認し、非平面の判定がどのレベルで変わるかを調べてみます。


#  ポリゴンメッシュの非平面チェック基準を確認する

#  出力されたポリゴンメッシュを選択し、編集モードに入って
#	メニュー > メッシュ >エラー診断選択 > 非平面
#  を選んで確認する

#  結果 : 平面と見なされる平面からのずれの許容値は bounding box size の 1/10000 くらい
#            しかし非平面の判定に上位マトリックスが影響している。



scene = xshade.scene()

Ls = []

L = []
L.append((1, 0.0006))
L.append((1, 0.0005))
L.append((1, 0.0004))
L.append((1, 0.0003))
L.append((1, 0.0002))
L.append((1, 0.0001))
L.append((1, 0.00009))
L.append((1, 0.00008))
Ls.append(L)

L = []
L.append((10, 0.006))
L.append((10, 0.005))
L.append((10, 0.004))
L.append((10, 0.003))
L.append((10, 0.002))
L.append((10, 0.001))
L.append((10, 0.0009))
L.append((10, 0.0008))
Ls.append(L)




scene.inhibit_update()

scene.exit_modify_mode()
scene.select_all()
scene.active_shape().transformation_matrix = ((1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1))

scene.create_part()
ob1 = scene.active_shape()
#ob1.cancel_transformation()

for L in Ls :
	scene.create_part()
	for data in L :
		y = float(data[0])/2
		x = y*2
		epsilon = data[1]
		
		scene.begin_creating()
		scene.begin_polygon_mesh(str(data))
	
		scene.append_polygon_mesh_vertex((0, y , epsilon))
		scene.append_polygon_mesh_vertex((0, -y, 0))
		scene.append_polygon_mesh_vertex((x, -y, 0))
		scene.append_polygon_mesh_vertex((x, y, 0))
	
		scene.append_polygon_mesh_vertex((0, -y, 0))
		scene.append_polygon_mesh_vertex((0, y, 0))
		scene.append_polygon_mesh_vertex((-x, y, 0))
		scene.append_polygon_mesh_vertex((-x, -y, 0))
			
		scene.append_polygon_mesh_edge(0, 1)
		scene.append_polygon_mesh_edge(1, 2)
		scene.append_polygon_mesh_edge(2, 3)
		scene.append_polygon_mesh_edge(3, 0)
		scene.append_polygon_mesh_edge(4, 5)
		scene.append_polygon_mesh_edge(5, 6)
		scene.append_polygon_mesh_edge(6, 7)
		scene.append_polygon_mesh_edge(7, 4)
	
		scene.end_polygon_mesh()
		scene.end_creating()
		
		ob = scene.active_shape()
		ob.name = str(data)
	ob.dad.activate()

ob1.activate()
#scene.reset_all_transformation()

ob1.copy_object(None)
ob2 = ob1.bro
ob2.activate()

scene.create_part()
ob3 = ob2.bro
ob3.transform(((10,0,0,0), (0,100,0,0),(0,0,100,0),(0,0,0,1)))

ob2.place_child(1)
ob2.activate()
scene.reset_all_transformation()

ob1.name = '単位マトリックス'
ob3.name = '非 - 単位マトリックス'
scene.active_shapes = (ob1, ob3)

scene.allow_update()
scene.force_update()




2)ポリゴンメッシュの頂点マージ cleanup_redundant_vertices()

下記のスクリプトでは、スケールの異なるポリゴンメッシュの一つの頂点だけを少しずらし、cleanup_redundant_vertices() を掛ける前後での頂点数をレポートします。

前半は上位パートにマトリックスがかかっていないケース
後半は掛かっているケースです。


#  ポリゴンメッシュの同一頂点座標判定基準を確認する

#  互いに微小の距離だけ離れた頂点を持つポリゴンメッシュに対して
#  同一頂点削除を cleanup_redundant_vertices() で行い、
#  同一頂点判定の基準を調べる。

#  結果 : 同一頂点座標であるか否かは 絶対値で判定している。
#            しかも判定に上位マトリックスが影響している。



scene = xshade.scene()

Ls = []

L = []
L.append((1, 0.1))
L.append((1, 0.01))
Ls.append(L)

L = []
L.append((10, 0.1))
L.append((10, 0.01))
Ls.append(L)

L = []
L.append((100, 0.1))
L.append((100, 0.01))
Ls.append(L)



obL = []

scene.inhibit_update()

scene.exit_modify_mode()
scene.select_all()
scene.active_shape().transformation_matrix = ((1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1))

scene.create_part()
ob1 = scene.active_shape()
#ob1.cancel_transformation()

for L in Ls :
	scene.create_part()
	for data in L :
		y = float(data[0])/2
		x = y*2
		epsilon = data[1]
		
		scene.begin_creating()
		scene.begin_polygon_mesh(str(data))
	
		scene.append_polygon_mesh_vertex((0, y , epsilon))
		scene.append_polygon_mesh_vertex((0, -y, 0))
		scene.append_polygon_mesh_vertex((x, -y, 0))
		scene.append_polygon_mesh_vertex((x, y, 0))
	
		scene.append_polygon_mesh_vertex((0, -y, 0))
		scene.append_polygon_mesh_vertex((0, y, 0))
		scene.append_polygon_mesh_vertex((-x, y, 0))
		scene.append_polygon_mesh_vertex((-x, -y, 0))
			
		scene.append_polygon_mesh_edge(0, 1)
		scene.append_polygon_mesh_edge(1, 2)
		scene.append_polygon_mesh_edge(2, 3)
		scene.append_polygon_mesh_edge(3, 0)
		scene.append_polygon_mesh_edge(4, 5)
		scene.append_polygon_mesh_edge(5, 6)
		scene.append_polygon_mesh_edge(6, 7)
		scene.append_polygon_mesh_edge(7, 4)
	
		scene.end_polygon_mesh()
		scene.end_creating()
		
		ob = scene.active_shape()
		ob.name = str(data)
		obL.append(ob)
	ob.dad.activate()

ob1.activate()
#scene.reset_all_transformation()

ob1.copy_object(None)
ob2 = ob1.bro
ob2.activate()

ob = ob2.son
while ob.has_bro :
	ob = ob.bro
	ob_ = ob.son
	while ob_.has_bro :
		ob_ = ob_.bro
		obL.append(ob_)

scene.create_part()
ob3 = ob2.bro
ob3.transform(((10,0,0,0), (0,100,0,0),(0,0,100,0),(0,0,0,1)))

ob2.place_child(1)
ob2.activate()
scene.reset_all_transformation()

ob1.name = '単位マトリックス'
ob3.name = '非 - 単位マトリックス'
scene.active_shapes = (ob1, ob3)

scene.allow_update()
scene.force_update()




count = 0
for ob in obL :
	n1 = ob.total_number_of_control_points
	ob.cleanup_redundant_vertices()
	n2 = ob.total_number_of_control_points
	
	print str(ob.name), '      ', n1, '   ', n2
	count = count + 1
	if count%2 == 0 :
		print ' '




3)線形状からポリゴンメッシュへの変換時の頂点マージ

下記のスクリプトでは、スケールの異なる線形状の一つの頂点だけを少しずらし、ポリゴンメッシュに変換してその頂点数をレポートします。

前半は上位パートにマトリックスがかかっていないケース
後半は掛かっているケースです。


#  線形状からポリゴンメッシュへの変換時の同一頂点座標判定基準を確認

#  互いに微小の距離だけ離れた頂点を持つ2つの線形状をポリゴンメッシュに変換し、
#  その際の同一頂点判定の基準を調べる。

#  結果 : 同一頂点座標であるか否かは 絶対値で判定している。
#            しかも判定に上位マトリックスが影響している。



scene = xshade.scene()

Ls = []

L = []
L.append((1, 0.1))
L.append((1, 0.01))
Ls.append(L)

L = []
L.append((10, 0.1))
L.append((10, 0.01))
Ls.append(L)

L = []
L.append((100, 0.1))
L.append((100, 0.01))
Ls.append(L)



obL = []

scene.inhibit_update()

scene.exit_modify_mode()
scene.select_all()
scene.active_shape().transformation_matrix = ((1,0,0,0),(0,1,0,0),(0,0,1,0),(0,0,0,1))

scene.create_part()
ob1 = scene.active_shape()
#ob1.cancel_transformation()

for L in Ls :
	scene.create_part()
	for data in L :
		y = float(data[0])/2
		x = y*2
		epsilon = data[1]

		scene.create_part(str(data))
		ob = scene.active_shape()
		obL.append (ob)
		
		scene.begin_creating()
		
		scene.begin_line(None, True)
		scene.append_point((0, y , epsilon))
		scene.append_point((0, -y, 0))
		scene.append_point((x, -y, 0))
		scene.append_point((x, y, 0))
		scene.end_line()
	
		scene.begin_line(None, True)
		scene.append_point((0, -y, 0))
		scene.append_point((0, y, 0))
		scene.append_point((-x, y, 0))
		scene.append_point((-x, -y, 0))
		scene.end_line()
		
		scene.end_creating()
		
		ob.activate()
	ob.dad.activate()

ob1.activate()
#scene.reset_all_transformation()

ob1.copy_object(None)
ob2 = ob1.bro
ob2.activate()

ob = ob2.son
while ob.has_bro :
	ob = ob.bro
	ob_ = ob.son
	while ob_.has_bro :
		ob_ = ob_.bro
		obL.append(ob_)

scene.create_part()
ob3 = ob2.bro
ob3.transform(((10,0,0,0), (0,100,0,0),(0,0,100,0),(0,0,0,1)))

ob2.place_child(1)
ob2.activate()
scene.reset_all_transformation()

ob1.name = '単位マトリックス'
ob3.name = '非 - 単位マトリックス'
scene.active_shapes = (ob1, ob3)




count = 0
for ob in obL :
	ob.activate()
	ob.convert_to_polygon_mesh_with_subdivision_level(0)
	ob = scene.active_shape()
		
	print str(ob.name), '      ', ob.total_number_of_control_points
	count = count + 1
	if count%2 == 0 :
		print ' '



scene.allow_update()
scene.force_update()
投稿数: 189
投稿日時: 2011-10-07 14:50
Re: 線形状(ベジェ曲線)の分割について教えてください
< 平行判定 - ゼロ判定 >


Bezier Patch 構築ロジックの中で 0 判定が bounding box 基準で行われていることについて述べましたが、ここでは 0 判定のもう一つの代表的な例として平行判定を考えてみます。



二つのベクトルの平行判定は、それらの外積で調べます。

  v1 = ( x1, y1, z1 )
  v2 = ( x2, y2, z2 )

  v = v1 × v2 ( 通常表記では、外積を × ( cross ) で表します )
   = ( x, y, z )

      x = y1*z2 - z1*y2
      y = z1*x2 - x1*z2
      z = x1*y2 - y1*x2

外積の意味については参考サイトをかかげておきます。

  http://yosshy.sansu.org/gaiseki.htm



で、この外積で与えられるベクトル v が

  | v | = 0 であれば、v1, V2 は平行( 互いの向きは問わず )

となります。



誤差を考慮した平行判定として色々な考え方があると思いますが、2例ばかりあげてみます。


A)長さで判定

  二つのベクトルの先端を結ぶ線分の長さを基準長 L とし、

  この基準線と二つのベクトルの起点との距離 r が

  基準長 L の epsilon 倍以内、あるいは delta 倍以上であれば、平行と見なす。

  注意 :
    片方のベクトル長さが他方に比べて極端に短いと、
    その向きに関わりなく平行と見なされてしまう可能性がある



B)角度で判定

  sin(θ) が 許容値以内であれば、平行と見なす。








Shade の Bezier 関連での平行判定の例として、ハンドルのリンク情報があげられます。

スクリプトで線形状を描くとき、あるアンカーポイントの inhandle とouthandle とが直線状に並んでいると Shade が判断すると、それらのハンドルは自動的にリンクされる仕様になっています。

リンク情報は xshade.scene.active_shape().control_point(XX).linked で取得できます。


この Shade によるハンドルの平行判定には 長さ基準の判定 が用いられています。

   in handle と outhandle の先端を結ぶ線分 L を基準長とし、

   anchor point と線分 L との距離 r が、

     r <= 0.00025*L あるいは r >= 1000*L

   であれば平行とみなし、linked = true とする

      r / L = | V1 × V2 | / L^2




次のスクリプトを作ってみました。

このスクリプトでは

  ハンドル長さを変化させ、直線状からのずれの量を変化させながらサンプル線形状を描き、
  当該アンカーポイントでの、linked 情報を調べるとともに、

  上記のしきい値を用いた 長さ判定 による平行判定と Shade が返す lnked が等しいことを示し、

  比較のためにもう一つの判定法 角度基準 による平行判定の結果もレポートします。



レポート結果から、

  長さ基準と角度基準とでは平行判定結果が異なる場合がある。( まあ、当たり前ですが )

  長さ基準を用いた Shade の linked 判定では、 片方のハンドル長さが極端に短い
  次のようなケースでも平行とみなされてしまうことがある。

     in handle  anchor point  out handle
     (-1, 0, 0)   (0, 0, 0)   (0, 0.0002, 0)

  これを避けるには、

     最小ハンドル長さ、あるいは、最小ハンドル長さ比( 短いハンドルと長いハンドルの長さの比 )
     等を指定して、平行判定の前にはじいてしまうロジックを追加する

     あるいは、角度基準の判定式を用いれる。




参考としてスクリプトの中に書いておいた 角度基準 の判定では、次の許容値を用いています。

  sin(θ) <= 0.001
  sin(θ) = | V1 × V2 | / ( | V1 | * | V2 | )

また、次の機能を付加してあります。

  許容最小ハンドル長さの指定
  許容ハンドル長さ比の指定
  方向を問わない平行判定, 同一方向の平行判定, 逆向きの平行判定 の選択



#  平行判定

#   inhandle と outhandle を直線状から僅かにずらして線形状を作成し、
#   control_point(X).linked で直線と見なしているか否かを調べ、
#   同時に 長さ基準と角度基準の平行判定式でも評価し、3つを比較する
#
#   長さ基準の平行判定式内の係数は、評価結果が Shade の linked と
#   同じになるように定めてある
#
#   角度基準の平行判定式では、許容最小ハンドル長さ、許容最小ハンドル長さ比、方向別の平行判定指定 を与えることができる。


#  結果 : linked 判定は
#	inhadle, outhandle 先端を結ぶ線分の長さを基準長とし、
#	その線分とアンカーポイントとの距離が
#	基準長の 0.00025 倍以下
#	あるいは 1000倍以上であれば
# 	直線状に並んでいると見なし linked = True を返す





#####  入力  #####
#  inhandle 座標, outhandle 座標, 回転

Ls = []

L = []
inh = (-0.5, 0, 0)
outh = (0.5, 0.0005001, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (-0.5, 0, 0)
outh = (0.5, 0.0004999, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)



L = []
inh = (-5, 0, 0)
outh = (5, 0.005001, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (-5, 0, 0)
outh = (5, 0.004999, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)



L = []
inh = (-2.5, 0, 0)
outh = (0.5, 0.0009001, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (-2.5, 0, 0)
outh = (0.5, 0.0008999, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)



L = []
inh = (-1, 0, 0)
outh = (0, 0.0003, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (-1, 0, 0)
outh = (0, 0.0002, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)



L = []
inh = (-0.5, 0, 0)
outh = (-0.1, 0, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (-0.5, 0, 0)
outh = (-0.5, 0, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)



L = []
inh = (-1, 0, 0)
outh = (-0.99999, 0.001, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (-1, 0, 0)
outh = (-0.99999, 0.00099999, 0)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)



e = 2.6469781e-23
f = 2.6469782e-23

L = []
inh = (0, 0, 0)
outh = (e, e, e)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))

inh = (0, 0, 0)
outh = (f, e, e)
rotation = [0, 0, 0]
L.append((inh, outh, rotation))
Ls.append(L)
#################



def sub_vec(v1, v2):		#  ベクトルの差
	return (v1[0]-v2[0], v1[1]-v2[1], v1[2]-v2[2])



def cross_vec(v1, v2) :		#  ベクトルの外積
	return[v1[1]*v2[2] - v1[2]*v2[1], v1[2]*v2[0] - v1[0]*v2[2], v1[0]*v2[1] - v1[1]*v2[0]]



def check_vector_length(v) :		#  ベクトル長さの 0 判定
#	delta = 0	#  python では double で計算されるので、0 で比較すると float で計算する Shade と結果が異なる場合がある
	delta = 2.10194792e-45
	return v[0]**2  + v[1]**2 + v[2]**2 > delta




def check_parallel_by_length(u, v) :	#  ベクトル u, v の平行判定 ( 長さ基準 )

	###  基準長に対する しきい値
	epsilon_1 = 0.00025
	epsilon_2 = 1000
	
	###  a ###	
	if (not check_vector_length(u)) and (not check_vector_length(v))  :  
		return False, 'ハンドル長さが 0 '		#  u, v 両方の長さが極端に短いケースは False と判定
	else :
		s = cross_vec(u, v)
		a = s[0]**2 + s[1]**2 + s[2]**2		
	
	###  b ###
	w = sub_vec(u, v)	
	if not check_vector_length(w) :	
		return True, 'inhandle, outhandle が等しい'	#  inhandle, outhandle がほとんど等しければ、平行
	else :
		b = w[0]**4 + w[1]**4 + w[2]**4
	
	###  r  ###	
	r = a/b
	if r <= epsilon_1**2 :
		return True, r
	elif r >= epsilon_2**2 :	#  これは inhandle, outhandle が同一方向で接近しているケースを想定
		return True, r
	else :
		return False, r
		



def check_parallel_by_angle(u, v, direction_, min_L_, min_ratio_) :	#  ベクトル u, v の平行判定 ( 角度基準 )
							#  direction  0 : 方向を問わない、1: 同一方向、-1: 逆方向
							#  min_L : 最小ベクトル長さ  これより短ければ False を返す
							#  min_raio : 最小ベクトル長さ比  ベクトル長さの比がこれを越えれば False を返す
	###  default 値
	if direction_ == None :
		direction = 0
	else :
		direction = direction_
		
	if min_L_ == None :
		min_L = 0
	else :
		min_L = mon_L_
		
	if min_L_ == None :
		min_ratio = 0		#  min_ratio = 0 ではベクトル長さ比のチェックを行わないことにします。
	else :
		min_ratio = mon_ratio_
		
	###  許容値
	allow = 0.001	#  交差角の sin 値
	
	###  a, b
	a = u[0]**2 + u[1]**2 + u[2]**2
	b = v[0]**2 + v[1]**2 + v[2]**2
	
	###  最小長さチェック
	if min_L != 0 :
		if a < min_L**2 or b < min_L**2 :
			return False, 'ハンドル長さが 0 '	#  いずれかのハンドルの長さが 0
		
	###  最小長さ比チェック
	if min_ratio != 0 :
		if a + b == 0 :
			return False, 'ハンドル長さが 0'
		elif a == 0 or b == 0 :
			return False, 'ハンドル長さ比が小さい'
		else :
			c = a/b
			if c > 1 :
				c = b/a
			if c < min_ratio :
				return False, 'ハンドル長さ比が小さい'

	###  平行判定
	if a*b == 0 :
		return False, 'ハンドル長さが短すぎる '
	else :
		s = cross_vec(u, v)
		area = s[0]**2 + s[1]**2 + s[2]**2
		r = area/(a*b)
		
	if direction == 0 :	#  方向を問わない
		return r <= allow**2, r
		
	elif direction == 1 :	#  同一方向のみ判定
		if u[0]*v[0] + u[1]*v[1] + u[2]*v[2] > 0 :
			return r <= allow**2, r
		else :
			return False, '方向ミスマッチ'
			
	elif direction == -1 :	#  逆方向のみ判定
		if u[0]*v[0] + u[1]*v[1] + u[2]*v[2] < 0 :
			return r <= allow**2, r
		else :
			return False, '方向ミスマッチ'
		


def main() :
	scene.inhibit_update()

	scene.create_part()
	ob0 = scene.active_shape()
	ob0.cancel_transformation()

	print ' '                        
	print 'linked,   長さ基準判定,  r / L,  ihandle,   outhandle'
	print '              角度基準判定,  sin(θ)'
	print ' '

	for L in Ls :
		scene.create_part()
	
		scene.begin_creating()
		ob = scene.active_shape()
		ob.disclosed = False
	
		for (inh, outh, rotation) in L :		
			scene.begin_line(None, False)
			scene.append_point((0, 0, 0), inh, outh, None, None)
			scene.append_point((1, 0, 0), inh, outh, None, None)
			scene.end_line()
		
			if rotation != [0, 0, 0] :
				xshade.scene().move_object([0, 0, 0], None, [0, 15, 30], None)

			shade_linked =  (scene.active_shape().control_point(0).linked == 1)		#  Shade の lined 判定
			linked, result = check_parallel_by_length(inh, outh)				#  平行判定 - 長さ基準 ( Shade 方式 )
			linked2, result2 = check_parallel_by_angle(inh, outh, None, None, None)	#  平行判定 - 角度基準
			
			###  レポート
			if isinstance(result, str) :
				print shade_linked, '    ', linked, '      ', result, '         ', inh, '   ', outh
			else :
				print shade_linked, '    ', linked, '      ', 'r / L = ', result**(0.5), '         ', inh, '   ', outh
				
			if isinstance(result2, str) :
				print '              ', linked2, '      ', result2
			else :
				print '               ', linked2, '      ', 'sin(θ) = ', result2**(0.5)
			
				
		scene.end_creating()
		ob.activate()
		print ' '
	
	ob0.activate()
	ob0.reset_transformation()

	scene.allow_update()
	scene.force_update()



scene = xshade.scene()
main()
投稿数: 106
投稿日時: 2011-10-01 23:29
Re: 線形状(ベジェ曲線)の分割について教えてください
重複関数があっても、単純上書きされるだけだと思いますので
問題ないかと。

ということで、PDFにして保存させていただきました。w
あざっす。
投稿数: 16
投稿日時: 2011-10-01 21:42
Re: 線形状(ベジェ曲線)の分割について教えてください
加藤さん、

楽しみにしてました^^
すごい情報量ですよねーー
わくわくしてきました♪

じっくりと、かみしめて読んでみたいと思います。


ちなみに関数の定義の重複はwindowsでも問題なかったです。
後の定義で上書きされている感じですね。
Pythonの仕様なんでしょうか?

投稿数: 189
投稿日時: 2011-10-01 19:17
Re: 線形状(ベジェ曲線)の分割について教えてください
あらら、最後のスクリプトの中で

def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])


def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])



って、同じものがダブルで入ってますねえ。

でも、エラーにはならないのか。

もし Windows で走らせてエラーになったら、お手数ですが、片方を削除して動かして下さい。
投稿数: 189
投稿日時: 2011-10-01 18:59
Re: 線形状(ベジェ曲線)の分割について教えてください
< Bezier Patch 補正  6 / 6 >


以上で Shade の Bezier Patch 構築方法が突き止められましたので、最後に自由曲面の指定する区間の Bezier Patch を DXF 経由と内部計算の二つで出力するスクリプトを書いておきます。

タブン、このロジックで 99.9% 合っていると思うのですが、ひょっとしたら未だ気づいていない ナニカ があるかもしれません。



もう一度、ロジックをまとめると、


[[ Shade の Bezier Patch の構築方 ]]

  1)16個の制御点の内、アンカーポイントとハンドル座標で直接求められない4つの制御点については
    A = H + ( h - p ) ルールで求める
   


  2)ハンドル長さの 0 判定を次の方法で行う

      アンカーポイントとハンドル座標で定義される bounding box を求め、
      その X, Y, Z 3方向の長さの和を基準長とし、
      基準長 × 1/10000 よりも短いハンドル長さは 0 とみなす


  3)一つのアンカーポイントについて、その交差する2つのハンドルの内、
    いずれか1つのハンドル長さが 0 と判定されたなら、Patch 補正を行う。
    補正は A = h0 + (B - h1)/2 ルールによる
   




最後のハンドル長さの 0 判定のところでは、かなり細かいところまで調べました。

言うまでもなく、図形ウィンド上で手作業でモデリングされた自由曲面に対してなら、マウス移動量がモニターの解像度に依存しますので、そんな細かい所まで考える必要はないでしょう。

でも、場合によっては僅かな誤差が結果に大きく影響する可能性も皆無ではありません。


ここで触れた「 長さの 0 判定 」あるいは「 同一座標判定 」は、一般ユーザーが意識する以上に重要で、判定ロジックはプログラム上で重要な位置を占めます。

このことに関する Shade の根底に流れる基本思想や、最近の困った傾向、などについて、また後日に述べたいと思います。




#  自由曲面の指定するブロックの Bezier Patch を出力する

#  DXF file 経由の Bezier Patch と 内部計算で求めた Bezier Patch の両方を出力する
#  内部計算には、ハンドル長さの0判定、Patch 補正 処理を含む。


#  script では直接 Bezier Patch データを入手できないので
#  一旦 scripts folder 内に DXF ファイルとして temporary file として
#  エクスポートし、それを読み込んで Patch データを入手する



###  注意!!!  #####################################
#  DXF エクスポートの際のプロパティー設定をスクリプト側から操作できないので、
#  事前にマニュアルで設定しておく必要がある。

#  メニュー > ファイル > エクスポート > DXF を選択し、
#  DXF エクスポート設定の欄の 「ベジェ曲面 」 にチェックを入れ、
#  OK ボタンを押して、ファイル出力ダイアログが出たら キャンセル を選ぶ。
#  一旦この操作をしておけば、Shade を再起動後も設定は保たれたままとなる。



####    入力     #################################################
#  対象とする Bezier Patch 区間 No.を (U, V ) で指定 ( 0基数 )
block = (0, 0)
############################################################



def sub_vec(va, vb) :		#  ベクトル減算
	return (va[0]-vb[0], va[1]-vb[1], va[2]-vb[2])
	

	
def check_object(ob, block) :	#  選択オブジェクトの種類と、その構成が指定する block No.とマッチするかチェック

	if ob.type != 2 or ob.part_type != 1 :
		print '自由曲面を選択して下さい。'
		return False
	elif ob.number_of_sons < 2 :
		print 'この自由曲面は面を構成しません。'
		return False
		
	if block[0] < 0 or block[1] < 0 :
		print 'block No. の指定が適切ではありません。'
		return False
		
	ob2 = ob.son.bro
	if ob2.closed :
		k = ob2.number_of_control_points - 1
	else :
		k = ob2.number_of_control_points - 2
	if k < block[0] :
		print 'block No. の指定が適切ではありません。'
		return False
		
	if ob.surface_closed :
		k = ob.number_of_sons - 1
	else :
		k = ob.number_of_sons - 2
	if k < block[1] :
		print 'block No. の指定が適切ではありません。'
		return False
		
	return True


		
			
	
def get_block_surface(ob, (u1, v1)) :	# block No. u1, v1 で与えられる区間のみを取り出した自由曲面を作成し、返す	

	###  u2, v2
	ob2 = ob.son.bro
	if ob2.closed and u1 == ob2.number_of_control_points - 1 :
		u2 = 0
	else :
		u2 = u1 + 1
		
	if ob.surface_closed and v1 == ob.number_of_sons - 1 :
		v2 = 0
	else :
		v2 = v1 + 1
		
	###  block = (u1, v1) で指定される区間のみを新たに自由曲面として作成する
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	scene.begin_creating()
	scene.begin_surface_part()
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob3 = ob2.bro

	scene.begin_line(None, False)
	
	ob3 = ob2.son
	for i in range(v1 + 1) :
		ob3 = ob3.bro
	(ap, out_H, lout_H) = (
		ob3.control_point(u1).position,
		ob3.control_point(u1).out_handle,
		ob3.control_point(u1).lateral_out_handle)
	scene.append_point(ap, None, out_H, None, lout_H)
	
	(ap, in_H, lout_H) = (
		ob3.control_point(u2).position,
		ob3.control_point(u2).in_handle,
		ob3.control_point(u2).lateral_out_handle)
	scene.append_point(ap, in_H, None, None, lout_H)

	scene.end_line()
	
	scene.begin_line(None, False)
	
	ob3 = ob2.son
	for i in range(v2 + 1) :
		ob3 = ob3.bro
	(ap, out_H, lin_H) = (
		ob3.control_point(u1).position,
		ob3.control_point(u1).out_handle,
		ob3.control_point(u1).lateral_in_handle)
	scene.append_point(ap, None, out_H, lin_H, None)
	
	(ap, in_H, lin_H) = (
		ob3.control_point(u2).position,
		ob3.control_point(u2).in_handle,
		ob3.control_point(u2).lateral_in_handle)
	scene.append_point(ap, in_H, None, lin_H, None)

	scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob2.remove()	#  ob の複製を削除
	
	return ob0	#  作成した自由曲面を返す
		



def get_bezier_patch_by_DXF(ob) :	#  自由曲面 ob の bezier Patch を DXF file 経由で返す
	import os
	
	ob.activate()
	name = ob.name	#  テキストで書かれた DXFfile 読み取り後のの文字列分割において、オブジェクト名称にスペースが含まれると、
	ob.name = ''		#  そこでも分割されてしまうため、デフォルト名称に変更しておいて DXF file をエクスポートする。
	
	path = xshade.shade().scripts_path + '/temp999.dxf'	#  temp file の path
	scene.save_DXF(path)					# scripts folder 内に temp file として DXF file を出力
	ob.name = name	#  オブジェクト名称復活
	if not os.path.exists :
		return None
	
	f = open(path)
	t = f.read()
	f.close()
	os.remove(path)	#  temporary file 削除	
	
	s = t.split()
	if (s[11] !='16') or (s[19] != 'VERTEX') :
		print 'DXF エクスポートの設定として、「ベジェ曲面」にチェックを入れておいて下さい。'
		return None
	else :
		bP = []
		k = 23
		for i in range(4) :
			v = []
			for j in range(4) :
				v.append((float(s[k]), float(s[k + 4]), -float(s[k + 2])))	#  DXF file の座標は傾いている
				k = k + 12
			bP.append(v)
		return bP





def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])




def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])



def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)




def get_standard_length(bP) :	#  ob の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if not ((i == 1 and j == 1) or (i == 1 and j == 2) or (i == 2 and j == 1) or (i == 2 and j == 2)) :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
	
	
	

def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	
#	ハンドル長さの0判定に has_handle( ) を用いない
#	has_uH = (
#		(ob2.control_point(0).has_out_handle, ob2.control_point(1).has_in_handle),
#		(ob2.control_point(2).has_out_handle, ob2.control_point(3).has_in_handle))
#	has_vH = (
#		(ob2.control_point(0).has_lateral_out_handle, ob2.control_point(1).has_lateral_out_handle),
#		(ob2.control_point(2).has_lateral_in_handle, ob2.control_point(3).has_lateral_in_handle))
	
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP

	


	
def write_bezier_patch(bP, name) :

	scene.begin_creating()
	scene.begin_surface_part()
	ob = scene.active_shape()
	ob.name = name
	ob.disclosed = False
	ob.cancel_transformation()	
	
	for i in range(4) :
		scene.begin_line(None, False)
		for j in range(4) :
			scene.append_point(bP[i][j], None, None, None, None)
		scene.end_line()
			
	scene.end_surface_part()
	scene.end_creating()
	
	ob.reset_transformation()
	
	return ob
	
	
	

def main(block) :
	ob = scene.active_shape()
	
	if not check_object(ob, block) :		#  選択オブジェクトの適合性を確認
		return
	
	scene.inhibit_update()
	
	scene.create_part()				#  出力全体を格納するパート
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	ob2 = get_block_surface(ob, block)		#  block = (u1, v1) で与えられる区間のみを取り出し、作成した temporary 自由曲面
	ob2.activate()
	
		
	
	bP = get_bezier_patch(ob2)			#  temporary 自由曲面の bezier Patch を 計算により取得
	ob3 = write_bezier_patch(bP, 'Bezier Patch')	#  Bezier Patch を出力	
		
	bP = get_bezier_patch_by_DXF(ob2)		#  temporary 自由曲面の bezier Patch を DXF file 経由で取得
	if bP != None :
		ob4 = write_bezier_patch(bP, 'Bezier Patch by DXF')	#  Bezier Patch を出力

	ob2.remove()				#  temporary 自由曲面削除			
	ob0.activate()
	scene.reset_all_transformation()
	
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main(block)
投稿数: 189
投稿日時: 2011-10-01 18:47
Re: 線形状(ベジェ曲線)の分割について教えてください
< Bezier Patch 補正  5 / 6 >


あるオブジェクトに対して、そのあるパーツの長さを求めるのに、

  その長さが微小の場合には 0 とみなした方が好ましい

  あるいは 0 とみなさないとロジックが破綻する

という状況は頻繁に遭遇します。

例えば、分母が 0 になるのをさけるために、考えられる諸々の誤差を考慮した上で、ある微小な値以下になる場合を避けたり、分岐処理したりする場合です。


この 0 判定はいささか面倒で、判定基準の定め方も、どのような処理を行おうとしているかによって変わることもあります。

この話については Bezier Patch 補正の話が終わった後に、改めて取り上げるつもりです。



で、Bezier Patch の組み立て時におけるハンドル長さの 0 判定ですが、最も一般的な判定基準は対象となるオブジェクトの bounding box size を基準とするものです。

つまり、求める長さがその bounding box size の XX % 以下ならば長さが 0 であると判定します。


bounding box size といっても 縦, 横, 高さ( X, Y, Z )3方向の長さがありますが、なにをもって基準長とするかは色々な考え方があります。

  X, Y, Z の最大値
  X, Y, Z の平均値
  box の対角方向の長さ = sqrt(X^2 + Y^2 + Z^2)
  X, Y, Z の和


また自由曲面の bounding box size の求め方にも
 
  a)アンカーポイントとハンドル座標から定義する簡易的なもの
  b)面の曲がり具合も考慮して、本当の意味での包括矩形領域を求めるもの

など、いくつか考えられます。

( ちなみに Bezier の内包性といって、必ず a) >= b) となります )



色々実験してみると、Shade の Bezier Patch の 0 判定は次のようになっています。

    当該 Bezier 区間の

      4つのアンカーポイント座標,
      4つのハンドル座標,
      4つの交差ハンドル座標

    を元に bounding box size を求め、

    bounding box size の「 X, Y, Z の和 」を基準長とし、

    基準長 の 1/10000 を 0 判定のしきい値とする


次のスクリプトでは has_handle( ) によるハンドル長さの判定ではなく、上記の 0 判定ロジックに置き換え、サンプル自由曲面に対して Bezer Patch を求めています。


スクリプトの結果として、3つのモデルに対して上記のロジックが正しいであろうことが示されています。

  モデル1, 2
    bounding box size = ( 1500, 1500, 1000)
    その和 = 4000 = 1500 + 1500 + 1000
    しきい長さ = 0.4 = 4000 × 1/10000
    
微小ハンドル長さ 0.4001 → ハンドル長さありと判定
微小ハンドル長さ 0.3999 → ハンドル長さ 0 と判定


  モデル3, 4
    bounding box size = ( 150000, 150000, 100000)
    その和 = 400000 = 150000 + 150000 + 100000
    しきい長さ = 40 = 400000 × 1/10000
    
微小ハンドル長さ 40.01 → ハンドル長さありと判定
微小ハンドル長さ 39.99 → ハンドル長さ 0 と判定


  モデル 5, 6
    bounding box size = ( 1457.1, 457.1, 2121.3)
    その和 = 4035.5 = 1457.1 + 457.1 + 2121.3
    しきい長さ = 0.40355 = 4035.5 × 1/10000
    
微小ハンドル長さ 0.4036 → ハンドル長さありと判定
微小ハンドル長さ 0.4035 → ハンドル長さ 0 と判定



なお、スクリプト のプロパティ xshade.scene().active_shape().bounding_box_size で得られる値はは上の方法で求める bounding box size とは値が異なります。


次回は Bezier Patch に係わるまとめです。




#  サンプル自由曲面の Bezier Patch を出力する

#  Patch 補正の必要性判定のための、ハンドル長さの0判定を
#  has_handle( ) を使用せず、バウンディングボックスサイズ基準にする



###  注意!!!  #####################################
#  DXF エクスポートの際のプロパティー設定をスクリプト側から操作できないので、
#  事前にマニュアルで設定しておく必要がある。

#  メニュー > ファイル > エクスポート > DXF を選択し、
#  DXF エクスポート設定の欄の 「ベジェ曲面 」 にチェックを入れ、
#  OK ボタンを押して、ファイル出力ダイアログが出たら キャンセル を選ぶ。
#  一旦この操作をしておけば、Shade を再起動後も設定は保たれたままとなる。



####    入力     #################################################
#  サンプル自由曲面
bLs = []

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (0.4001, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = '長さ 0.4001 のハンドル'
bLs.append((p1, p2, p3, p4, name, (0, 0, 0)))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (0.39999, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = '長さ 0.39999 のハンドル'
bLs.append((p1, p2, p3, p4, name, (0, 0, 0)))


p1= ((0, 100000, 0), None, (0, 100000, 50000), None, (40.01, 100000, 0))
p2 = ((0, 0, 100000), (0, 50000, 100000), None, None, (50000, 0, 100000))
p3 = ((50000, 100000, -50000), None, (100000, 100000, -50000), None, None)
p4 = ((150000, 0, -50000), (150000, 50000, -50000), None, (150000, 0, 0), None)
name = '長さ 40.01 のハンドル'
bLs.append((p1, p2, p3, p4, name, (0, 0, 0)))

p1= ((0, 100000, 0), None, (0, 100000, 50000), None, (39.99, 100000, 0))
p2 = ((0, 0, 100000), (0, 50000, 100000), None, None, (50000, 0, 100000))
p3 = ((50000, 100000, -50000), None, (100000, 100000, -50000), None, None)
p4 = ((150000, 0, -50000), (150000, 50000, -50000), None, (150000, 0, 0), None)
name = '長さ 39.99 のハンドル'
bLs.append((p1, p2, p3, p4, name, (0, 0, 0)))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (0.4036, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = '長さ 0.4036 のハンドル'
bLs.append((p1, p2, p3, p4, name, (0, 45, 45)))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (0.4035, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = '長さ 0.4035 のハンドル'
bLs.append((p1, p2, p3, p4, name, (0, 45, 45)))

############################################################



def sub_vec(va, vb) :		#  ベクトル減算
	return (va[0]-vb[0], va[1]-vb[1], va[2]-vb[2])

		
	
def write_sample_surface(bL) :	# block No. u1, v1 で与えられる区間のみを取り出した自由曲面を作成し、返す	
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	scene.begin_line(None, False)
	scene.append_point(bL[0][0], bL[0][1], bL[0][2], bL[0][3], bL[0][4])
	scene.append_point(bL[1][0], bL[1][1], bL[1][2], bL[1][3], bL[1][4])
	scene.end_line()
	
	scene.begin_line(None, False)
	scene.append_point(bL[2][0], bL[2][1], bL[2][2], bL[2][3], bL[2][4])
	scene.append_point(bL[3][0], bL[3][1], bL[3][2], bL[3][3], bL[3][4])
	scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	if bL[4] != '' :
		ob.name = bL[4]
		
	if bL[5] != (0, 0, 0) :
		scene.move_object(None, None, bL[5], None)
		scene.reset_transformation()
	
	
	return ob 	#  作成した自由曲面を返す
		



def get_bezier_patch_by_DXF(ob) :	#  自由曲面 ob の bezier Patch を DXF file 経由で返す
	import os
	
	ob.activate()
	name = ob.name	#  テキストで書かれた DXFfile 読み取り後のの文字列分割において、オブジェクト名称にスペースが含まれると、
	ob.name = ''		#  そこでも分割されてしまうため、デフォルト名称に変更しておいて DXF file をエクスポートする。
	
	path = xshade.shade().scripts_path + '/temp999.dxf'	#  temp file の path
	scene.save_DXF(path)					# scripts folder 内に temp file として DXF file を出力
	ob.name = name	#  オブジェクト名称復活	
	if not os.path.exists :
		return None
	
	f = open(path)
	t = f.read()
	f.close()
	os.remove(path)	#  temporary file 削除	
	
	s = t.split()
	if (s[11] !='16') or (s[19] != 'VERTEX') :
		print 'DXF エクスポートの設定として、「ベジェ曲面」にチェックを入れておいて下さい。'
		return None
	else :
		bP = []
		k = 23
		for i in range(4) :
			v = []
			for j in range(4) :
				v.append((float(s[k]), float(s[k + 4]), -float(s[k + 2])))	#  DXF file の座標は傾いている
				k = k + 12
			bP.append(v)
		return bP





def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])



def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)




def get_standard_length(bP) :	#  ob の bounding box size から、基準長を求める
	
	b1 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	b2 = [bP[0][0][0], bP[0][0][1], bP[0][0][2]]
	for i in range(4) :
		for j in range(4) :
			if not ((i == 1 and j == 1) or (i == 1 and j == 2) or (i == 2 and j == 1) or (i == 2 and j == 2)) :
				for k in range(3) :
					b1[k] = max(b1[k], bP[i][j][k])
					b2[k] = min(b2[k], bP[i][j][k])
				
	return (b1[0] - b2[0]) + (b1[1] - b2[1]) + (b1[2] - b2[2])	#  bounding box size の和を基準長とする
		
	

def check_handle(p1, p2, std_L) :	#  ハンドル長さの有無を基準長ベースに判定する
	v = sub_vec(p2, p1)
	L = v[0]**2 + v[1]**2 + v[2]**2	#  v の長さの2乗
	epsilon = 0.0001			#  ハンドル長さが 基準長 × epsilon 以下ならば、ハンドル長さがないとみなす
	
	return L/(std_L**2) > epsilon**2



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	
#	ハンドル長さの0判定に has_handle( ) を用いない
#	has_uH = (
#		(ob2.control_point(0).has_out_handle, ob2.control_point(1).has_in_handle),
#		(ob2.control_point(2).has_out_handle, ob2.control_point(3).has_in_handle))
#	has_vH = (
#		(ob2.control_point(0).has_lateral_out_handle, ob2.control_point(1).has_lateral_out_handle),
#		(ob2.control_point(2).has_lateral_in_handle, ob2.control_point(3).has_lateral_in_handle))
	
	std_L = get_standard_length(bP)	#  bounding box size 基準の基準長を求める
	has_uH = (
		(check_handle(bP[0][0], bP[0][1], std_L),
		check_handle(bP[0][3], bP[0][2], std_L)),
		(check_handle(bP[3][0], bP[3][1], std_L),
		check_handle(bP[3][3], bP[3][2], std_L)))
	has_vH = (
		(check_handle(bP[0][0], bP[1][0], std_L),
		check_handle(bP[0][3], bP[1][3], std_L)),
		(check_handle(bP[3][0], bP[2][0], std_L),
		check_handle(bP[3][3], bP[2][3], std_L)))

	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	


	
def write_bezier_patch(bP, name) :

	scene.begin_creating()
	scene.begin_surface_part()
	ob = scene.active_shape()
	ob.name = name
	ob.disclosed = False
	ob.cancel_transformation()	
	
	for i in range(4) :
		scene.begin_line(None, False)
		for j in range(4) :
			scene.append_point(bP[i][j], None, None, None, None)
		scene.end_line()
			
	scene.end_surface_part()
	scene.end_creating()
	
	ob.reset_transformation()
	
	return ob
	
	
	
def main() :

	scene.inhibit_update()
	
	scene.create_part('ハンドル長さの0判定を別途考慮')
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	for bL in bLs :
		scene.create_part()
		ob = write_sample_surface(bL)
		bP = get_bezier_patch(ob)			#  temporary 自由曲面の bezier Patch を 計算により取得
		ob2 = write_bezier_patch(bP, 'Bezier Patch')	#  Bezier Patch を出力
		bP = get_bezier_patch_by_DXF(ob)		#  temporary 自由曲面の bezier Patch を DXF file 経由で取得
		if bP != None :
			ob3 = write_bezier_patch(bP, 'Bezier Patch by DXF')	#  Bezier Patch を出力
		scene.active_shape().dad.activate()

	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main()
投稿数: 189
投稿日時: 2011-10-01 18:32
Re: 線形状(ベジェ曲線)の分割について教えてください
< Bezier Patch 補正  4 / 6 >


次のスクリプトを実行してみます。

このスクリプトでは、微小長さのハンドルを持つサンプル自由曲面を作り、その Bezier Patch を DXF 経由と、ここまで解析してきたロジックで計算したものとを出力し、比較します。


結果は、

  Shade は微小のハンドル長さに対しては、「 ハンドル長さが 0 である 」と見なして Bezier Patch を組み立てている。

  スクリプトの has_handle( ) でハンドル長さの有無を判定するのと、Shade が長さが無いと見なす基準は異なっている。


この Shade のハンドル長さ判定の基準がどのような仕組みになっているのか、次回で調べます。


#  サンプル自由曲面の Bezier Patch を出力する

#  微小長さのハンドルによる Patch 補正の有無を確認



###  注意!!!  #####################################
#  DXF エクスポートの際のプロパティー設定をスクリプト側から操作できないので、
#  事前にマニュアルで設定しておく必要がある。

#  メニュー > ファイル > エクスポート > DXF を選択し、
#  DXF エクスポート設定の欄の 「ベジェ曲面 」 にチェックを入れ、
#  OK ボタンを押して、ファイル出力ダイアログが出たら キャンセル を選ぶ。
#  一旦この操作をしておけば、Shade を再起動後も設定は保たれたままとなる。



####    入力     #################################################
#  サンプル自由曲面
bLs = []

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (1, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = '長さ 1 のハンドル'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, (0.1, 1000, 0))
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = '長さ 0.1 のハンドル'
bLs.append((p1, p2, p3, p4, name))



############################################################



		
			
	
def write_sample_surface(bL) :	# block No. u1, v1 で与えられる区間のみを取り出した自由曲面を作成し、返す	
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	scene.begin_line(None, False)
	scene.append_point(bL[0][0], bL[0][1], bL[0][2], bL[0][3], bL[0][4])
	scene.append_point(bL[1][0], bL[1][1], bL[1][2], bL[1][3], bL[1][4])
	scene.end_line()
	
	scene.begin_line(None, False)
	scene.append_point(bL[2][0], bL[2][1], bL[2][2], bL[2][3], bL[2][4])
	scene.append_point(bL[3][0], bL[3][1], bL[3][2], bL[3][3], bL[3][4])
	scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	if bL[4] != '' :
		ob.name = bL[4]
	
	return ob 	#  作成した自由曲面を返す
		



def get_bezier_patch_by_DXF(ob) :	#  自由曲面 ob の bezier Patch を DXF file 経由で返す
	import os
	
	ob.activate()
	name = ob.name	#  テキストで書かれた DXFfile 読み取り後のの文字列分割において、オブジェクト名称にスペースが含まれると、
	ob.name = ''		#  そこでも分割されてしまうため、デフォルト名称に変更しておいて DXF file をエクスポートする。
	
	path = xshade.shade().scripts_path + '/temp999.dxf'	#  temp file の path
	scene.save_DXF(path)					# scripts folder 内に temp file として DXF file を出力
	ob.name = name	#  オブジェクト名称復活	
	if not os.path.exists :
		return None
	
	f = open(path)
	t = f.read()
	f.close()
	os.remove(path)	#  temporary file 削除	
	
	s = t.split()
	if (s[11] !='16') or (s[19] != 'VERTEX') :
		print 'DXF エクスポートの設定として、「ベジェ曲面」にチェックを入れておいて下さい。'
		return None
	else :
		bP = []
		k = 23
		for i in range(4) :
			v = []
			for j in range(4) :
				v.append((float(s[k]), float(s[k + 4]), -float(s[k + 2])))	#  DXF file の座標は傾いている
				k = k + 12
			bP.append(v)
		return bP





def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])



def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(ob2.control_point(0).has_out_handle, ob2.control_point(1).has_in_handle),
		(ob2.control_point(2).has_out_handle, ob2.control_point(3).has_in_handle))
	has_vH = (
		(ob2.control_point(0).has_lateral_out_handle, ob2.control_point(1).has_lateral_out_handle),
		(ob2.control_point(2).has_lateral_in_handle, ob2.control_point(3).has_lateral_in_handle))
	
	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	


	
def write_bezier_patch(bP, name) :

	scene.begin_creating()
	scene.begin_surface_part()
	ob = scene.active_shape()
	ob.name = name
	ob.disclosed = False
	ob.cancel_transformation()	
	
	for i in range(4) :
		scene.begin_line(None, False)
		for j in range(4) :
			scene.append_point(bP[i][j], None, None, None, None)
		scene.end_line()
			
	scene.end_surface_part()
	scene.end_creating()
	
	ob.reset_transformation()
	
	return ob
	
	
	
def main() :

	scene.inhibit_update()
	
	scene.create_part('微小長さのハンドルによる Patch 補正')
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	for bL in bLs :
		scene.create_part()
		ob = write_sample_surface(bL)
		bP = get_bezier_patch(ob)			#  temporary 自由曲面の bezier Patch を 計算により取得
		ob2 = write_bezier_patch(bP, 'Bezier Patch')	#  Bezier Patch を出力
		bP = get_bezier_patch_by_DXF(ob)		#  temporary 自由曲面の bezier Patch を DXF file 経由で取得
		if bP != None :
			ob3 = write_bezier_patch(bP, 'Bezier Patch by DXF')	#  Bezier Patch を出力
		scene.active_shape().dad.activate()

	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main()
投稿数: 189
投稿日時: 2011-10-01 18:29
Re: 線形状(ベジェ曲線)の分割について教えてください
< Bezier Patch 補正  3 / 6 >


では、Patch 補正の内容はどのようになっているかというと、

  A = h0 + (B - h1)/2


もし、B も Patch 補正対象となっているなら、補正前の B 座標値を用いて A を求めます。






前回検証した5種類のサンプル形状に対する検証スクリプトに、この Patch 補正のロジックを組み込んでみたのが、次のスクリプトです。

このスクリプトでは自由曲面のハンドル長さの有無を has_handle( ) でチェックし、補正が必要と判断されたなら、上記の補正を加えて Patch を出力します。

このスクリプトでは DXF 経由で出力された Patch と、計算で求めた Patch がはれて一致しました。


ところで、もし何らかの計算に基づいて自由曲面を作成し、プラグインやスクリプトでて出力した場合、その数値計算過程には必ず数値の丸め誤差が生じるので、ハンドル長さ 0 になるはずのものが、誤差で 微小の長さのハンドルを持ってしまうことが考えられます。

次回は、そのような微小長さのハンドルを持つ自由曲面に対して、これまで考えてきた Patch 形成のロジックがそのまま適用できるか否かを調べます。

答えは No です。



#  サンプル自由曲面の Bezier Patch を出力する

#  計算で求める Bezier Patch に補正処理を組み込む



###  注意!!!  #####################################
#  DXF エクスポートの際のプロパティー設定をスクリプト側から操作できないので、
#  事前にマニュアルで設定しておく必要がある。

#  メニュー > ファイル > エクスポート > DXF を選択し、
#  DXF エクスポート設定の欄の 「ベジェ曲面 」 にチェックを入れ、
#  OK ボタンを押して、ファイル出力ダイアログが出たら キャンセル を選ぶ。
#  一旦この操作をしておけば、Shade を再起動後も設定は保たれたままとなる。



####    入力     #################################################
#  サンプル自由曲面
bLs = []

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((0, 1000, 0), None, (500, 1000, 0), None, None)
p4 = ((1000, 0, 0), (1000, 500, 0), None, (1000, 0, 500), None)
name = 'case - 1'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 2'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -250), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 3'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, None, None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -250), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 4'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, None)
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -250), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 5'
bLs.append((p1, p2, p3, p4, name))

############################################################



		
			
	
def write_sample_surface(bL) :	# block No. u1, v1 で与えられる区間のみを取り出した自由曲面を作成し、返す	
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	scene.begin_line(None, False)
	scene.append_point(bL[0][0], bL[0][1], bL[0][2], bL[0][3], bL[0][4])
	scene.append_point(bL[1][0], bL[1][1], bL[1][2], bL[1][3], bL[1][4])
	scene.end_line()
	
	scene.begin_line(None, False)
	scene.append_point(bL[2][0], bL[2][1], bL[2][2], bL[2][3], bL[2][4])
	scene.append_point(bL[3][0], bL[3][1], bL[3][2], bL[3][3], bL[3][4])
	scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	if bL[4] != '' :
		ob.name = bL[4]
	
	return ob 	#  作成した自由曲面を返す
		



def get_bezier_patch_by_DXF(ob) :	#  自由曲面 ob の bezier Patch を DXF file 経由で返す
	import os
	
	ob.activate()
	name = ob.name	#  テキストで書かれた DXFfile 読み取り後のの文字列分割において、オブジェクト名称にスペースが含まれると、
	ob.name = ''		#  そこでも分割されてしまうため、デフォルト名称に変更しておいて DXF file をエクスポートする。
	
	path = xshade.shade().scripts_path + '/temp999.dxf'	#  temp file の path
	scene.save_DXF(path)					# scripts folder 内に temp file として DXF file を出力
	ob.name = name	#  オブジェクト名称復活
	if not os.path.exists :
		return None
	
	f = open(path)
	t = f.read()
	f.close()
	os.remove(path)	#  temporary file 削除	
	
	s = t.split()
	if (s[11] !='16') or (s[19] != 'VERTEX') :
		print 'DXF エクスポートの設定として、「ベジェ曲面」にチェックを入れておいて下さい。'
		return None
	else :
		bP = []
		k = 23
		for i in range(4) :
			v = []
			for j in range(4) :
				v.append((float(s[k]), float(s[k + 4]), -float(s[k + 2])))	#  DXF file の座標は傾いている
				k = k + 12
			bP.append(v)
		return bP





def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])



def correct_patch(p1, p2, p3) :
	return (p1[0] + (p2[0] - p3[0])/2, p1[1] + (p2[1] - p3[1])/2, p1[2] + (p2[2] - p3[2])/2)



def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	
	###  ハンドル長さ判定  ###
	has_uH = (
		(ob2.control_point(0).has_out_handle, ob2.control_point(1).has_in_handle),
		(ob2.control_point(2).has_out_handle, ob2.control_point(3).has_in_handle))
	has_vH = (
		(ob2.control_point(0).has_lateral_out_handle, ob2.control_point(1).has_lateral_out_handle),
		(ob2.control_point(2).has_lateral_in_handle, ob2.control_point(3).has_lateral_in_handle))
	
	###  Bezier Patch 補正処理  ###
	bP11 = bP[1][1]
	bP12 = bP[1][2]
	bP21 = bP[2][1]
	bP22 = bP[2][2]
	
	if (not has_uH[0][0] and has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[1][0], bP21, bP[2][0])
	elif (has_uH[0][0] and not has_vH[0][0]) :
		bP[1][1] = correct_patch(bP[0][1], bP12, bP[0][2])
		
	if (not has_uH[0][1] and has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[1][3], bP22, bP[2][3])
	elif (has_uH[0][1] and not has_vH[0][1]) :
		bP[1][2] = correct_patch(bP[0][2], bP11, bP[0][1])
	
	if (not has_uH[1][0] and has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[2][0], bP11, bP[1][0])
	elif (has_uH[1][0] and not has_vH[1][0]) :
		bP[2][1] = correct_patch(bP[3][1], bP22, bP[3][2])
		
	if (not has_uH[1][1] and has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[2][3], bP12, bP[1][3])
	elif (has_uH[1][1] and not has_vH[1][1]) :
		bP[2][2] = correct_patch(bP[3][2], bP21, bP[3][1])

	ob2.remove()
	
	return bP
	


	
def write_bezier_patch(bP, name) :

	scene.begin_creating()
	scene.begin_surface_part()
	ob = scene.active_shape()
	ob.name = name
	ob.disclosed = False
	ob.cancel_transformation()	
	
	for i in range(4) :
		scene.begin_line(None, False)
		for j in range(4) :
			scene.append_point(bP[i][j], None, None, None, None)
		scene.end_line()
			
	scene.end_surface_part()
	scene.end_creating()
	
	ob.reset_transformation()
	
	return ob
	
	
	
def main() :

	scene.inhibit_update()
	
	scene.create_part('Patch 補正を考慮')
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	for bL in bLs :
		scene.create_part()
		ob = write_sample_surface(bL)
		bP = get_bezier_patch(ob)			#  temporary 自由曲面の bezier Patch を 計算により取得
		ob2 = write_bezier_patch(bP, 'Bezier Patch')	#  Bezier Patch を出力
		bP = get_bezier_patch_by_DXF(ob)		#  temporary 自由曲面の bezier Patch を DXF file 経由で取得
		if bP != None :
			ob3 = write_bezier_patch(bP, 'Bezier Patch by DXF')	#  Bezier Patch を出力
		scene.active_shape().dad.activate()

	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main()
投稿数: 189
投稿日時: 2011-10-01 18:18
Re: 線形状(ベジェ曲線)の分割について教えてください
< Bezier Patch 補正 2 / 6 >


どのようなケースでPatch 補正が行われているのかを以下のスクリプトで調べました。

このスクリプトでは5種類のサンプル形状を作成し、それぞれについて DXF 経由と計算による2種類の Bezier Patch を出力し、比較します。


検証の結果、Patch 補正の有無は次のようになっています。

  1)自由曲面の一つのアンカーポイントにおいて、交差する U, V 2つのハンドルのうち、
    一方の長さが 0 の場合、そのアンカーポイント周りに補正がなされる。

  2)2つのハンドルの両方共が長さがある場合、あるいは、両方とも長さが 0 の場合は
    補正は行われない。

補正の内容については次回、







#  サンプル自由曲面の Bezier Patch を出力する

#  Bezier Patch 補正 の適用条件を確認


###  注意!!!  #####################################
#  DXF エクスポートの際のプロパティー設定をスクリプト側から操作できないので、
#  事前にマニュアルで設定しておく必要がある。

#  メニュー > ファイル > エクスポート > DXF を選択し、
#  DXF エクスポート設定の欄の 「ベジェ曲面 」 にチェックを入れ、
#  OK ボタンを押して、ファイル出力ダイアログが出たら キャンセル を選ぶ。
#  一旦この操作をしておけば、Shade を再起動後も設定は保たれたままとなる。



####    入力     #################################################
#  サンプル自由曲面
bLs = []

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((0, 1000, 0), None, (500, 1000, 0), None, None)
p4 = ((1000, 0, 0), (1000, 500, 0), None, (1000, 0, 500), None)
name = 'case - 1'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), None, None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 2'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -250), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 3'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, None, None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, (500, 0, 1000))
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -250), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 4'
bLs.append((p1, p2, p3, p4, name))

p1= ((0, 1000, 0), None, (0, 1000, 500), None, None)
p2 = ((0, 0, 1000), (0, 500, 1000), None, None, None)
p3 = ((500, 1000, -500), None, (1000, 1000, -500), (500, 1000, -250), None)
p4 = ((1500, 0, -500), (1500, 500, -500), None, (1500, 0, 0), None)
name = 'case - 5'
bLs.append((p1, p2, p3, p4, name))

############################################################



		
			
	
def write_sample_surface(bL) :	# block No. u1, v1 で与えられる区間のみを取り出した自由曲面を作成し、返す	
	
	scene.begin_creating()
	scene.begin_surface_part()
	
	scene.begin_line(None, False)
	scene.append_point(bL[0][0], bL[0][1], bL[0][2], bL[0][3], bL[0][4])
	scene.append_point(bL[1][0], bL[1][1], bL[1][2], bL[1][3], bL[1][4])
	scene.end_line()
	
	scene.begin_line(None, False)
	scene.append_point(bL[2][0], bL[2][1], bL[2][2], bL[2][3], bL[2][4])
	scene.append_point(bL[3][0], bL[3][1], bL[3][2], bL[3][3], bL[3][4])
	scene.end_line()
	
	scene.end_surface_part()
	scene.end_creating()
	
	ob = scene.active_shape()
	ob.disclosed = False
	if bL[4] != '' :
		ob.name = bL[4]
	
	return ob 	#  作成した自由曲面を返す
		



def get_bezier_patch_by_DXF(ob) :	#  自由曲面 ob の bezier Patch を DXF file 経由で返す
	import os
	
	ob.activate()
	name = ob.name	#  テキストで書かれた DXFfile 読み取り後のの文字列分割において、オブジェクト名称にスペースが含まれると、
	ob.name = ''		#  そこでも分割されてしまうため、デフォルト名称に変更しておいて DXF file をエクスポートする。
	
	path = xshade.shade().scripts_path + '/temp999.dxf'	#  temp file の path
	scene.save_DXF(path)					# scripts folder 内に temp file として DXF file を出力
	ob.name = name	#  オブジェクト名称復活
	
	if not os.path.exists :
		return None
	
	f = open(path)
	t = f.read()
	f.close()
	os.remove(path)	#  temporary file 削除	
	
	s = t.split()
	if (s[11] !='16') or (s[19] != 'VERTEX') :
		print 'DXF エクスポートの設定として、「ベジェ曲面」にチェックを入れておいて下さい。'
		return None
	else :
		bP = []
		k = 23
		for i in range(4) :
			v = []
			for j in range(4) :
				v.append((float(s[k]), float(s[k + 4]), -float(s[k + 2])))	#  DXF file の座標は傾いている
				k = k + 12
			bP.append(v)
		return bP





def get_bezier_patch_1(p1, p2, p3) :
	return (p2[0] + p3[0] - p1[0], p2[1] + p3[1] - p1[1], p2[2] + p3[2] - p1[2])




def get_bezier_patch(ob) :		#  自由曲面 ob の bezier Patch を 計算し返す
	
	m = ob.local_to_world_matrix
	ob.copy_object(m)		#  座標値取得のため、local_to_world_matrix を乗じた自由曲面のコピーを作成
	ob2 = ob.bro
	ob2.reset_transformation()	#  自由曲面内の線形状の座標値取得のため、この処置を施しておく
	
	bP = []
	bP.append([
		ob2.control_point(0).position,
		ob2.control_point(0).out_handle,
		ob2.control_point(1).in_handle,
		ob2.control_point(1).position])
	bP.append([
		ob2.control_point(0).lateral_out_handle,
		None,
		None,
		ob2.control_point(1).lateral_out_handle])
	bP.append([
		ob2.control_point(2).lateral_in_handle,
		None,
		None,
		ob2.control_point(3).lateral_in_handle])
	bP.append([
		ob2.control_point(2).position,
		ob2.control_point(2).out_handle,
		ob2.control_point(3).in_handle,
		ob2.control_point(3).position])

	bP[1][1] = get_bezier_patch_1(bP[0][0], bP[0][1], bP[1][0])
	bP[1][2] = get_bezier_patch_1(bP[0][3], bP[0][2], bP[1][3])
	bP[2][1] = get_bezier_patch_1(bP[3][0], bP[3][1], bP[2][0])
	bP[2][2] = get_bezier_patch_1(bP[3][3], bP[3][2], bP[2][3])
	
	ob2.remove()
	
	return bP
	


	
def write_bezier_patch(bP, name) :

	scene.begin_creating()
	scene.begin_surface_part()
	ob = scene.active_shape()
	ob.name = name
	ob.disclosed = False
	ob.cancel_transformation()	
	
	for i in range(4) :
		scene.begin_line(None, False)
		for j in range(4) :
			scene.append_point(bP[i][j], None, None, None, None)
		scene.end_line()
			
	scene.end_surface_part()
	scene.end_creating()
	
	ob.reset_transformation()
	
	return ob
	
	
	
def main() :

	scene.inhibit_update()
	
	scene.create_part('ハンドルのない Bezier Patch')
	ob0 = scene.active_shape()
	ob0.cancel_transformation()
	
	for bL in bLs :
		scene.create_part()
		ob = write_sample_surface(bL)
		bP = get_bezier_patch(ob)			#  temporary 自由曲面の bezier Patch を 計算により取得
		ob2 = write_bezier_patch(bP, 'Bezier Patch')	#  Bezier Patch を出力
		bP = get_bezier_patch_by_DXF(ob)		#  temporary 自由曲面の bezier Patch を DXF file 経由で取得
		if bP != None :
			ob3 = write_bezier_patch(bP, 'Bezier Patch by DXF')	#  Bezier Patch を出力
		scene.active_shape().dad.activate()

	ob0.activate()
	scene.reset_all_transformation()
		
	scene.allow_update()
	scene.force_update()
	

	
	
scene = xshade.scene()
main()
スレッド表示 | 古いものから 投稿するには登録が必要です 前のスレッド | 次のスレッド | トップ

最近の投稿

フォーラム スレッド 返信 閲覧 最終投稿
Free Talk DNAの2重らせんの水素結合部位の作成 0 11757 2016-08-01 21:37 Benthos
Free Talk パート内の名前を一括返還 2 14266 2016-03-07 12:21 画像投稿機
Dev Forum イームズシェルチェアーの作成 2 14179 2015-11-25 14:44 CR7
Free Talk MOVファイルについて 17 35284 2014-12-29 17:14 momokuma
Dev Forum 2種類の液体アニメーションを作る方法 0 14259 2014-11-13 10:42 mejapan
Free Talk 面取りについて 0 13373 2014-11-08 15:18 MoonChild
Free Talk 丸太を結ぶ縄の作成について 1 19582 2014-09-18 22:33 kenslab
Free Talk パーティクルフィジックスのメタパーティクルについて 0 13835 2014-09-03 20:40 penta
Free Talk データの保存に関して 2 13664 2014-08-18 01:24 sierra
Free Talk Shade 3D ver14での、ポリゴンメッシュへの変換以上終了 1 14106 2014-04-23 12:04 MASA_