Geometrie/Racionální algoritmus de Boor

Z Wikiknih

Přejít na: navigace, hledání

Obsah

[editovat] De Boor algoritmus

Bod na racionální B-spline křivce můžeme stejně jako u racionálních Bézierových křivek spočítat hned dvěma způsoby.

1. Pomocí algoritmu de Boor pro racionální křivky:

Zvolíme-li u\in\left\langle u_i,u_{i+1}\right\rangle, pak bod na racionální B-spline křivce vypočteme opakovanou lineární interpolací:

P^j_i\left(u\right)=\left(1-\alpha^j_i\right)\cdot\frac{\omega^{j-1}_{i-1}}{\omega^j_i}\cdot P^{j-1}_{i-1}+\alpha^j_i\cdot\frac{\omega^{j-1}_i}{\omega^j_i}\cdot P^{j-1}_i,

kde

\omega^j_i\left(u\right)=\left(1-\alpha^j_i\right)\cdot\omega^{j-1}_{i-1}+\alpha^j_i\cdot\omega^{j-1}_i

a

\alpha^j_i=\frac{u-u_i}{u_{i+k-j+1}}, i=1-k,\dots,1 a P^0_i=P_i.

Bod P^k_1\left(u\right)=P\left(u\right) je hledaný bod racionální B-spline křivky.

Toto rozšíření algoritmu de Boor vzniklo obdobným způsobem jako rozšíření algoritmu de Casteljau pro Bézierovy křivky.

2. Převodem na neracionální křivku a zpětným promítnutím:

Převedeme racionální B-spline na neracionální (o jednu dimenzi výše), provedeme de Boor algoritmus pro neracionální křivku a výsledný bod převedeme (promítneme) zpět.


[editovat] Algoritmizace

[editovat] ComputeKnotVector(int n, int k)

Spočítá uzlový vektor.

Parametry:

  • n - počet kontrolních bodů mínus 1
  • k - stupeň de Boor bázové funkce
void ComputeKnotVector(int n, int k)
{
 parametrization = new ParameterCollection();
 for(int i=0; i<=n+k+1; i++)
 {
  if(i<=k) parametrization.Insert(i,0);
   else
  if(i>n) parametrization.Insert(i,n-k+1);
   else
  parametrization.Insert(i,i-k);
 }
}

[editovat] Vector GetPoint(double t)

Přetížená metoda třídy Curve. Spočítá a vrátí bod na křivce.

Parametry:

  • t - parametr výpočtu
public override Vector GetPoint(double t) 
{
 ParameterCollection u = parametrization;
 int n = ctrlPoly.Count;
 Vector[] points = new Vector[n];
 double[] weight = new double[n];
 for (int i=0;i<n;i++)
 {
  points[i] = ((RichPoint)(ctrlPoly[i])).Locate;
  weight[i] = ((RacionalRichPoint)(ctrlPoly[i])).Weight;
 }
 int l=0;
 while(t>u[l]) l++;
 l--;
 if(l<Degree) l=Degree;
 for(int j=0;j<Degree;j++)
 {
  for(int i=l-Degree+1+j;i<=l;i++)
  {
   double alfa = (t - u[i])/(u[i+Degree-j]-u[i]);
   points[i]=(1-alfa)*points[i-1]+alfa*points[i];
  }
 }
 return points[l];
}