樣條曲線長度--數值積分

轉載自：https://www.cnblogs.com/flytrace/p/8413255.html

文章目錄

• • • 1.根據長度求坐标思路
    • 2.用數值積分實作s的積分
    • 3.牛頓疊代求解s對應的t
    • 4.直接推導二階貝塞爾曲線的長度積分解析式

前幾天在做錄影機軌道時，解決勻速定長運動問題，剛好了解了題述問題。

首先做完這個曲線上勻速運動的程式設計實作後，是十分愉悅的。這個問題裡，兩個科學史上的偉人，牛頓和高斯，相繼大顯身手，怎能不讓人膜拜？我之前并不了解所涉及的數值分析方法，是以覺得很神奇。

以下為簡潔我會使用向量的表示方法，使用大寫數字代表向量。以catmullrom樣條曲線(spline)為例，對于一個關于t的多次樣條曲線:

Q ( t ) = A t 3 + B t 2 + C t + D Q(t) = At^{3} + Bt^{2} + Ct + D Q(t)=At3+Bt2+Ct+D

1.根據長度求坐标思路

一般來說實際中， t t t的範圍都是[0, 1]。當 t t t在[0, 1]中滑動時， Q ( t ) Q(t) Q(t)的集合即為一條平滑曲線。注意到當 t t t以勻速移動時， Q ( t ) Q(t) Q(t)在沿曲線方向上并不是勻速移動。要做到這一點，顯然曲線長度 s s s是需要考慮的。假設在 t t t點處我們有曲線長度，（根據 s s s長度，計算控制量 t t t，再計算坐标 Q ( t ) Q(t) Q(t)）。

s = L ( t ) s = L(t) s=L(t)

那麼s長度所對應的t即為上面方程的反函數:

t = L − 1 ( s ) t = L^{-1}(s) t=L−1(s)

反帶入 Q ( t ) Q(t) Q(t)，得

Q ( t ) = Q ( L − 1 ( s ) ) = Z ( s ) Q(t) = Q(L^{-1}(s)) = Z(s) Q(t)=Q(L−1(s))=Z(s)

s s s的區間顯然是[0, 曲線總長度]，于是我們可以對 s s s做勻速遞增，進而得到 Q ( t ) Q(t) Q(t)的位置， 進而完成勻速插值。

難點在于求 S = L ( t ) S=L(t) S=L(t)，這實際上是一個求積分的問題。很遺憾，一般來說很難得到一個解析的積分公式的，實際中是通過數值計算方法得到積分值。

s s s通過數值方法得到，也很難求逆得到 t t t，可以通過牛頓疊代，或二分查找得到。

2.用數值積分實作s的積分

數值積分使用Guass-lengendre積分來實作，通過標明的幾個積分點及權重，它神奇而快速的得到精度很高的積分值。我們直接來重點吧，該方法需要一張n階權重表，比如5階權重表:

及一個計算公式:

∫ a b f ( x ) d x = b − a 2 ∑ i = 1 n w i f ( b − a 2 x i + b + a 2 ) \int_{a}^{b}f(x)dx = \frac{b-a}{2}\sum_{i=1}^{n}{w_{i}f(\frac{b-a}{2}x_{i}+\frac{b+a}{2})} ∫ab​f(x)dx=2b−a​i=1∑n​wi​f(2b−a​xi​+2b+a​)

對于樣條曲線，我們有:

d s d t = ∣ ∣ Q ′ ( t ) ∣ ∣ d s = ∣ ∣ Q ′ ( t ) ∣ ∣ d t s = ∫ 0 t ∣ ∣ Q ′ ( x ) ∣ ∣ d x \frac{ds}{dt} = ||Q'(t)|| \\ ds = ||Q'(t)||dt \\ s = \int_{0}^{t}||Q'(x)||dx dtds​=∣∣Q′(t)∣∣ds=∣∣Q′(t)∣∣dts=∫0t​∣∣Q′(x)∣∣dx

• Q ′ ( t ) Q'(t) Q′(t)為 Q ( t ) Q(t) Q(t)的導

  Q ′ ( t ) = 3 A t 2 + 2 B t + C Q'(t) = 3At^{2} + 2Bt + C Q′(t)=3At2+2Bt+C

• ∣ ∣ Q ′ ( t ) ∣ ∣ ||Q'(t)|| ∣∣Q′(t)∣∣為該點處的長度:

  ∣ ∣ Q ′ ( t ) ∣ ∣ = Q x ′ ( t ) 2 + Q y ′ ( t ) 2 ||Q'(t)|| = \sqrt{Q_x'(t)^2 + Q_y'(t)^2} ∣∣Q′(t)∣∣=Qx′​(t)2+Qy′​(t)2

  ​

• 具體對于我們的公式，有:

  s = L ( t ) = ∫ 0 t ∣ ∣ Q ′ ( x ) ∣ ∣ d x = t 2 ∑ i = 1 n w i ∣ ∣ Q ′ ( t 2 x i + t 2 ) ∣ ∣ s= L(t) = \int_{0}^{t}||Q'(x)||dx = \frac{t}{2}\sum_{i=1}^{n}{w_{i}||Q'(\frac{t}{2}x_{i}+\frac{t}{2})||} s=L(t)=∫0t​∣∣Q′(x)∣∣dx=2t​i=1∑n​wi​∣∣Q′(2t​xi​+2t​)∣∣

将表中的值依次帶入，疊代計算即可得到 t t t點處的積分值 s s s。比我之前快速實作時所采用的線性分割方法快了100倍有沒有？特别是該方程的普适性，對于多階多項式都适用而且精度非常高。有種魔性在這裡，高斯是神。

3.牛頓疊代求解s對應的t

現在 s = L ( t ) s=L(t) s=L(t)對我們來說是已知了，那麼如何求 s s s所對應的 t t t呢?這實際上是求方程:

L ( t ) − s = 0 L(t) - s = 0 L(t)−s=0

的根。因為 L ( t ) L(t) L(t)沒有解析表達，我們将使用牛頓疊代法:

b = a − f ( a ) f ′ ( a ) b= a - \frac{f(a)}{f'(a)} b=a−f′(a)f(a)​

将得到的b的值作為a再帶入上述公式，疊代幾次即可逼近 s = f ( x ) s = f(x) s=f(x)的根。具體到我們的公式：

b = a − L ( a ) − s L ′ ( a ) b= a - \frac{L(a) - s}{L'(a)} b=a−L′(a)L(a)−s​

對于特定的 s s s，我們首先要拟定一個初值 a a a，因為樣條曲線性質良好，而且一般說來實際中不會允許2點間出現陡峭的曲線，一定會多加一些點令2點間曲線形狀為凸或凹的，是以我們可以近似認為 t t t在區間[0, 1]的比例，近似等于 s s s與總長度的比例。而顯然總長為L(1),通過前面的高斯積分已經可以算出，于是初始值 a a a我們可拟定為：

a = s L ( 1 ) a = \frac{s}{L(1)} a=L(1)s​

反複疊代(3-4次的精度就不錯了)，得出最終 t t t值，再把該值反帶入 Q ( t ) Q(t) Q(t)中，得到對于給定s值所對應的位置，至此我們完成了樣條曲線的重參數化。

4.直接推導二階貝塞爾曲線的長度積分解析式

https://blog.csdn.net/LANGZI7758521/article/details/52101672