接線・法線ベクトル
vlax-curve-getFirstDeriv 関数や vlax-curve-getSecondDeriv 関数が精度や境界値での動作が不安定といった少し残念なものでしたので、使える関数を使用して近似で接線ベクトルと法線ベクトルを計算する関数を作ってみます。
単位接線ベクトルと法線ベクトルを求める関数は、次の curve-tangent 関数と curve-normal 関数になります。カーブ上の位置は【パラメータ】で指定します。始点からの長さで指定するバージョン curve-tangent:dist 関数と curve-normal:dist 関数も使えます。法線ベクトルは、調べた点でカーブに接する円の中心を指すベクトルとなります。このベクトルの大きさは「曲率半径」と呼ばれ、その逆数が「曲率」です。直線の場合の法線ベクトルは理論的には無限大の大きさを持ちますが、ここでは nil を返します。
やはり制限があります。微小区間でのカーブ上の点の変化を調べて接線ベクトルの近似としています。あまり考えられないことですが、計算誤差に対応して微小区間がカーブの長さより大きくなった場合は nil を返します。カーブが連続して微分可能という前提ですので、ライトウェイトポリラインの折れ点付近では論理的には正しくなりません。2D/3D ポリラインの場合は getPointAt～ 系関数の曲線の解像度が低いため、単位接線ベクトルの精度は低く法線ベクトルは取得できません。
;;;=====================================
;;; support function in curve function
;;;=====================================
(defun curve:vector:size (vector1)
(sqrt (vector:InnerProduct vector1 vector1))
)
(defun curve:vector:unit (vector1 / size)
(if vector1
(progn (setq size (curve:vector:size vector1))
(mapcar (function (lambda (e) (/ e size))) vector1)
)
)
)
(defun curve:point:dist (ename dist / clength)
(setq clength (vlax-curve-getDistAtParam
ename
(vlax-curve-getEndParam ename)
)
)
(if (equal dist clength (* clength (expt 10.0 -10)))
(vlax-curve-getEndPoint ename)
(vlax-curve-getPointAtDist ename dist)
)
)
(defun curve:differential:set-p0/p1-statement
(dist clength delta-s/2 tolerance)
(cond ((equal 0.0 dist tolerance)
(if (vlax-curve-isClosed ename)
(setq p0-statement '(- clength delta-s/2)
p1-statement 'delta-s/2
)
(setq p0-statement '0.0
p1-statement '(* 2.0 delta-s/2)
)
)
)
((equal clength dist tolerance)
(if (vlax-curve-isClosed ename)
(setq p0-statement '(- clength delta-s/2)
p1-statement 'delta-s/2
)
(setq p0-statement '(- clength (* 2.0 delta-s/2))
p1-statement 'clength
)
)
)
(T
(setq p0-statement '(- dist delta-s/2)
p1-statement '(+ dist delta-s/2)
)
)
)
)
;;;=====================================
;;; (curve-tangent ename param)
;;; (curve-tangent:dist ename dist)
;;;=====================================
(defun curve-tangent:dist:approximation (delta-s/2 / result)
(if (< delta-s/2 (/ clength 2.0))
(if (vl-every 'zerop
(setq result (vector:sub
(curve:point:dist ename (eval p1-statement))
(curve:point:dist ename (eval p0-statement))
)
)
)
(curve-tangent:dist:approximation (* delta-s/2 10.0))
result
)
)
)
(defun curve-tangent:dist
(ename dist / clength tolerance delta-s/2 p0-statement p1-statement)
(setq clength (vlax-curve-getDistAtParam
ename
(vlax-curve-getEndParam ename)
)
tolerance (* clength (expt 10.0 -10))
delta-s/2 (* clength (expt 10.1 -10))
)
(curve:differential:set-p0/p1-statement dist clength delta-s/2 tolerance)
(curve:vector:unit (curve-tangent:dist:approximation delta-s/2))
)
(defun curve-tangent (ename param / tangent)
(curve-tangent:dist ename (vlax-curve-getDistAtParam ename param))
)
;;;=====================================
;;; (curve-normal ename param)
;;; (curve-normal:dist ename dist)
;;;=====================================
(defun curve-normal:dist:approximation
(delta-s/2 / normal curvature result)
(if (< delta-s/2 (/ clength 20.0))
(if (equal (setq curvature (curve:vector:size
(setq normal (vector:div
(vector:sub
(curve-tangent:dist
ename
(eval p1-statement)
)
(curve-tangent:dist
ename
(eval p0-statement)
)
)
(* 2.0 delta-s/2)
)
)
)
)
0.0
(expt 10.0 -10)
)
(curve-normal:dist:approximation (* delta-s/2 10.0))
(vector:scale normal (/ 1.0 curvature curvature))
)
)
)
(defun curve-normal:dist
(ename dist / clength tolerance delta-s/2 p0-statement p1-statement)
(setq clength (vlax-curve-getDistAtParam
ename
(vlax-curve-getEndParam ename)
)
tolerance (* clength (expt 10.0 -4))
delta-s/2 (* clength (expt 10.1 -4))
)
(curve:differential:set-p0/p1-statement dist clength delta-s/2 tolerance)
(curve-normal:dist:approximation delta-s/2)
)
(defun curve-normal (ename param)
(curve-normal:dist ename (vlax-curve-getDistAtParam ename param))
)
この接線と法線を求める関数を使ったプログラムを三つほど書いてみます。
次の一つ目の例は、図形を選択すると、長さを10等分した位置の単位接線ベクトルと法線ベクトルをグリーンとマジェンタで描画します。
(defun curve:point:dist (ename dist / clength)
(setq clength (vlax-curve-getDistAtParam
ename
(vlax-curve-getEndParam ename)
)
)
(if (equal dist clength (* clength (expt 10.0 -10)))
(vlax-curve-getEndPoint ename)
(vlax-curve-getPointAtDist ename dist)
)
)
(defun drawVector (origin v color /) ;_ origin,v(UCS)
(grvecs (list color origin (vector:add origin v)))
)
(defun ParamList:sub (index)
(if (<= index n)
(cons (+ start (* step index)) (ParamList:sub (1+ index)))
)
)
(defun ParamList (start end n / step)
(setq step (/ (float (- end start)) n))
(ParamList:sub 0)
)
(while
(setq ename (car (entsel "\n図形を選択（選択なしで終了）[最後(L)]：")))
(if (setq endParam (vlax-curve-getEndParam ename))
(foreach dist (ParamList 0.0
(vlax-curve-getDistAtParam ename endParam)
10.0
)
(if (setq normal (curve-normal:dist ename dist))
(drawVector
(trans (curve:point:dist ename dist) acWorld acUCS)
(trans normal acWorld acUCS T)
6
)
)
(drawVector
(trans (curve:point:dist ename dist) acWorld acUCS)
(trans (curve-tangent:dist ename dist) acWorld acUCS T)
3
)
)
)
)
次の二つ目の例は、図形を選択すると、選択した点で接する円を描きます。
(defun range (s e /)
(if (<= s e)
(cons s (range (1+ s) e))
)
)
(defun PolarPolygon (center radius offset numberOfVertex / step)
(setq step (/ (* 2 PI) numberOfVertex))
(mapcar (function
(lambda (index) (polar center (+ offset (* step index)) radius))
)
(range 0 (1- numberOfVertex))
)
)
(defun grvecs:DrawData:sub (plist / pt1 pt2)
(if (and plist
(setq pt1 (car plist))
(setq pt2 (if (cadr plist)
(cadr plist)
(if closed
startPoint
nil
)
)
)
)
(cons (list color pt1 pt2) (grvecs:DrawData:sub (cdr plist)))
nil
)
)
(defun grvecs:DrawData (plist color closed / startPoint)
(setq startPoint (car plist))
(apply 'append (grvecs:DrawData:sub plist))
)
(while (setq selInfo (entsel "\n図形を選択（選択なしで終了）[最後(L)]："))
(setq ename (car selInfo)
point (osnap (cadr selInfo) "NEA")
param (vlax-curve-getParamAtPoint ename (trans point acUCS acWorld))
)
(if (setq normal (curve-normal ename param))
(progn
(setq normal (trans normal acWorld acUCS T)
amatrix (make-transformation-matrix
(vector:add point normal)
(trans (curve-tangent ename param) acWorld acUCS T)
normal
)
)
(grvecs
(grvecs:DrawData
(mapcar (function (lambda (point) (transform amatrix point nil)))
(PolarPolygon '(0.0 0.0 0.0) (vector:size normal) 0.0 32)
)
6
T
)
)
)
)
)
最後の三つ目の例は、選択した図形の曲率を画面に描きます。図形サイズと画面表示によって、コードの中間辺りにある変数 scale 値を調整してください。
(defun curve:point:dist (ename dist / clength)
(setq clength (vlax-curve-getDistAtParam
ename
(vlax-curve-getEndParam ename)
)
)
(if (equal dist clength (* clength (expt 10.0 -10)))
(vlax-curve-getEndPoint ename)
(vlax-curve-getPointAtDist ename dist)
)
)
(defun grvecs:DrawData:sub (plist / pt1 pt2)
(if (and plist
(setq pt1 (car plist))
(setq pt2 (if (cadr plist)
(cadr plist)
(if closed
startPoint
nil
)
)
)
)
(cons (list color pt1 pt2) (grvecs:DrawData:sub (cdr plist)))
nil
)
)
(defun grvecs:DrawData (plist color closed / startPoint)
(setq startPoint (car plist))
(apply 'append (grvecs:DrawData:sub plist))
)
(defun ParamList:sub (index)
(if (<= index n)
(cons (+ start (* step index)) (ParamList:sub (1+ index)))
)
)
(defun ParamList (start end n / step)
(setq step (/ (float (- end start)) n))
(ParamList:sub 0)
)
(setq scale 1000.0) ;_ change scale by view size
(while
(setq ename (car (entsel "\n図形を選択（選択なしで終了）[最後(L)]：")))
(if (setq endParam (vlax-curve-getEndParam ename))
(grvecs
(grvecs:DrawData
(mapcar
(function
(lambda (dist / normal)
(if (setq normal (curve-normal:dist ename dist))
(progn (setq normal (trans normal acWorld acUCS T))
(vector:add
(trans (curve:point:dist ename dist) acWorld acUCS)
(vector:scale
(vector:invert (vector:unit normal))
(/ scale (vector:size normal))
)
)
)
(trans (curve:point:dist ename dist) acWorld acUCS)
)
)
)
(ParamList 0.0 (vlax-curve-getDistAtParam ename endParam) 1000.0)
)
6
nil
)
)
)
)
