向量
在編程上面,向量用(X,Y,Z)表示,也就是(123.156,194.156,215,00)
唉!那么和點的數據結構是一樣的.
它主要的目的是累加方向和計算夾角.
記住一個原則就行:它表示原點(0,0,0)到(X,Y,Z)的方向和長度.
那么兩個點怎么轉換成向量呢？
數學上就是b點-a點==(x2-x1,y2-y1,z2-z1),而cad提供了一個方法pta.GetVectorTo(ptb)
而法向量是什么?
X軸Y軸形成的平面(XOY面)的法向量就是Z軸.
還有的人喜歡把向量和點用標記區分:(x,y,z,標記)
那么他們運算是豎式計算:
(1,2,3,1) 點
+(4,5,6,0) 向量
---------------
=(5,7,9,1) 點
我擼了一個PointV類更好的說明了4個數字表達的優美,其中涉及了后續的矩陣齊次坐標等操作.
而單位向量是什么?
cad是GetNormal,但是我用數學實現了一下,
看下面點乘中的 public static Point3d GetUnitNormal(this Point3d ob_vector) 了解得知 😃
點乘的意義
下圖所示,通過向量OA,向量OB,求OP的距離和P點坐標.
現象就是A點投影到OB向量,呈現一個90度直角和P點,而且OA是任意角度的,不需要和x軸重合.
雖然說B站有視頻教程,但是直接看代碼貌似更方便理解.
再深入,通過A,B,C求E點,就是 E == (AC點乘AB).至於為什么要算E點,可以看 數學篇 cad.net 葛立恆凸包算法和面積最小包圍盒.
測試命令:
/// <summary>
/// 測試點乘代碼,點1=o,點2=a,點3=b
/// </summary>
[CommandMethod("test_DotProductOrCross")]
public void test_DotProductOrCross()
{
var dm = Acap.DocumentManager;
var ed = dm.MdiActiveDocument.Editor;
var ppo = new PromptPointOptions("");
var pts = new List<Point3d>();
for (int i = 0; i < 3; i++)
{
ppo.Message = $"{Environment.NewLine}點位置{i}:";
var pot = ed.GetPoint(ppo);
if (pot.Status != PromptStatus.OK)
return;
pts.Add(pot.Value);
}
var pickPoint = MathHelper.DotProduct(pts[0], pts[1], pts[2]);
ed.WriteMessage($"{Environment.NewLine}點乘P點是:" + pickPoint);
var pickLength = MathHelper.DotProductLength(pts[0], pts[1], pts[2]);
ed.WriteMessage($"{Environment.NewLine}點乘OP長度是:" + pickLength);
var ve1 = new Vector3d(1, 3, 4);
var ve2 = new Vector3d(2, 5, 9);
var ss1 = MathHelper.CrossNormal(ve1, ve2);
ed.WriteMessage($"{Environment.NewLine}叉乘測試****法向量是:" + ss1); //(7,-1,-1)
var aa = MathHelper.CrossAclockwise(pts[0], pts[1], pts[2]);
ed.WriteMessage($"{Environment.NewLine}叉乘測試*****逆時針:" + aa.ToString());
}
點積函數:
/// <summary>
/// 點積,求坐標
/// </summary>
/// <param name="o">原點</param>
/// <param name="a"></param>
/// <param name="b"></param>
/// <returns>坐標點</returns>
public static Point3d DotProduct(Point3d o, Point3d a, Point3d b)
{
//點乘就是將oa向量投射到ob向量上面,求得op坐標(也就是呈現90度角的坐標)
var oa = o.GetVectorTo(a);
var obUnit = o.GetVectorTo(b).GetUnitNormal();
var dot = (oa.X * obUnit.X) + (oa.Y * obUnit.Y) + (oa.Z * obUnit.Z);
var p = obUnit * dot;
//如果O不在坐標原點,則需要平移
return new Point3d(p.X + o.X, p.Y + o.Y, p.Z + o.Z);
}
/// <summary>
/// 點積,求長度
/// </summary>
/// <param name="o">原點</param>
/// <param name="a"></param>
/// <param name="b"></param>
/// <returns>長度</returns>
public static double DotProductLength(Point3d o, Point3d a, Point3d b)
{
var p = DotProduct(o, a, b);
return o.GetDistanceTo(p);
}
/// <summary>
/// 獲取單位向量,仿照向量的GetNormal的數學實現,了解原理,用了Point3d代替向量
/// </summary>
/// <param name="ob_vector"></param>
/// <returns></returns>
//https://www.bilibili.com/video/BV1qb411M7wL?from=search&seid=9280697047969917119
public static Point3d GetUnitNormal(this Point3d ob_vector)
{
#if true2
var ob_length = Point3d.Origin.GetDistanceTo(ob_vector);
return ob_vector / ob_length;
#else
//兩點的距離,向量模長數學版 √(x²+y²+z²)
var ob_length = Math.Sqrt(ob_vector.X * ob_vector.X +
ob_vector.Y * ob_vector.Y +
ob_vector.Z * ob_vector.Z);
return new Point3d(ob_vector.X / ob_length, ob_vector.Y / ob_length, ob_vector.Z / ob_length)
#endif
}
叉乘的意義
叉乘有三個意義!
二維叉乘求面積:
A:返回值是數值,有正負,表示繞原點四象限的位置變換,也就是有向面積,面積屬性可用來求解凸包的單峰函數.
見圖理解一下交叉相乘求面積的數學知識(就是交叉相乘再相減),下方所有Cross函數也是基於此制作.
求平行四邊形面積本來就是底乘以高只是換成坐標的求解方式,也可以參考知乎
代碼同下.
二維叉乘求方向:
B:返回值是數值,有正負,如果>0就是逆時針.
/// <summary>
/// 叉積,二維叉乘計算
/// </summary>
/// <param name="a">傳參是向量,表示原點是0,0</param>
/// <param name="b">傳參是向量,表示原點是0,0</param>
/// <returns>其模為a與b構成的平行四邊形面積</returns>
public static double Cross(Vector2d a, Vector2d b)
{
return a.X * b.Y - a.Y * b.X;
}
/// <summary>
/// 叉積,二維叉乘計算
/// </summary>
/// <param name="o">原點</param>
/// <param name="a">oa向量</param>
/// <param name="b">ob向量,此為判斷點</param>
/// <returns>返回值有正負,表示繞原點四象限的位置變換,也就是有向面積</returns>
public static double Cross(Point2d o, Point2d a, Point2d b)
{
return Cross(o.GetVectorTo(a),o.GetVectorTo(b));
}
/// <summary>
/// 叉積,逆時針方向為真
/// </summary>
/// <param name="o">直線點1</param>
/// <param name="a">直線點2</param>
/// <param name="b">判斷點</param>
/// <returns>b點在oa的逆時針<see cref="true"/></returns>
public static bool CrossAclockwise(Point2d o, Point2d a, Point2d b)
{
return Cross(o, a, b) > -1e-6;//浮點數容差考慮
}
三維叉乘:
C:返回值是向量.例如(X軸叉乘Y軸==Z軸),如果叉乘的不是XY軸,而是兩個三維向量(並不需要90度),那么就可以求出一個平面的法向量.
同樣叫做A叉乘B,左右手叉乘的方式剛好方向相反.圖來自,必要的視頻參考
左右手坐標系相互轉換:
如果將左a重合到右a,將左b重合到右b,那么左c和右c的關系剛好是正負的關系,結果*-1即可將左c變成右c.
(說明這個僅僅是為了,鏡像導致cad圖元的210組碼(normal)為負數而做出關聯修改方式)
三維叉乘(X軸*Y軸==Z軸)
/// <summary>
/// 叉積,求法向量
/// </summary>
/// <param name="a">向量a</param>
/// <param name="b">向量b</param>
/// <returns>右手坐標系系法向量</returns>
public static Vector3d CrossNormal(Vector3d a, Vector3d b)
{
//叉乘:依次用手指蓋住每列,交叉相乘再相減,注意主副順序
//(a.X a.Y a.Z)
//(b.X b.Y b.Z)
return new Vector3d(a.Y * b.Z - b.Y * a.Z, //主-副(左上到右下是主,左下到右上是副)
a.Z * b.X - b.Z * a.X, //副-主
a.X * b.Y - b.X * a.Y); //主-副
}
矩陣乘法
首先要學習一下矩陣乘法,主要掌握相乘再相加.
但是編程應用在於調用方法,只需要了解:
構造一個xx矩陣 * 點矩陣 == 變換后的矩陣.
齊次坐標
由於向量乘法需要滿足a矩陣行和b矩陣列個數相同才能進行矩陣乘法,所以需要引入齊次坐標.
假設使用2x2的矩陣,是沒有辦法描述平移操作的,只有引入3x3矩陣形式,才能統一描述二維中的平移、旋轉、縮放操作. 同理必須使用4x4的矩陣才能統一描述三維的變換.
若不然則可嘗試將平移矩陣的結尾1去掉,你會發現不滿足矩陣乘法了.
平移矩陣
透視矩陣
這里還有一個好玩的透視矩陣
旋轉矩陣
有了矩陣乘法的思想了,所以應該去看看一些人家寫好的方法了,我這里就沒必要重復敘述了,
他寫得太好了: 旋轉矩陣
然后就可以去隨便找一個c#矩陣代碼,或者使用math.net庫,抄一下跑一下,懂得原理后再抄代碼心里舒服多啦
此處代碼是我做的PointV類中使用到的.
public partial class Transform
{
/// <summary>
/// 旋轉矩陣
/// </summary>
/// <param name="angle">角度</param>
/// <param name="axis">任意旋轉軸</param>
/// <returns></returns>
public static Matrix GetRotateTransform(double angle, PointV axis)
{
angle = -angle;//保證是逆時針
var cos = Math.Cos(angle);
var sin = Math.Sin(angle);
var cosMinus = 1 - cos;
axis = axis.GetUnitNormal();
var u = axis.X;
var v = axis.Y;
var w = axis.Z;
var pOut = new double[4, 4];
{
pOut[0, 0] = cos + u * u * cosMinus;
pOut[0, 1] = u * v * cosMinus + w * sin;
pOut[0, 2] = u * w * cosMinus - v * sin;
pOut[0, 3] = 0;
pOut[1, 0] = u * v * cosMinus - w * sin;
pOut[1, 1] = cos + v * v * cosMinus;
pOut[1, 2] = w * v * cosMinus + u * sin;
pOut[1, 3] = 0;
pOut[2, 0] = u * w * cosMinus + v * sin;
pOut[2, 1] = v * w * cosMinus - u * sin;
pOut[2, 2] = cos + w * w * cosMinus;
pOut[2, 3] = 0;
pOut[3, 0] = 0;
pOut[3, 1] = 0;
pOut[3, 2] = 0;
pOut[3, 3] = 1;
}
return new Matrix(pOut);
}
/// <summary>
/// 應用旋轉矩陣
/// </summary>
/// <param name="pts">點集</param>
/// <param name="rotateMat">旋轉矩陣</param>
/// <returns></returns>
public static PointV[] WarpRotate(PointV[] pts, Matrix rotateMat)
{
#if DEBUG
if (rotateMat.Rows != rotateMat.Cols || rotateMat.Cols != 4)
throw new ArgumentNullException("WarpRotate矩陣大小不對");
#endif
var outPts = new List<PointV>();
foreach (var ptItem in pts)
outPts.Add(SetRotateMat(ptItem, rotateMat));
return outPts.ToArray();
}
/// <summary>
/// 應用旋轉矩陣
/// </summary>
/// <param name="ptItem">點</param>
/// <param name="rotateMat">旋轉矩陣</param>
/// <returns></returns>
public static PointV WarpRotate(PointV ptItem, Matrix rotateMat)
{
#if DEBUG
if (rotateMat.Rows != rotateMat.Cols || rotateMat.Cols != 4)
throw new ArgumentNullException("WarpRotate矩陣大小不對");
#endif
return SetRotateMat(ptItem, rotateMat);
}
static PointV SetRotateMat(PointV ptItem, Matrix rotateMat)
{
var ptRo = rotateMat * new Matrix(ptItem.ToArray());//左乘
return new PointV(ptRo.ToArray());
}
}
/////////////////////////////調用例子//////////////////////////////
/// <summary>
/// 旋轉
/// </summary>
/// <param name="angle">角度</param>
/// <param name="vector">任意轉軸</param>
/// <param name="centerPoint">繞點</param>
/// <returns></returns>
public PointV RotateBy(double angle, PointV vector, PointV centerPoint)
{
var roMat = Transform.GetRotateTransform(angle, vector);
var pt = this - centerPoint;//繞點旁邊平移到原點旁邊
pt = Transform.WarpRotate(pt, roMat);
return pt + centerPoint;//原點旁邊平移到繞點旁邊
}
二維旋轉
如果不寫成矩陣形式,寫成函數形式呢？貌似更方便入門和理解.
我發現人家的圖不是很好,自己做了一個.
public static Point2d Rotation(Point2d pt, double al)
{
var xa = pt.X;
var ya = pt.Y;
var cos = Math.Cos(al);
var sin = Math.Sin(al);
var xb = xa * cos + ya * sin;
var yb = ya * cos - xa * sin;
return new Point2d(xb, yb);
}
有了二維旋轉之后,你就能知道點是如何在原點上圍繞Z軸旋轉了.
因為Z軸在XOY面上有無數個,所以必須繞原點,否則平移到原點再繞,最后再平移回去.
而二維擴展到三維,無非就是點Z的位置是0,其他軸旋轉就是其他軸的點是0,
這樣就可以明白了XYZ軸的旋轉了,而任意軸旋轉看下面的三維旋轉.
三維旋轉
那么三維的點怎么從用戶坐標系到世界坐標系？或者世界轉用戶？
首先我們要知道這個點所在的坐標系(用戶坐標系,物體坐標系,相對坐標系),
拿到這個坐標系的Z軸(提及過的任意軸)
- 重合用戶原點到世界原點(就是求向量,向量就是原點開始)
- 重合用戶Z軸到世界Z軸,使得剩下世界Z軸可以旋轉(不就可以套用了二維旋轉了嗎).
- 記錄下XYZ三個軸的轉動角度.
變換必然會依照順序,出現XYZ按照夾角alx,aly,alz轉動,
逆向則為ZYX按照夾角-alz,-aly,-alx轉動.
這樣即可變換兩個坐標系. - 平移回去用戶原點.
步驟2就是核心:
重合Z軸只需要旋轉XY兩個軸即可,如圖所示: 原點和P點形成了一個向量(我習慣稱作用戶z軸)
2.1. 把P點投影到YOZ面.獲取oq與世界Z軸的夾角alx.
2.2. 通過alx旋轉世界X軸后,P點將會變換為r點,r點在XOZ平面上.獲取or與世界Z軸角度aly.
2.3. 接着獲取alz,看代碼就曉得了.
public static partial class MathHelper
{
/// 原理見:http://www.cnblogs.com/graphics/archive/2012/08/10/2627458.html
/// 以及:http://www.doc88.com/p-786491590188.html
/// <summary>
/// 輸入一個法向量(用戶Z軸),獲取它與世界坐標X軸和Y軸的夾角,
/// 旋轉這兩個夾角可以重合世界坐標的Z軸
/// </summary>
/// <param name="userZxis">用戶坐標系的Z軸</param>
/// <param name="alx">輸出X軸夾角</param>
/// <param name="aly">輸出Y軸夾角</param>
public static void ToWcsAngles(this Vector3d userZxis, out double alx, out double aly)
{
//處理精度
double X = Math.Abs(userZxis.X) < 1e-10 ? 0 : userZxis.X;
double Y = Math.Abs(userZxis.Y) < 1e-10 ? 0 : userZxis.Y;
double Z = Math.Abs(userZxis.Z) < 1e-10 ? 0 : userZxis.Z;
//YOZ面==旋轉X軸..投影的
var oq = Point3d.Origin.GetVectorTo(new Point3d(0, Y, Z));
alx = Vector3d.ZAxis.GetAngleTo(oq);
alx = oq.Y > 0 ? alx : Pi2 - alx;
alx = Math.Abs(Pi2 - alx) < 1e-10 ? 0 : alx;
//XOZ面==旋轉Y軸..旋轉的
var userZ = Math.Pow(Y * Y + Z * Z, 0.5);
var or = Point3d.Origin.GetVectorTo(new Point3d(X, 0, userZ));
aly = Vector3d.ZAxis.GetAngleTo(or);
aly = or.X < 0 ? aly : Pi2 - aly;
aly = Math.Abs(Pi2 - aly) < 1e-10 ? 0 : aly;
}
/// <summary>
/// X軸到向量的弧度,cad的獲取的弧度是1PI,所以轉換為2PI(上小,下大)
/// </summary>
/// <param name="ve">向量</param>
/// <returns>弧度</returns>
public static double GetAngle2XAxis(this Vector2d ve)
{
double alz = Vector2d.XAxis.GetAngleTo(ve);//觀察方向不要設置
alz = ve.Y > 0 ? alz : Pi2 - alz; //逆時針為正,如果-負值控制正反
alz = Math.Abs(Pi2 - alz) < 1e-10 ? 0 : alz;
return alz;
}
/// <summary>
/// X軸到向量的弧度,cad的獲取的弧度是1PI,所以轉換為2PI(上小,下大)
/// </summary>
public static double GetAngle2XAxis(this Point2d startPoint, Point2d endtPoint)
{
return startPoint.GetVectorTo(endtPoint).GetAngle2XAxis();
}
/// <summary>
/// 三個軸角度
/// 重合兩個坐標系的Z軸,求出三個軸旋轉角度,后續利用旋轉角度正負和次序來變換用戶和世界坐標系
/// </summary>
/// <param name="ucs">用戶坐標系</param>
/// <param name="alx">坐標系間的X軸旋轉角度</param>
/// <param name="aly">坐標系間的Y軸旋轉角度</param>
/// <param name="alz">坐標系間的Z軸旋轉角度</param>
public static void ToWcsAngles(this CoordinateSystem3d ucs, out double alx, out double aly, out double alz)
{
//ucs可能帶有新原點,設置到0,世界坐標系原點重合
var ucs2o = new CoordinateSystem3d(Point3d.Origin, ucs.Xaxis, ucs.Yaxis);
//XY軸通過叉乘求得Z軸,但是這個類幫我求了
ucs2o.Zaxis.ToWcsAngles(out alx, out aly);
//使用戶X軸與世界XOY面共面,求出Z軸旋轉角
var newXa = ucs2o.Xaxis.RotateBy(alx, Vector3d.XAxis)
.RotateBy(aly, Vector3d.YAxis);
alz = -newXa.GetAngle2XAxis();
}
private static Point3d Wcs2Ucs(this Point3d pt, CoordinateSystem3d ucs)
{
ucs.ToWcsAngles(out double alx, out double aly, out double alz);
pt = new Point3d(pt.X - ucs.Origin.X,
pt.Y - ucs.Origin.Y,
pt.Z - ucs.Origin.Z);
pt = pt.RotateBy(alx, Vector3d.XAxis, Point3d.Origin)
.RotateBy(aly, Vector3d.YAxis, Point3d.Origin)
.RotateBy(alz, Vector3d.ZAxis, Point3d.Origin);
return pt;
}
private static Point3d Ucs2Wcs(this Point3d pt, CoordinateSystem3d ucs)
{
ucs.ToWcsAngles(out double alx, out double aly, out double alz);
//把pt直接放在世界坐標系上,此時無任何重合.
//進行逆旋轉,此時向量之間重合.
pt = pt.RotateBy(-alz, Vector3d.ZAxis, Point3d.Origin)
.RotateBy(-aly, Vector3d.YAxis, Point3d.Origin)
.RotateBy(-alx, Vector3d.XAxis, Point3d.Origin);
//平移之后就是點和點重合,此點為世界坐標系.
pt = new Point3d(pt.X + ucs.Origin.X,
pt.Y + ucs.Origin.Y,
pt.Z + ucs.Origin.Z);
return pt;
}
/// <summary>
/// 把一個點從一個坐標系統變換到另外一個坐標系統
/// </summary>
/// <param name="userPt">來源點</param>
/// <param name="source">來源坐標系</param>
/// <param name="target">目標坐標系</param>
/// <returns>目標坐標系的點</returns>
public static Point3d Transform(this Point3d userPt, CoordinateSystem3d source, CoordinateSystem3d target)
{
//世界坐標是獨特的
var wcs = Matrix3d.Identity.CoordinateSystem3d;
Point3d pt = Point3d.Origin;
if (Equals(source, target))
pt = userPt;
//al的角度是一樣的,旋轉方式取決於正負號
else if (!Equals(source, wcs) && !Equals(target, wcs))//用戶轉用戶
pt = userPt.Ucs2Wcs(source).Wcs2Ucs(target);
if (Equals(target, wcs))
pt = userPt.Ucs2Wcs(source);
else if (!Equals(target, wcs))
pt = userPt.Wcs2Ucs(target);
return pt;
}
/// <summary>
/// 獲取坐標系統三維
/// </summary>
/// <param name="ed"></param>
/// <param name="wcs"></param>
/// <param name="ucs"></param>
public static void GetWcsUcs(this Editor ed, out CoordinateSystem3d wcs, out CoordinateSystem3d ucs)
{
Matrix3d Ide_wcs = Matrix3d.Identity;//獲取世界坐標系
wcs = Ide_wcs.CoordinateSystem3d;
Matrix3d used_ucs = ed.CurrentUserCoordinateSystem;//當前用戶坐標系
ucs = used_ucs.CoordinateSystem3d;
}
}
凸度
參考
mac-lee,它已經說得很好了
bulgeconversion-參考鏈接1
bulgeconversion-參考鏈接2
已知圓弧的起點端點和凸度計算圓心
說明
A: 此函數將把用於定義Acad圓弧轉多邊形弧段(頂點,凸度).
B: 還可以從 凸度==0 判斷點是否在直線上.
返回:凸度可能是正的,也可能是負的,這取決於通過這三個點的弧線是順時針路徑還是逆時針.
為了盡量避免歧義,在交流的時候不要使用中點和終點,改為腰點和尾點.
lisp的處理方式:
3-Points to Bulge縮寫是 3P to B..感覺怪怪吼!
;; 3-Points to Bulge - Lee Mac
(defun LM:3p->bulge ( pt1 pt2 pt3 )
((lambda ( a ) (/ (sin a) (cos a))) (/ (+ (- pi (angle pt2 pt1)) (angle pt2 pt3)) 2))
)
c#的處理方式:
/// <summary>
/// 凸度的測試命令,點1=圓弧頭點,點2=圓弧腰點,點3=圓弧尾點
/// </summary>
[CommandMethod("test_BulgeCenter")]
public void test_BulgeCenter()
{
var dm = Acap.DocumentManager;
Editor ed = dm.MdiActiveDocument.Editor;
var ppo = new PromptPointOptions("");
var pts = new List<Point3d>();
for (int i = 0; i < 3; i++)
{
ppo.Message = $"{Environment.NewLine}點位置{i}:";
var pot = ed.GetPoint(ppo);
if (pot.Status != PromptStatus.OK)
return;
pts.Add(pot.Value);
}
var getBulge = MathHelper.GetArcBulge(pts[0], pts[1], pts[2]);
var center = MathHelper.GetArcBulgeCenter(pts[0].ToPoint2d(), pts[2].ToPoint2d(), getBulge);
ed.WriteMessage($"{Environment.NewLine}凸度是:{getBulge}");
ed.WriteMessage($"{Environment.NewLine}圓心點是:" + center.DebugString());
var arcLength = MathHelper.GetArcLength(pts[0].ToPoint2d(), pts[2].ToPoint2d(), getBulge);
ed.WriteMessage($"{Environment.NewLine}弧長是:{arcLength}");
}
計算凸度
/// http://www.lee-mac.com/bulgeconversion.html
/// <summary>
/// 求凸度,判斷三點是否一條直線上
/// </summary>
/// <param name="arc1">圓弧起點</param>
/// <param name="arc2">圓弧腰點</param>
/// <param name="arc3">圓弧尾點</param>
/// <returns>逆時針為正,順時針為負</returns>
public static double GetArcBulge(Point2d arc1, Point2d arc2, Point2d arc3)
{
double dStartAngle = GetAngle2XAxis(arc2, arc1);
double dEndAngle = GetAngle2XAxis(arc2, arc3);
//求的P1P2與P1P3夾角
var talAngle = (Math.PI - dStartAngle + dEndAngle) / 2;
//凸度==拱高/半弦長==拱高比值/半弦長比值
//有了比值就不需要拿到拱高值和半弦長值了,因為接下來是相除得凸度
double bulge = Math.Sin(talAngle) / Math.Cos(talAngle);
//處理精度
if (bulge > 0.9999 && bulge < 1.0001)
bulge = 1;
else if (bulge < -0.9999 && bulge > -1.0001)
bulge = -1;
else if (Math.Abs(bulge) < 1e-10)
bulge = 0;
return bulge;
}
/// <summary>
/// 求凸度,判斷三點是否一條直線上
/// </summary>
/// <param name="arc" >圓弧</ param>
/// <returns></returns>
public static double GetArcBulge(this Arc arc)
{
//還有一種求凸度方法是tan(圓心角 / 4),但是丟失了方向,
//再用叉乘來判斷腰點方向,從而凸度是否 * -1
double bulge = Math.Tan(arc.TotalAngle / 4); //凸度是圓心角的四分之一的正切
Point3d midpt = arc.GetPointAtDist(arc.GetDistAtPoint(arc.EndPoint) / 2);//圓弧的腰點
Vector3d vmid = midpt - arc.StartPoint; //起點到腰點的向量
Vector3d vend = arc.EndPoint - arc.StartPoint; //起點到尾點的向量
Vector3d vcross = vmid.CrossProduct(vend); //叉乘求正負
//根據右手定則,腰點向量在尾點向量右側,則叉乘向量Z值為正,圓弧為逆時針
if (vcross.Z < 0)
bulge *= -1;
//處理精度
if (bulge > 0.9999 && bulge < 1.0001)
bulge = 1;
else if (bulge < -0.9999 && bulge > -1.0001)
bulge = -1;
else if (Math.Abs(bulge) < 1e-10)
bulge = 0;
return bulge;
}
凸度求圓心
/// http://bbs.xdcad.net/thread-722387-1-1.html
/// <summary>
/// 凸度求圓心
/// </summary>
/// <param name="arc1">圓弧頭點</param>
/// <param name="arc3">圓弧尾點</param>
/// <param name="bulge">凸度</param>
/// <returns>圓心</returns>
public static Point2d GetArcBulgeCenter(Point2d arc1, Point2d arc3, double bulge)
{
if (bulge == 0)
throw new ArgumentNullException("凸度為0,此線是平的");
var x1 = arc1.X;
var y1 = arc1.Y;
var x2 = arc3.X;
var y2 = arc3.Y;
var b = (1 / bulge - bulge) / 2;
var x = (x1 + x2 - b * (y2 - y1)) / 2;
var y = (y1 + y2 + b * (x2 - x1)) / 2;
return new Point2d(x, y);
}
凸度求弧長
/// <summary>
/// 凸度求弧長
/// </summary>
/// <param name="arc1">圓弧頭點</param>
/// <param name="arc3">圓弧尾點</param>
/// <param name="bulge">凸度</param>
/// <returns></returns>
public static double GetLength(Point2d arc1, Point2d arc3, double bulge)
{
var bowLength = arc1.DistanceTo(arc3); //弦長
var bowLength2 = bowLength / 2; //半弦
var archHeight = Math.Abs(bulge) * bowLength2; //拱高==凸度*半弦
//根據三角函數: (弦長/2)²+(半徑-拱高)²=半徑²
//再根據:完全平方公式變形: (a+b)²=a²+2ab+b²、(a-b)²=a²-2ab+b²
var r = (bowLength2 * bowLength2 + archHeight * archHeight) / (2 * archHeight); //半徑
//求圓心角:一半圓心角[(對邊==半弦長,斜邊==半徑)]; *2就是完整圓心角
var asin = Math.Asin(bowLength2 / r) * 2; //反正弦
//弧長公式: 弧長=絕對值(圓心弧度)*半徑...就是單位*比例..縮放過程
var arcLength = asin * r;
return arcLength;
}
反三角函數
由於我經常忘記三角函數和反三角函數之間的關系...敲了一點cad的測試命令.
[CommandMethod("test_InverseTrigonometric")]
public void test_InverseTrigonometric()
{
var doc = Acap.DocumentManager.MdiActiveDocument;
var db = doc.Database;
var ed = doc.Editor;
ed.WriteMessage($"{Environment.NewLine}反三角函數");
var ppo = new PromptPointOptions("")
{
AllowNone = false,
};
var pts = new List<Point3d>();
for (int i = 0; i < 3; i++)
{
ppo.Message = $"{Environment.NewLine}選擇三角形角點{i}";
var gp = ed.GetPoint(ppo);
if (gp.Status != PromptStatus.OK)
return;
pts.Add(gp.Value);
}
//由於某些實際情況缺少一條邊,此時需要求夾角,則使用反三角函數
var 鄰邊 = pts[0].DistanceTo(pts[1]);
var 對邊 = pts[1].DistanceTo(pts[2]);
var 斜邊 = pts[2].DistanceTo(pts[0]);
//Math.Cos(angle) == 鄰邊 / 斜邊;
//Math.Sin(angle) == 對邊 / 斜邊;
//Math.Tan(angle) == 對邊 / 鄰邊;
var acos = Math.Acos(鄰邊 / 斜邊); //反余弦
var asin = Math.Asin(對邊 / 斜邊); //反正弦
var atan = Math.Atan(對邊 / 鄰邊); //反正切
var atan2 = Math.Atan2(56, 30) * 180 / Math.PI;//(Y,X)與(0,0 X軸)的夾角
ed.WriteMessage($"{Environment.NewLine}角度a是:{asin},{atan},{acos}");
}
它們的實現方式是泰勒展開
高精度反三角函數的實現
導數應用(光線反射)
計算方式是:
一階導數y'=dy/dx=df(x)/dx
二階導數y''=dy'/dx=[d(dy/dx)]/dx=d²y/dx²=d²f(x)/dx²
...
看一點偽代碼,形象的理解一下:
//一階導數
double GetFirstDerivative(double dy)
{
double dx = 1e-10;//趨近於0
return dy / dx;
}
double 二階導數(double dy)
{
return GetFirstDerivative(GetFirstDerivative(dy));
}
double 三階導數(double dy)
{
return GetFirstDerivative(GetFirstDerivative(GetFirstDerivative(dy)));
}
//四階五階六階...無限嵌套就是無限階...
cad實現
起初我根據群友認為的二階導數就是法線,實際上求了一個特殊解:只有圓弧的二階導數等於法線.
而要擴展到曲線的法線,則通用解是切線旋轉90度.
using Autodesk.AutoCAD.EditorInput;
using Autodesk.AutoCAD.Geometry;
using Autodesk.AutoCAD.Runtime;
using Autodesk.AutoCAD.DatabaseServices;
using Acap = Autodesk.AutoCAD.ApplicationServices.Application;
using System.Collections.Generic;
using System;
namespace JoinBox
{
public partial class CmdTest
{
[CommandMethod("test_Derivative")]
public void test_Derivative()
{
var dm = Acap.DocumentManager;
var doc = dm.MdiActiveDocument;
var db = doc.Database;
var ed = doc.Editor;
//兩點確定射線入射角
var pts = new List<Point3d>();
var ppo = new PromptPointOptions("");
for (int i = 0; i < 2; i++)
{
ppo.Message = $"{Environment.NewLine}測試導數,點位置{i + 1}:";
var pot = ed.GetPoint(ppo);
if (pot.Status != PromptStatus.OK)
return;
pts.Add(pot.Value);
}
var ray0 = new Ray
{
BasePoint = pts[0],
SecondPoint = pts[1],
};
db.Action(tr => {
//獲取屏幕內的圖元,再獲取射中的
var ids = SelectViewportObjectId(ed, tr);
if (ids == null)
return;
foreach (var id in ids)
{
var ent = id.ToEntity(tr);
if (ent is Curve cur)
{
if (ent is not Arc)//測試時候過濾掉無用注釋用
continue;
var ptc = new Point3dCollection();
ray0.IntersectWith(cur, Intersect.OnBothOperands, ptc.Collection, (int)IntPtr.Zero, (int)IntPtr.Zero); //得出的所有交點
if (ptc.Count == 0)
continue;
var intPt = ptc[0]; //交點可能是3d點
// 一階導數是切線
Vector3d yt1 = cur.GetFirstDerivative(intPt);
var ray1 = new Ray
{
BasePoint = intPt,
SecondPoint = new Point3d(yt1.X + intPt.X, yt1.Y + intPt.Y, 0),
};
ray1.ColorIndex = 1;
tr.AddEntityToMsPs(db, ray1);
#if true3
// 只有圓弧的二階導數是法線
Vector3d yt2 = cur.GetSecondDerivative(intPt);
var ray2 = new Ray
{
BasePoint = intPt,
SecondPoint = new Point3d(yt2.X + intPt.X, yt2.Y + intPt.Y, 0),
};
#else
// 曲線的法線通用解是:切線是旋轉90度
var ray2 = ray1.Clone() as Ray;
ray2.EntityRotate(Math.PI / 2, Vector3d.ZAxis, intPt);
#endif
ray2.ColorIndex = 2;
tr.AddEntityToMsPs(db, ray2);
// 射出光線
var ray3 = new Ray
{
BasePoint = intPt,
SecondPoint = pts[0],
};
ray3.EntityMirror(ray2.BasePoint, ray2.SecondPoint, out Entity ray3mr);
ray3mr.ColorIndex = 3;
tr.AddEntityToMsPs(db, ray3mr);
}
}
});
}
/// <summary>
/// 獲取當前屏幕的所有圖元
/// </summary>
/// <param name="ed"></param>
/// <param name="tr"></param>
/// <returns></returns>
ObjectId[] SelectViewportObjectId(Editor ed, Transaction tr)
{
var view = ed.GetCurrentView();
var cpt = view.CenterPoint;
var w = view.Width / 2;
var h = view.Height / 2;
//先忽略視口原點角度等信息
var min = new Point3d(cpt.X - w, cpt.Y - h, 0);
var max = new Point3d(cpt.X + w, cpt.Y + h, 0);
var ssget = ed.SelectCrossingWindow(min, max);
if (ssget.Status != PromptStatus.OK)
return null;
return ssget.Value.GetObjectIds();
}
}
}
public static partial class EntityEdit
{
// 鏡像
public static Matrix3d EntityMirror(this Entity ent,
Point3d point1, Point3d point2, out Entity entity)
{
//計算變換矩陣
Matrix3d mt = Matrix3d.Mirroring(new Line3d(point1, point2));//鏡像矩陣
entity = ent.GetTransformedCopy(mt); //這里可能要替換成克隆
return mt;
}
}
public static partial class EntityAdd
{
/// <summary>
/// 將圖形添加到數據庫的當前空間中
/// </summary>
/// <param name="db">圖形數據庫</param>
/// <param name="ent">圖形對象</param>
/// <returns>圖形的ObjectId</returns>
public static ObjectId AddEntityToMsPs(this Transaction tr, Database db, Entity ent)
{
ObjectId entId;
//在位編輯的時候自動加到塊內了,而不是模型空間,因為會把臨時的數據拷貝回去塊表記錄上面
//出現eLockViolation 是因 CommandFlags.Session | CommandFlags.Redraw 又沒有鎖文檔
var btRec = tr.GetObject(db.CurrentSpaceId, OpenMode.ForWrite) as BlockTableRecord;
//加圖元到塊表記錄,設定返回值
entId = btRec.AppendEntity(ent);
//更新數據
tr.AddNewlyCreatedDBObject(ent, true);
btRec.DowngradeOpen();
btRec.Dispose();
return entId;
}
}
子函數:
本博客引用: 委托db.Action cad.net 委托無返回值寫法
本博客引用: ptc.Collection 二度封裝Point3dCollection
二階導數
那曲線的二階導數有什么用呢？
它表達斜率的變化率,所以越凸的地方值越大.
看圖上綠色的部分,長度越長表示值越大,表示你油門踩得越深,然后也表達加速度了.
所以運動控制上面有用.
而對橢圓用二階導數,會得到橢圓中心點,而不是兩個焦點.
橢圓實際上是任意軸觀察的圓形,這樣就求出圓心了.
三階導數
它是變化率的變化率,它有什么用呢？
找拐點
相關閱讀
(完)