Shade3D 公式

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


#1

**16 / 05 / 17**    script を 四カ所 修正 / 追加

16 / 05 / 18    script を 一カ所 修正

16 / 06 / 05    script を 一カ所 修正

          修正 / 追加箇所については 末尾 27 - 4, 5 Script 修正 の項を参照

16 / 06 / 14, 17     script に v2_prev に係わる記述を追加 / 修正

          追加内容は 29 掃引の断面配置への追加事項 参照


関連記事 :
  • 08 Bezier line 等分割
  • 25 Custom 掃引
  • 26 滑らかな掃引( 1/2 )

27 - 1  曲線の曲がりに応じた分割( その1 )


曲線の曲がりに応じた分割を行おうとする時、**08 Bezier line 等分割** で述べた曲線の長さを N 分割する方法と同じ手法が考えられます。
  1. 当該 bezier 区間の 全曲率 を Simpson 法による数値積分で求め、分割数 N を定める

  2. Newton 法により、全曲率を N 分割する bezier parameter を求める


曲線が緩やかであればこの方法で可能ですが、曲がりが急激に変化する部分では Newton 法では解の収束が極めて困難になります。

全曲率は曲率を積分したものですが、曲がりが急激に変化する所では曲率の変化が大きすぎるため、単純な数値積分の Simpson 法では変化の大きさに追従しきれないのが原因です。




[ script 27 - 1 ] ではサンプル線形状の曲率を bezier parameter 0.25 ~ 0.37 の区間内で10等分してレポートします。

parameter のわずかな変化に対して、曲率には 1000倍 もの開きがあります。

もっと高度な数値積分を用いればこのような厳しい条件下でも収束が得られますが、それでも収束速度は遅くなり、ある程度の正確さを得ようとすれば全体としての演算コストは高めになってしまいます。


結果:

          [ table 26 - 1 ]


[ script 27 - 1 ]
import labo
from vec3 import *
from matrix import *

scene = xshade.scene()

bz = [vec3(-5000, 0, 0), vec3( -5000, 15000, 0), vec3(-10000, -2500, 0), vec3(15000, 2500, 0)]

scene.begin_creating()
scene.begin_line(None, False)
scene.append_point(bz[0], None, bz[1], None, None)
scene.append_point(bz[3], bz[2], None, None, None)
scene.end_line()
scene.end_creating()


bz = matrix(True, bz)
M = labo.bezier_m()

for i in range(11) :
	t = 0.25 + 0.12*float(i)/10
	Mt1 = labo.bezier_t_derivative(t)
	Mt2 = labo.bezier_t_second_derivative(t)
	v1 = Mt1*M*bz
	v2 = Mt2*M*bz
	v = v1*v2
		
	print t, '       ', v.abs()/(v1.abs())**3

#2

27 - 2  曲線の曲がりに応じた分割( その2 )


滑らかな掃引の作成では、**正確な N 分割を得る必要はなく**、曲がり具合に対して **ある程度の均一さで分割できればいい** と考えることにします。

そこで、全曲率を求めるのに数値積分は用いず、全曲率の可視化で述べた 単位球上の接線単位ベクトル先端軌跡の長さ として求める方法を用います。( [ fig 26 - 5 ] 参照 )

また N 分割 bezier parameter を求めるのに Newton 法 ではなく、単純な 線形補間 を用います。


          [ fig 26 - 5 ]



大まかな処理の流れは次のようになります。
**1)当該 bezier 区間の屈曲点( 変曲点 )を求め、その点で分割( 26 - 5 参照 )**
  • 屈曲点の求め方は既に 26 - 6 にて、関数 bezier_inflection_point( ) として作成済み( [ script 26 - 3 ] 参照 )

**2)1)で分割された bezier 区間を parameter 基準で微小区間に 最適分割**
  • 線形補間でもある程度の正確さを求めるなら、全曲率を求めるための分割はかなり細かくなる

  • そこで分割数を定数とはせず、曲がり具合に応じた最適分割を行い、演算コストの削減を図る

  • 平面曲線の場合、予め1)の変曲点での分割を行わなければ、精度を保てない。

  • 平面曲線の場合、変曲点での分割のみで正しい全曲率が求まるが、線形補間で N 分割点を求めるので、精度上ある程度の分割数は必要


**3)最適分割は次の要領による**
  • 最初に8分割し、各区間のおおよその全曲率を簡易的に求める

  • 簡易的な全曲率から分割数を求め、再分割

  • 以上を再帰関数で処理し、分割区間の全曲率が1.5度相当以下になるまで繰り返す


**4)当該 Bezier 区間の全曲率 = 微小各区間での全曲率の和 の大きさから分割数 N を決定**
  • N = ( 全曲率 - 0.0001) / しきい値 とする。

  • 0.0001 を減ずるのは、丸め誤差による分割数の不均一を避けるため


**5)3)で求めた全曲率と parameter リストから N 分割を与える bezier parameter を線形補間で求め 分割**

[ script 27 - 2 ] は選択された線形状を屈曲点で分割した後、指定するしきい角度に応じて分割する script です。

Shade 上で分割を行うのではなく、script 内で分割したものを出力します。


[ script 27 - 2 ]

import labo
from vec3 import *
from matrix import *

	
		
##   入力  ##########################################
level = 1./3			#  しきい角度 = level * π
####################################################




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

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

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

	A = Q[2]*Q[3]
	B = Q[1]*Q[3]
	C = Q[1]*Q[2]
	
	bb = labo.vec3_bounding_box_size(Q)
	bb.pop()
	bb2 = sorted(bb)
	for i in range(3) :
		if bb2[0] == bb[i] :
			a = A[i]
			b = B[i]
			c = C[i]
			t = []
			if a != 0 :
				d = b**2 - 4*a*c
				if d > 0 :
					tt = (-b - d**0.5)/(2*a)
					if tt > 0 and tt < 1 :
						t.append(tt)
					tt = (-b + d**0.5)/(2*a)
					if tt >= 0 and tt <= 1 :
						t.append(tt)
				elif d == 0 :
					tt = -b/(2*a)
					if tt > 0 and tt < 1 :
						t.append(tt)
			elif b != 0 :
				tt = -c/b
				if tt > 0 and tt < 1 :
					t.append(tt)
			
			if len(t) == 1 :
				v1 = labo.bezier_line_tangent(bz, 0)
				v2 = labo.bezier_line_tangent(bz, 1)
				v0 = labo.bezier_line_tangent(bz, t[0])
				if v0.dot(v1) >= 0.9999 or  v0.dot(v2) >= 0.9999 :
					t.pop()
					
			elif len(t) == 2 :
				t.sort()	
				v1 = labo.bezier_line_tangent(bz, 0)
				v2 = labo.bezier_line_tangent(bz, 1)
				v0 = labo.bezier_line_tangent(bz, t[1])
				if v0.dot(v2) >= 0.9999 :
					del t[1]
				v0 = labo.bezier_line_tangent(bz, t[0])
				if v0.dot(v1) >= 0.9999 :
					del t[0]
			return t
				
	
	
#  bezier_divide_by_total_curvature( bz as vec3 list, threshold as float) as float list
#
#		bz : 		bezier control point 座標を格納した vec3 list
#		threshold :	分割しきい値
#
#		bezier line bz を分割する bezier parameter を返す。
#		分割がない場合は、空 list が返される	

def bezier_divide_by_total_curvature(bz, threshold) :
	paramL = []								#  分割点の bezier parameter を格納する list
	
	v1 = labo.bezier_line_tangent(bz, 0)
	if v1.abs2() == 0 :
		return paramL						#  重点による長さのない区間
	
    #  bz に最適分割を行って cL[], pL[] にデータを格納する	
	pL = [0.]								#  各区間を分割する bezer parameter を格納する list
	cL = [0.]								#  各区間までの全曲率を格納する list
	v2 = labo.bezier_line_tangent(bz, 1)
	
	a = v1.dot(v2)
	if abs(a) > 0.99999994 :
		w = bz[3] - bz[0]
		if w.abs2() == 0 :
			return paramL					#  重点とみなす
		elif a < - 0.99999994 :
			w.norm()
			if abs(w.dot(v1)) > 0.99999994 :
				return paramL				#  重点と同じ扱いにする
	
	make_total_curvature_list(bz, 0., v1, 1., v2, pL, cL, 8)		#  最適分割を行って cL[], pL[] にデータを格納する 再帰関数
	theta = cL[len(cL) - 1]										#  bz の全曲率

	if theta > threshold :								#  bz の全曲率がしきい値を越えている
		n = int((theta - 0.00001)/threshold) + 1		#  分割数 ( 丸め誤差による分割数の不均一防止のため 0.00001 を減ずる )
		dc = theta/n				
		paramL = [0.]									#  分割点の bezier parameter を格納 (最初は 0. )
		j = 1
		for i in range(1, n) :
			c = i*dc
			while c > cL[j] :
				j += 1
			paramL.append(pL[j] - (pL[j] - pL[j - 1])*(cL[j] - c)/(cL[j ] - cL[j - 1]))	#  分割点の bezier parameter
		paramL.append(1.)																#  分割点の bezier parameter list の最後に 1 を追加
		
	return paramL
	

		
		
#  make_total_curvature_list(bz as vec3 list, t1 as float, v1 as vec3, t2 as float, v2 as vec3, pL as float list, cL as float list, n as int)
#
#		bz の parameter 区間 t1~t2 を最適サイズに分割して、cL[], pL[] にデータを格納する再帰関数
#
#		bz : 		bezier control point 座標を格納した vec3 list
#		t1, v1 :		計算対象区間の 始端 bezier parameter と その点での接線ベクトル	
#		t2, v2 :		計算対象区間の 終端 bezier parameter と その点での接線ベクトル
#		pL :			各区間を分割する bezer parameter を格納
#		cL :			分割された各区間までの全曲率を格納
#		n :			分割数 ( bezier parameter t1~t2 を n 分割する )

def make_total_curvature_list(bz, t1, v1, t2, v2, pL, cL, n) :
	from math import acos
	
	u = v1
	s1 = t1
	ds = (t2 - t1)/n								#  parameter  増分
	
	for i in range(1, n + 1) :						#  分割数 n で分割
		s2 = s1 + ds

		if i < n :
			v = labo.bezier_line_tangent(bz, s2)	# parameter t での接線ベクトル
		else :
			v = v2
			
		cosT = u.dot(v)
		
		if cosT >= 0.99966 :						#  0.99966 = cos(1.5度)
			cosT = max(-1, min(1, cosT))
			pL.append(s2)
			cL.append(cL[len(cL) - 1] + acos(cosT))	
		else :
			theta =  acos(cosT)
			m = max(2, 1 + int(theta*180/pi/1.5))					#  分割数
			make_total_curvature_list(bz, s1, u, s2, v, pL, cL, m)	#  再分割
			
		u = v
		s1 = s2
		
				
	
			

#  make_bezier_line(bZs as vec3 list, ob as shape, name as string)
#		bZs を線形状として出力
#		bZs = [bz, bz, bz, …]   bz = [vec3, vec3, vec3, vec3]
#		ob :		original 線形状
#		name :	出力する線形状の名称

def make_bezier_line(bZs, ob, name) :
	n = len(bZs)
	p = []
	inH = []
	outH = []

	closed = ob.closed
	if not closed :
		m = ob.number_of_control_points
		mx = matrix(ob.local_to_world_matrix)
		inH_ = mx*vec3(ob.control_point(0).in_handle)
		outH_ = mx*vec3(ob.control_point(m - 1).out_handle)
	
	if closed :
		inH.append(bZs[n - 1][2])
	else :
		inH.append(inH_)
	p.append(bZs[0][0])
	outH.append(bZs[0][1])
	
	for i in range(1, n) :
		inH.append(bZs[i - 1][2])
		p.append(bZs[i][0])
		outH.append(bZs[i][1])
		
	if not closed :
		inH.append(bZs[n - 1][2])
		p.append(bZs[n - 1][3])
		outH.append(outH_)
		
	labo.make_line(xshade, p, inH, outH, ob.closed, name)



	
	
from math import pi
	
scene = xshade.scene()
ob = scene.active_shape()


#  線形状全データ取得			bZs = [bz, bz, bz, …]   bz = [vec3, vec3, vec3, vec3] transposed
bZs0 = labo.get_bezier(xshade, -1)


scene.create_part()			#   出力線形状を格納する part
ob1 = scene.active_shape()


#  bZs1,  屈曲点にポイント追加
bZs1 = []
for bz in bZs0 :
	t = bezier_inflection_point(bz)
	n = len(t)
	if n == 0 :				#  屈曲点なし
		bZs1.append(bz)
	elif n == 1 :			#  屈曲点が一カ所あり
		bZs1.append(labo.bezier_line_subdivision(bz, 0, t[0]))
		bZs1.append(labo.bezier_line_subdivision(bz, t[0], 1))
	else :					#  屈曲点が二カ所あり
		bZs1.append(labo.bezier_line_subdivision(bz, 0, t[0]))
		bZs1.append(labo.bezier_line_subdivision(bz, t[0], t[1]))
		bZs1.append(labo.bezier_line_subdivision(bz, t[1], 1))
		
make_bezier_line(bZs1, ob, '屈曲点で分割')

	
#  bZs2,	 level に従って bezier 区間を分割
bZs2 = []
for bz in bZs1 :
	paramL = bezier_divide_by_total_curvature(bz, level*pi)
	n = len(paramL)
	if n == 0 : 
		bZs2.append(bz)
	else :
		for i in range(n - 1) :
			bZs2.append(labo.bezier_line_subdivision(bz, paramL[i], paramL[i + 1]))
		
make_bezier_line(bZs2, ob, 'level で分割')


ob1.activate()

#3

27 - 3  Script


以上から関数 **custom_sweep( )** に滑らかな掃引を作成する機能として [ script 27 - 2 ] を組み込んだ関数 **c_sweep( )** を作成、しきい角度などを指定するユーザーインターフェースを追加して [ script 27 - 3 ] として下に記します。
ただし、一カ所 [ script 27 - 2 ] とは異なっている所があります。

c_sweep( ) では屈曲点での接線ベクトルと始端あるいは終端での接線ベクトルがほとんど同じに近ければ、その屈曲点では分割しません。


[ 条件 27 - 1 ]

     


これはほとんど直線と言ってもいいような区間が発生するのを防ぐためです。

          [ fig 27 - 1 ]


**出力例**

          [ fig 27 - 2 ]



[ script 27 - 3 ]

16 / 05 / 17 script を一部修正

16 / 05 / 18 script を一部修正

16 / 06 / 05 script を一部修正

          修正箇所については 末尾 27 - 4, 5 Script 修正 の項を参照

16 / 06 / 14, 17 script に v2_prev に係わる記述を追加 / 修正

          追加内容は 29 掃引の断面配置への追加事項 参照


custom_sweep( ) からの変更や追加部分には行末に ############ を付しています。

import labo
from matrix import *
from quaternion import *

		

#  sweep_section_quaternion(bz as vec3 matrix ) as quaternion
#
#	bezier data bz から掃引時の断面形状の姿勢を与える quaternion を返す

def sweep_section_quaternion(bz) :
	v1 = labo.bezier_line_tangent(bz, 0)			#  outhandle 側接線ベクトル
	v2 = labo.bezier_line_tangent(bz, 1)			#  inhandle 側接線ベクトル

	if (v1.abs2() == 0) or (v2.abs2() == 0) :		#   重点
		return quaternion()
		
	elif abs(v1.dot(v2)) >= 0.99999994 :			#  区間 bz の outhandle と inhandle が 0度平行
		return quaternion()	
		
	else :
		t1 = 0
		Q = quaternion()
		for i in range(16) :						#  bezier parameter を16等分して細かく回転を16回繰り返す
			t2 = float(i + 1)/16
			v2 = labo.bezier_line_tangent(bz, t2)
			Q = quaternion().slerp(v1, v2)*Q
			t1 = t2
			v1 = v2
	
		return Q
			
		
			
				
#  sweep_plane_normal(p0 as vec3, p1 as vec3, p2 as vec3, p3 as vec3) as vec3
#
#	掃引中心線の bezier data bz (p0, p1, p2) から、掃引体交差方向の handle 長さを求めるための面法線を返す
#	基準となる anchor point に handle  が出ていることが前提
#
#		p0, p1, p2, p3 :		outhandle 側を求める場合は	bz[0], bz[1], bz[2], bz[3]
#							inhandle 側を求める場合は		bz[3], bz[2], bz[1], bz[0]

def sweep_plane_normal(p0, p1, p2, p3) :
	from math import pi, acos, sin, atan
	
	v1 = p1 - p0					#  v1.abs() != 0 が前提
	v2 = p2 - p1
	v1.norm()
	v2.norm()
	
	if v2.abs2() == 0 :
		v3 = p3 - p1
		v3.norm()
		if v3.abs2() == 0 :
			return v1
		else :
			vn = v1 + v3
			vn.norm()
			vp = v1 + vn
			vp.norm()
			return vp
		
	else :
		a = v1.dot(v2)
		if abs(a) >= 0.99999994 :
			return v1
		
		elif a >= 0 :
			vn = v1 + v2
			vn.norm()
			
			sinA = a
			cosA = (1 - sinA**2)**0.5
			beta = (pi - acos(a))/2
			r = (1 - sinA)/3
			
			axis = v1*v2
			axis.norm()
			v = v1 - v2
			v.norm()
			u = v/sin(beta) + v1*r/cosA
			u.norm()
			phi = acos(u.dot(v))

		else :
			axis = v1*v2
			phi = atan(4./3) - pi/4
			Q = quaternion().rotation(axis, pi/4)
			vn = Q*v1

		v3 = p3 - p2
		v3.norm()
		if v3.abs2() == 0 :
			return vn
		else :
			d = v1.dot(v2)
			w = -(v1-2*d*v2)
			w.norm()
			e = v3.dot(w)
			Q = quaternion().rotation(axis, e*phi)
			vp = Q*vn
			vp.norm()
			return vp

	
	
	
#  sweep_handle(w as vec3, u as vec3,  Q as vec3, r as float, p as vec3, vp as vec3, corner as bool) as vec3
#
#	handle 座標を返す

def sweep_handle(w, u, Q, r, p, vp, corner) :
	#     L :  点 w を通り、u 方向に進む直線
	#	PL0 : 点 Q を通り面法線 u なる面
	#	PL1 : 点 p を通り面法線 vp なる面
	#	精度確保のため r = 1 とはせず、r = bounding box size / 16 としてある

	a = u.dot(w - Q)			#  点 w と面 PL0 との符号付き距離		
	d1 = vp.dot(w - p)			#  点 w と 面 PL1 との符号付き距離
	d2 = vp.dot(w + r*u - p)	#  点 w + r*u と 面 PL1 との符号付き距離
	d = r*d1/(d1 - d2)			#  点 w と、直線 L と 面 PL1 との交点までの 符号付き距離
	if not corner :
		d = d + a				#  handle の符号つき長さ
	if d < 0 :					#  食い込みが発生する
		d = -d/12
	return w + d*u
		

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

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

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

	A = Q[2]*Q[3]
	B = Q[1]*Q[3]
	C = Q[1]*Q[2]
	
	bb = labo.vec3_bounding_box_size(Q)
	bb.pop()
	bb2 = sorted(bb)
	for i in range(3) :
		if bb2[0] == bb[i] :
			a = A[i]
			b = B[i]
			c = C[i]
			t = []
			if a != 0 :
				d = b**2 - 4*a*c
				if d > 0 :
					tt = (-b - d**0.5)/(2*a)
					if tt > 0 and tt < 1 :
						t.append(tt)
					tt = (-b + d**0.5)/(2*a)
					if tt >= 0 and tt <= 1 :
						t.append(tt)
				elif d == 0 :
					tt = -b/(2*a)
					if tt > 0 and tt < 1 :
						t.append(tt)
			elif b != 0 :
				tt = -c/b
				if tt > 0 and tt < 1 :
					t.append(tt)
			
			if len(t) == 1 :
				v1 = labo.bezier_line_tangent(bz, 0)
				v2 = labo.bezier_line_tangent(bz, 1)
				v0 = labo.bezier_line_tangent(bz, t[0])
				if v0.dot(v1) >= 0.9999 or  v0.dot(v2) >= 0.9999 :	#  屈曲点での接線ベクトルが、始端 / 終端での接線ベクトルと
					t.pop()											#  ほとんど同じであれば、その屈曲点では分割しない
					
			elif len(t) == 2 :
				t.sort()	
				v1 = labo.bezier_line_tangent(bz, 0)
				v2 = labo.bezier_line_tangent(bz, 1)
				v0 = labo.bezier_line_tangent(bz, t[1])
				if v0.dot(v2) >= 0.9999 :							#  屈曲点での接線ベクトルが、終端での接線ベクトルと
					del t[1]										#  ほとんど同じであれば、その屈曲点では分割しない
				v0 = labo.bezier_line_tangent(bz, t[0])
				if v0.dot(v1) >= 0.9999 :							#  屈曲点での接線ベクトルが、始端での接線ベクトルと
					del t[0]										#  ほとんど同じであれば、その屈曲点では分割しない
			return t
				
	

#  bezier_divide_by_total_curvature( bz as vec3 list, threshold as float) as float list
#
#		bz : 		bezier control point 座標を格納した vec3 list
#		threshold :	分割しきい値
#
#		bezier line bz を分割する bezier parameter を返す。
#		分割がない場合は、空 list が返される	

def bezier_divide_by_total_curvature(bz, threshold) :
	paramL = []								#  分割点の bezier parameter を格納する list
	
	v1 = labo.bezier_line_tangent(bz, 0)
	if v1.abs2() == 0 :
		return paramL						#  重点による長さのない区間
	
 	#  bz に最適分割を行って cL[], pL[] にデータを格納する	
	pL = [0.]								#  各区間を分割する bezer parameter を格納する list
	cL = [0.]								#  各区間までの全曲率を格納する list
	v2 = labo.bezier_line_tangent(bz, 1)
	
	a = v1.dot(v2)
	if abs(a) > 0.99999994 :
		w = bz[3] - bz[0]
		if w.abs2() == 0 :
			return paramL					#  重点とみなす
		elif a < - 0.99999994 :
			w.norm()
			if abs(w.dot(v1)) > 0.99999994 :
				return paramL				#  重点と同じ扱いにする
	
	make_total_curvature_list(bz, 0., v1, 1., v2, pL, cL, 8)		#  最適分割を行って cL[], pL[] にデータを格納する 再帰関数
	theta = cL[len(cL) - 1]										#  bz の全曲率
	
	if theta > threshold :								#  bz の全曲率がしきい値を越えている
		n = int((theta - 0.00001)/threshold) + 1		#  分割数 ( 丸め誤差による分割数の不均一防止のため 0.00001 を減ずる )
		dc = theta/n				
		paramL = [0.]									#  分割点の bezier parameter を格納 (最初は 0. )
		j = 1
		for i in range(1, n) :
			c = i*dc
			while c > cL[j] :
				j += 1
			paramL.append(pL[j] - (pL[j] - pL[j - 1])*(cL[j] - c)/(cL[j ] - cL[j - 1]))	#  分割点の bezier parameter
		paramL.append(1.)															#  分割点の bezier parameter list の最後に 1 を追加
		
	return paramL
	

		
		
#  make_total_curvature_list(bz as vec3 list, t1 as float, v1 as vec3, t2 as float, v2 as vec3, pL as float list, cL as float list, n as int)
#
#		bz の parameter 区間 t1~t2 を最適サイズに分割して、cL[], pL[] にデータを格納する再帰関数
#
#		bz : 		bezier control point 座標を格納した vec3 list
#		t1, v1 :		計算対象区間の 始端 bezier parameter と その点での接線ベクトル	
#		t2, v2 :		計算対象区間の 終端 bezier parameter と その点での接線ベクトル
#		pL :			各区間を分割する bezer parameter を格納
#		cL :			分割された各区間までの全曲率を格納
#		n :			分割数 ( bezier parameter t1~t2 を n 分割する )

def make_total_curvature_list(bz, t1, v1, t2, v2, pL, cL, n) :
	from math import acos
	
	u = v1
	s1 = t1
	ds = (t2 - t1)/n								#  parameter  増分
	
	for i in range(1, n + 1) :						#  分割数 n で分割
		s2 = s1 + ds

		if i < n :
			v = labo.bezier_line_tangent(bz, s2)	# parameter t での接線ベクトル
		else :
			v = v2
			
		cosT = max(-1, min(1, u.dot(v)))
		
		if cosT >= 0.99966 :						#  0.99966 = cos(1.5度)
			pL.append(s2)
			cL.append(cL[len(cL) - 1] + acos(cosT))	
		else :
			theta =  acos(cosT)
			m = max(2, 1 + int(theta*180/pi/1.5))					#  分割数
			make_total_curvature_list(bz, s1, u, s2, v, pL, cL, m)	#  再分割
			
		u = v
		s1 = s2
		
		
		

#  c_sweep(ob1 as shape, ob2 as shape, threshold as float, limitter1 as float = 10., limitter2 as float = 3.) as shape
#
#	掃引体を Shade 上に作成し、その shape を返す
#		ob1 :	掃引中心線形状
#		ob2 :	掃引断面形状
#		threshold:	分割しきい値
#		limitter1 :	直線コーナーの scale liitter	
#		limitter2 :	曲線コーナーの scale liitter	

def c_sweep(ob1, ob2, threshold, limitter1 = 10., limitter2 = 3.) :
	if (ob1 == None) or (ob2 == None) :
		return
	if (not isinstance(ob1, xshade.line)) or (not isinstance(ob2, xshade.line)) :
		return
		
	#  check							############
	if threshold < 0 :					############
		return							############
	if limitter1 < 1 or limitter2 < 1 :	############
		return							############
	

	#  掃引中心線データ取得 nop, closed, bZ
	ob1.activate()
	nop = ob1.number_of_control_points		#  ポイント数
	if nop < 2 :
		return
	closed = ob1.closed						#  開閉情報
	bZs = labo.get_bezier(xshade, -1)		#  全区間の bezier data [matrix, matrix, … ]
	
	
	# mx1, mx2, mx22, ps, pe
	mx1 = matrix(ob2.world_to_local_matrix)
	mx2 = matrix(ob2.local_to_world_matrix)
	if not closed :
		mx22 = matrix(ob1.local_to_world_matrix)
		ps = vec3(ob1.control_point(0).in_handle)*mx22			#  掃引体の交差方向 handle 作成時に使用	############
		pe = vec3(ob1.control_point(nop - 1).out_handle)*mx22	#  掃引体の交差方向 handle 作成時に使用	############
		
	#  r							
	bb = labo.vec3_bounding_box_size(bZs)														############
	r = bb[3]/16												#  sweep_handle( ) の引数に使用  	############	
	
		
	#  掃引中心線 の handle 情報取得
	has_in_handle = []			#  掃引中心線の inhandle の有無のリスト
	has_out_handle = []			#  掃引中心線の outhandle の有無のリスト
	linked = []					#  掃引中心線の hanndle link リスト
	for i in range(nop) :
		has_in_handle.append(ob1.control_point(i).has_in_handle)
		has_out_handle.append(ob1.control_point(i).has_out_handle)
		linked.append(ob1.control_point(i).linked)

		
	#  bZs を分割		###  追加  ####################################
	if threshold > 0 :

		#  屈曲点で分割
		bZs1 = []
		insert_list = []
		k = 0
		for bz in bZs :
			k += 1
			t = bezier_inflection_point(bz)
			n = len(t)
			if n == 0 :				#  屈曲点なし
				bZs1.append(bz)
			elif n == 1 :				#  屈曲点が一カ所あり
				bZs1.append(labo.bezier_line_subdivision(bz, 0, t[0]))
				bZs1.append(labo.bezier_line_subdivision(bz, t[0], 1))
				insert_list.append(k)
				nop += 1
			else :					#  屈曲点が二カ所あり
				bZs1.append(labo.bezier_line_subdivision(bz, 0, t[0]))
				bZs1.append(labo.bezier_line_subdivision(bz, t[0], t[1]))
				bZs1.append(labo.bezier_line_subdivision(bz, t[1], 1))
				insert_list.append(k)
				insert_list.append(k)
				nop += 2

		if len(insert_list) > 0 :
			insert_list.reverse()
			for k in insert_list :
				has_in_handle.insert(k, True)
				has_out_handle.insert(k, True)
				linked.insert(k, True)

		
		#  threshold に従って分割
		bZs2 = []
		insert_list = []
		k = 0
		for bz in bZs1 :
			k += 1
			paramL = bezier_divide_by_total_curvature(bz, threshold)
			n = len(paramL)
			if n == 0 : 
				bZs2.append(bz)
			else :
				for i in range(n - 1) :
					bZs2.append(labo.bezier_line_subdivision(bz, paramL[i], paramL[i + 1]))
					insert_list.append(k)
					nop += 1
				insert_list.pop()
				nop += -1
					
		if len(insert_list) > 0 :
			insert_list.reverse()
			for k in insert_list :
				has_in_handle.insert(k, True)
				has_out_handle.insert(k, True)
				linked.insert(k, True)
				
		bZs = bZs2			###   追加 ここまで #########################################
											
											
	
	#  ob0 : 掃引体を格納する自由曲面
	ob2.activate()
	scene.create_surface_part(None)
	ob0 = scene.active_shape()				#   掃引体を格納する自由曲面
	ob0.surface_closed = closed
	ob2.place_child(1)

	
	#  掃引断面を配置
	p1 = bZs[0][0]
	Mt1 = matrix().translate(-p1[0], -p1[1], -p1[2])
	Q1 = quaternion()
	n = nop - 1
	v2_prev = None
	cornerL = [False]
		
	for i in range(n) :
		p2 = bZs[i][3]
		v2 = labo.bezier_line_tangent(bZs[i], 1)
		Mt2 = matrix().translate(p2[0], p2[1], p2[2])
		Q2 = sweep_section_quaternion(bZs[i])
		Q2 = Q2*Q1
		Mr = Q2.matrix()
		corner = False
		
		if i == n - 1 and not closed :
			M = mx2*Mt1*Mr*Mt2*mx1
			
		else :
			v3 = labo.bezier_line_tangent(bZs[i + 1], 0)
			if v3.abs2() == 0 :							#  重点
				M = mx2*Mt1*Mr*Mt2*mx1
				if v2_prev == None :					#  最初の重点
					v2_prev = v2
				
			elif v2.abs2() == 0 :							#  最後の重点
				if v2_prev == None or v2_prev.abs2() == 0 :	#  始点が重点から始まっている
					M = mx2*Mt1*Mr*Mt2*mx1
				else :
					Qc = quaternion().slerp(v2_prev, v3, 1, False)	
					if Qc != None :	
						Q2 = Qc*Q2
						Mc = Qc.matrix()
						M = mx2*Mt1*Mr*Mc*Mt2*mx1
					else :
						M = mx2*Mt1*Mr*Mt2*mx1
				v2_prev = None
				
			else :
				if abs(v2.dot(v3)) >= 0.99999994 :	#  掃引中心線の point が corner point  ではない
					M = mx2*Mt1*Mr*Mt2*mx1
				
				else :								#  掃引中心線の point が corner point 
					w = v2 + v3
					w.norm()
					corner = True
					v1 = labo.bezier_line_tangent(bZs[i], 0)
					v4 = labo.bezier_line_tangent(bZs[i + 1], 1)

					if v1.dot(v2) + v3.dot(v4) > 1.98 :	#  直線コーナー
						limitter = limitter1
					else :								#  曲線コーナー
						limitter = limitter2
					Q3 = quaternion().slerp(v2, w)
					Mq = Q3.matrix()
					
					v = v2 - v3
					cosT = w.dot(v2)
					Q4 = quaternion().slerp(v, vec3(0, 1, 0))
					Mqs1 = Q4.matrix()
					Q5 = Q4.inverse()
					Mqs2 = Q5.matrix()
					Ms = matrix().scale(1, min(limitter, 1/cosT), 1)	
					M = mx2*Mt1*Mr*Mq*Mqs1*Ms*Mqs2*Mt2*mx1

					Q6 = quaternion().slerp(v2, v3)
					Q2 = Q6*Q2						#   次の point での掃引断面配置に必要

	
		ob2.copy_object(M)
		ob = ob2.bro
		ob.place_brother(i)
	
		Q1 = Q2
		cornerL.append(corner)

	
	#  掃引体の交差方向 handle をセット
#	bb = labo.vec3_bounding_box_size(bZ)
#	r = bb[3]/16												#  sweep_handle( ) の引数に使用
	
	n = ob2.number_of_control_points
	ob0.switch()
	ob = ob0.son
	for i in range(n) :
		ob = ob.bro												#  掃引体内部の交差方向線形状
		
		for j in range(nop) :
			w = vec3(ob.control_point(j).position)*mx2			#  交差方向線形状の point 座標
			
			if not has_in_handle[j] :
				inH = w											#  inhandle 座標
			else :
				if (j == 0) and not closed :
					inH = w + ps - bZs[0][0]					#  inhandle 座標		############
				else :
					if j == 0 :
						k = nop - 1
					else :
						k = j - 1
					vp = sweep_plane_normal(bZs[k][3], bZs[k][2], bZs[k][1], bZs[k][0])
					u = bZs[k][2] - bZs[k][3] 
					u.norm()
					p = bZs[k][2]
					Q = bZs[k][3]
					inH = sweep_handle(w, u, Q, r, p, vp, cornerL[j])		#  inhandle 座標
					
	
			if not has_out_handle[j] :
				outH = w										#  outhandle 座標
			else :
				if (j == nop - 1) and not closed :
					outH = w + pe - bZs[nop - 2][3]				#  outhandle 座標	############
				else :
					vp = sweep_plane_normal(bZs[j][0], bZs[j][1], bZs[j][2], bZs[j][3])
					u = bZs[j][1] - bZs[j][0] 
					u.norm()
					p = bZs[j][1]
					Q = bZs[j][0]
					outH = sweep_handle(w, u, Q, r, p, vp, cornerL[j])		#  outhandle 座標
						

			ob.control_point(j).in_handle = inH*mx1
			ob.control_point(j).out_handle = outH*mx1					
			ob.control_point(j).linked = linked[j]
	
		
	ob0.switch()
	ob0.activate()
	ob0.disclosed = False
	
	return ob0
	

	
	
#  Shade  上で 掃引中心線 , 掃引断面 の順で選択して実行
scene = xshade.scene()

#  dialog uuid は Online GUID Generator により生成		https://www.guidgenerator.com
dialog = xshade.create_dialog_with_uuid('1f8f4d7b-4634-4219-8567-431224ecf60d')

dialog.begin_tab_group('分割')
level_ = dialog.append_selection('level/分割せず/90度/60度/45度/30度/22.5度/15度')
dialog.end_tab_group()

dialog.begin_tab_group('limitter')
limitter1_ = dialog.append_float('直線コーナー', '( defalt 10 )')
limitter2_ = dialog.append_float('曲線コーナー', '( defalt 3 )')
dialog.end_tab_group()

dialog.set_value(4, 10)
dialog.set_value(5, 3)

if dialog.ask('滑らかな掃引'):
	level = dialog.get_value(level_)
	from math import pi
	if level == 0 :
		threshold = 0.
	elif level == 1 :	
		threshold = pi/2	
	elif level == 2 :		
		threshold = pi/3
	elif level == 3 :	
		threshold = pi/4
	elif level == 4 :	
		threshold = pi/6
	elif level == 5 :	
		threshold = pi/8
	elif level == 6 :	
		threshold = pi/12
		
	limitter1 = dialog.get_value(limitter1_)
	limitter2 = dialog.get_value(limitter2_)

	scene.exit_modify_mode()
	[ob1, ob2] = scene.active_shapes
	ob2.copy_object(None)
	ob3 = ob2.bro
	ob4 = c_sweep(ob1, ob3, threshold, limitter1, limitter2)
	if ob4 != None :
		ob4.activate()    
	else :
		scene.active_shapes = [ob1, ob2]

#4

27 - 4  Script 修正


16 / 07 / 10 削除

#5

27 - 5  Script 修正( 16 / 06 / 05 )


[ script 27 - 3 ] に対して、下記一カ所 追加 / 変更 を行いました
**def bezier_divide_by_total_curvature(bz, threshold) :**
**旧**
#  bz に最適分割を行って cL[], pL[] にデータを格納する	
pL = [0.]								#  各区間を分割する bezer parameter を格納する list
cL = [0.]								#  各区間までの全曲率を格納する list
v2 = labo.bezier_line_tangent(bz, 1)
	
w = bz[3] - bz[0]										
if (w.abs2() == 0) and(abs(v1.dot(v2)) > 0.99999994) :		
	return paramL						#  重点とみなす
	
make_total_curvature_list(bz, 0., v1, 1., v2, pL, cL, 8)	#  最適分割を行って cL[], pL[] にデータを格納する 再帰関数
theta = cL[len(cL) - 1]										#  bz の全曲率

**追加 / 変更**
#  bz に最適分割を行って cL[], pL[] にデータを格納する	
	pL = [0.]								#  各区間を分割する bezer parameter を格納する list
	cL = [0.]								#  各区間までの全曲率を格納する list
	v2 = labo.bezier_line_tangent(bz, 1)
	
	a = v1.dot(v2)
	if abs(a) > 0.99999994 :
		w = bz[3] - bz[0]
		if w.abs2() == 0 :
			return paramL					#  重点とみなす
		elif a < - 0.99999994 :
			w.norm()
			if abs(w.dot(v1)) > 0.99999994 :
				return paramL				#  重点と同じ扱いにする
	
	make_total_curvature_list(bz, 0., v1, 1., v2, pL, cL, 8)		#  最適分割を行って cL[], pL[] にデータを格納する 再帰関数
	theta = cL[len(cL) - 1]										#  bz の全曲率