根據貝塞爾曲線上的點反算t值

時間 2019-11-06

標籤根據曲線简体版

原文原文鏈接

這是一個項目中遇到的實際需求。場景是一個智能倉庫管理系統，場景裏面有直線和曲線構成的環穿軌道。環穿軌道上面會有小車運動，後臺推進小車的兩個點位A和B，其中A和B都會在軌道上面，前端須要根據這兩個推送點，自動播放小車從A點沿軌道到B點的動畫。下面是項目截圖：前端

項目中使用的是二次貝塞爾曲線，因此本文也主要以二次貝塞爾曲線爲講解重點。要實現上述動畫，須要首先肯定A點和B點在曲線上面的比例值t_a和t_b程序員

最終的需求變成：「根據貝塞爾曲線上的點反算t值」。大概有如下幾種方法。現假設貝塞爾曲線上的點爲點P（後續會用到該點）。canvas

分片迭代

分片迭代是一種近似的方法。咱們知道，二次貝塞爾曲線的公式以下：B(t) = (1-t)2 P0 + 2t(1-t) P1 + t2 * P2其中： $t in $[0,1]，P0爲二次貝塞爾曲線的起始點，P1爲控制點，P2爲終止點。架構

若是你對於上面的知識點不是很熟悉，建議學習貝塞爾曲線相關知識。推薦學習本人的專欄Canvas高級進階, 裏面有專門的章節對貝塞爾曲線進行了全面詳細的講解。本文也是從該專欄的文章中摘錄並適當改編而成的。函數

從以上公式，咱們能夠獲得，對於任意給定的比例值t，能夠求出對應該比例值的點B(t)。分片迭代思路是：如今加設把範圍[0,1]平均分紅N（好比100）等份，造成一系列的比例值t，對於每個t值，求取對應的點B(t) ，而後讓點B(t)和已知在貝塞爾曲線上的點P進行比較，若是點B(t)和點P之間的直線距離在必定的偏差範圍以內，則認爲B(t)等於P，而此時的t值，就是咱們要求的t值。如下是主要代碼：性能

function computeT(p0,p1,p2,p) {
  var t = 0;
  for(var i = 0;i < 1000;i ++){
      var point = getPointOnQuadraticCurve(p0,p1,p2,t);//根據二次貝塞爾曲線公式求B(t)，其中point = B(t)
      if(distance(point,p) < 0.01){ // 判斷point和p點的距離是否在特定偏差以內
        return t;
      }
      t+= 0.001;
  }
  return null;
}複製代碼

上述分片迭代的方法，思路最簡單，最直觀。在精度要求不高的狀況下是能夠知足的。而在精度要求高的時候，即代碼中的「特定偏差」值要很小，可能會出現函數返回值爲null的狀況，在精度要求高的時候要可以計算出值，就要增長迭代次數，此時會極大增長性能消耗。好比上面代碼的迭代次數可能會變成10000甚至10000。學習

迭代方法一樣適用於三次貝塞爾曲線和更加高階的貝塞爾曲線。優化

分片迭代優化版本

上面提到在精度要求高的狀況下，要獲得正確結果，要極大的增長迭代次數，形成性能的極大消耗。有沒有辦法既提升精度，又不大量增長迭代次數呢？通過筆者的思考，發現是能夠的。想一想假設要求的t值在0.5附近，那麼咱們只須要在0.5附近加大分片的數量，而不須要在其餘地方（0.1~0.4，0.6~1.0）增長分片的數量。應此升級版本的思路就是，先用比較粗的分片初步肯定t值的一個大體範圍，再在該範圍之類，比較細的分片肯定t值。注意這是個遞歸的過程，若是在第二次比較細的分片狀況下，仍然不能肯定t值，那麼就肯定一個t值的更小分範圍；重複上面過程，直到找到t值爲止。大體步驟以下：動畫

首先，經過一個小的迭代次數進行分片迭代；
在迭代的過程當中若是找到了符合的比例值t，直接返回；
在迭代的過程當中同時記錄離目標點P最近的t值，若是上一步未找到符合的t值，則進行下一步操做。
上一步找到了離目標點P最近的t值，在t值的附近(t - step,t + step)(其中step爲上一次分片的步進值)進行分片迭代查找，在迭代的過程當中若是找到了符合的比例值t，直接返回。
若是沒找到，重複上面的不斷縮小範圍並加大分片精度的過程。直到找到t值爲止。

下面是示例代碼：職業規劃

function computeT(p0, p1, p2, p,startT = 0,endT = 1) {
  var t = startT;
  var minDistance = Infinity,
      minDistanceT = null;
  var step = (endT - startT) / 100;
  for (var i = 0; i < 100; i++) {
    var point = getPointOnQuadraticCurve(p0, p1, p2, t);
    var dst = distance(point,p);
    if (dst < minDistance) {
      minDistance = dst;
      minDistanceT = t;
    }
    if (dst < 0.0001) {
      return t;
    }
    t += step;
  }
  return computeT(p0, p1, p2, p, minDistanceT - step,minDistanceT + step);
}複製代碼

以上過程雖然增長了必定的迭代次數，可是是常量級別的增長，而非數量級別的增長，因此會極大提升性能。好比目標t值在0.5附近，第一次經過100次迭代能夠肯定t值的範圍在0.4 ~ 0.6之間；而後進行第二次迭代，第二次迭代這次數仍然爲100次，假設肯定t值的範圍在0.51 ~ 0.53之間；而後進行第三次迭代，第三次迭代這次數仍然爲100次，此時能夠獲取t值爲0.516，能夠看出最多值迭代了300次。假設總共通過第N次迭代，每次迭代次數爲M，才找到t值，那麼總共的迭代次數是N * M。

該迭代方法一樣適用於三次貝塞爾曲線和更加高階的貝塞爾曲線。並且相對於未優化的版本，該方法的性能好了不少。是適合全部貝塞爾曲線的比較好的反算t值的方法。

二分法

二分法的思路是：

首先肯定一個起始t值和結束t值t₀和t₁，初始值t₀ = 0，t₁ = 1。
取t₀和t₁的中間值t_m = (t₀+t₁)/2
經過t_m計算出點P_m，若是P_m和目標點P之間的距離在偏差值範圍以內，則t_m爲須要計算的目標t值。
若是上一步P_m和目標點P之間的距離不在偏差值範圍以內，則判斷P_m和目標點P的先後順序，若是P_m在目標點P的前面，則把t_m賦值給t₁；不然把t_m賦值給t₀。
重複以上步驟直到找到合適的t_m值。

上述步驟有一個難點：如何判斷P_m和目標點P的先後順序？對於二次貝塞爾曲線，以下圖所示：其中，P0爲起始點，P2爲終止點，P1爲控制點。二次貝塞爾曲線有以下特色：線段(P1，P0)、(P1，P2)和曲線相切，這也就意味着曲線必定在三角形(P0，P1，P2)以內，並且二次貝塞爾曲線自己不會自身相交，全部咱們能夠有以下結論，

對於曲線上面的點A，直線(P1,A)和線段(P0,P1)相交於點a；對於曲線上面的點B，直線(P1,B)和線段(P0,P1)相交於點b。點A和點B的前後順序與點a和點b的前後順序是一致的，而直線上面的點（a和b）的先後順序是容易判斷的。也就是說若是點a在點b的前面，則點A也在點B的前面，反之亦然。以下圖所示：有了以上的結論，咱們就找到了判斷P_m和目標點P的先後順序的方法。

若是你對上述結論不熟悉，建議學習貝塞爾曲線的相關知識，推薦學習本人的專欄Canvas高級進階, 裏面有專門的章節對貝塞爾曲線進行了全面詳細的講解。本文也是從該專欄的文章中摘錄並適當改編而成的。

有了這個方法，加上前面描述的二分查找的步驟，能夠得出示例代碼以下：

function computeT2(p0,p1,p2,p,startT = 0,endT = 1) {
   var halfT  = (startT + endT) / 2;
   var halfPoint = getPointOnQuadraticCurve(p0,p1,p2,halfT);
    if(distance(halfPoint,p) < 0.0001){
      return halfT;
   }
   //求交點：
   var inter1 = segmentsIntr(p0,p2,p1,p);
   var inter2 = segmentsIntr(p0,p2,p1,halfPoint);
  
   var r1 = interpolationRate(p0,inter1,p2),
       r2 = interpolationRate(p0,inter2,p2);
   
   if(r1 > r2){
       startT = halfT;  
   }else {
     endT = halfT;
   }
   return computeT2(p0,p1,p2,p,startT,endT);
}複製代碼

解方程

前面說過，貝塞爾曲線的公式以下：B(t) = (1-t)2 P0 + 2t(1-t) P1 + t2 * P2其中： $t in $[0,1]，P0爲二次貝塞爾曲線的起始點，P1爲控制點，P2爲終止點。分別表示成x和y的方程，則能夠表示以下：

x_P = (1-t)² x_P0 + 2t(1-t) x_P1 + t² * x_P2
y_P = (1-t)² y_P0 + 2t(1-t) y_P1 + t² * y_P2

實際上就是兩個變量t的二次元方程，取上面任意一個方程，帶入相關的值解方程，方程的解即爲咱們要求的目標t值。

整理方程: x_P = (1-t)² x_P0 + 2t(1-t) x_P1 + t² * x_P2，能夠得出二次方程以下：(x_P2 + x_P0 - 2 x_P1 ) t² + 2(x_P1 - x_P0) t + (x_P0 - x_P) = 0。咱們已知二次方程的： at² + b t + c = 0的解爲：

若是a = 0，則解爲 -c/b
若是a != 0，解以下圖所示：

應此令：

a = (x_P2 + x_P0 - 2 * x_P1 )
b = 2*(x_P1 - x_P0)
c = (x_P0 - x_P)
能夠方便求出方程的解。

須要注意的是，二次方程的解可能會有兩個。若是求出的解有兩個怎麼辦呢。首先咱們知道貝塞爾曲線的t值的範圍是$t in $[0,1]，因此若是有兩個解：

其中一個再也不[0,1]的範圍以內，則另一個解就是目標t值。（注意不可能兩個都不在[0,1]範圍以內，由於咱們知道，目標點P在貝塞爾曲線上面）。
若是兩個解都在[0,1]範圍以內，那就把兩個解再帶入貝塞爾曲線的公式，分別求出兩個B(t)點，那個離目標點P近，就取那個解。

下面是示例代碼,其中函數equation2用於解曲線的方程：

function computeT(p0,p1,p2,p) {
    let interpolationx =  (p1.x - p0.x) / (p2.x - p0.x);
    let tt;
    if(interpolationx >= 0 && interpolationx <= 1){
      let ty = equation2(p0.y,p1.y,p2.y,p.y);
      return ty;
    }else{
      tt = equation2(p0.x,p1.x,p2.x,p.x);
      if(tt.tt1){
        var pointTest = getPointOnQuadraticCurve(p0,p1,p2,tt.tt1);
        if(distance(pointTest,p) < 0.01){
          return tt.tt1;
        }else{
          return tt.tt2;
        }
      }else{
        return tt;
      }
    }
}

function equation2(z0,z1,z2,zp){ // z0、z1，z2表明P0、P一、P2的x座標值或者y座標值，zp表明目標點P的x座標值或者y座標值
    var a = z0 - z1 * 2 + z2,
      b = 2*(z1 - z0),
      c = z0 - zp;
  var tt = null;
  if(a == 0 && b != 0){
    tt = - c / b;
  } else {
    var sq = Math.sqrt( b * b - 4 * a * c );
    var tt1 = (sq - b)/ (2 * a),
        tt2 = (-sq - b) / (2 * a);
    // console.log("tt1,tt2:",tt1,tt2); 
    if((tt1 <= 1 && tt1>= 0) && (tt2 <= 1 && tt2>= 0)){
      return {tt1,tt2};
    }else if(tt1 <= 1 && tt1>= 0){
      tt = tt1;
    }else {
      tt = tt2;
    }
  }
  return tt;
}複製代碼