Code-Eli  0.3.6
bezier.hpp
Go to the documentation of this file.
1 /*********************************************************************************
2 * Copyright (c) 2013 David D. Marshall <ddmarsha@calpoly.edu>
3 *
4 * All rights reserved. This program and the accompanying materials
5 * are made available under the terms of the Eclipse Public License v1.0
6 * which accompanies this distribution, and is available at
7 * http://www.eclipse.org/legal/epl-v10.html
8 *
9 * Contributors:
10 * David D. Marshall - initial code and implementation
11 ********************************************************************************/
12 
13 #ifndef eli_geom_curve_bezier_h
14 #define eli_geom_curve_bezier_h
15 
16 #include <iostream>
17 #include <vector>
18 
19 #include "eli/code_eli.hpp"
20 
21 #include "eli/util/tolerance.hpp"
22 
25 
32 
33 // TODO: MOVE THESE TO utility namespace
34 namespace eli
35 {
36  namespace geom
37  {
38  namespace curve
39  {
40  namespace internal
41  {
42  template<typename Derived1, typename Derived2, typename PointType>
43  void build_fit_Ab(Eigen::MatrixBase<Derived1> &A,
44  Eigen::MatrixBase<Derived2> &b,
45  std::vector<typename Derived1::Scalar> &t,
46  const std::vector<PointType, Eigen::aligned_allocator<PointType> > &pts,
47  const typename Derived1::Index &n, const size_t &dim)
48  {
49  typedef Eigen::Matrix<typename Derived1::Scalar, Eigen::Dynamic, Eigen::Dynamic> mat_type;
50  typedef Eigen::Matrix<typename Derived1::Scalar, Eigen::Dynamic, 1> col_type;
51 
52  typename Derived1::Index i, sz(n+1), npts(pts.size()), nrows(std::max(npts, sz));
53  mat_type N, T(nrows,sz);
54 
55  // resize the return matrices
56  A.derived().resize(nrows, sz);
57  b.derived().resize(nrows, dim);
58 
59  // build the least squares terms
61  for (i=0; i<npts; ++i)
62  {
63  col_type Tvec;
64  eli::geom::utility::bezier_T(Tvec, t[i], n);
65  T.row(i)=Tvec.transpose();
66  b.row(i)=pts[i];
67  }
68  for (i=npts; i<nrows; ++i)
69  {
70  T.row(i).setZero();
71  b.row(i).setZero();
72  }
73 
74  A=T*N;
75  }
76 
77  template<typename index_type1, typename index_type2, typename index_type3>
78  index_type1 determine_n(const index_type1 &deg_in, const index_type2 &nconstrs, const index_type3 &npts)
79  {
80  index_type1 n;
81 
82  // if have an under-determined system then what should be done? For now, reducing
83  // the order of resulting curve
84  if (deg_in+1>static_cast<index_type1>(npts+nconstrs))
85  {
86  n=npts+nconstrs-1;
87  std::cerr << "deg_in (" << deg_in << ") too low. Order should be " << n << std::endl;
88  assert(false);
89  }
90  else
91  n=deg_in;
92 
93  return n;
94  }
95  }
96  }
97  }
98 }
99 
100 namespace eli
101 {
102  namespace geom
103  {
104  namespace curve
105  {
106  // TODO: Integrate the tol__ class into this class to replace open_flag and any other place
107  // where numerical error might affect an equivalence comparison
108  template<typename data__, unsigned short dim__, typename tol__=eli::util::tolerance<data__> >
109  class bezier
110  {
111  public:
112  typedef unsigned short dimension_type;
113  typedef data__ data_type;
114  typedef Eigen::Matrix<data_type, 1, dim__> point_type;
115  typedef point_type control_point_type;
116  typedef typename point_type::Index index_type;
118  typedef tol__ tolerance_type;
119  typedef Eigen::Matrix<data_type, dim__, dim__> rotation_matrix_type;
121 
122  public:
123  bezier() : B(1, dim__) {}
124  bezier(const index_type &n) : B((n<=0)?(1):(n+1), dim__) {}
126  ~bezier() {}
127 
129  {
130  if (this != &bc)
131  {
132  B=bc.B;
133  }
134  return *this;
135  }
136 
138  {
139  if (this == &bc)
140  return true;
141  if ((B.rows()!=B.rows()) || (B.cols()!=B.cols()))
142  return false;
143  if (B!=bc.B)
144  return false;
145  return true;
146  }
147 
149  {
150  return !operator==(bc);
151  }
152 
154  {
155  tolerance_type tol;
156 
157  if (this==&bc)
158  return true;
159 
160  if ((B.rows()!=B.rows()) || (B.cols()!=B.cols()))
161  return false;
162 
163  for (index_type i=0; i<=degree(); ++i)
164  {
165  if (!tol.approximately_equal(get_control_point(i), bc.get_control_point(i)))
166  {
167  return false;
168  }
169  }
170 
171  return true;
172  }
173 
175  {
176  if (this==&bc)
177  return 0;
178 
179  data_type d;
181  return d;
182  }
183 
184  void clear() {resize(0);}
185 
186  void resize(const index_type &t_dim)
187  {
188  B.resize(t_dim+1, dim__);
189  }
190 
191  index_type degree() const
192  {
193  return B.rows()-1;
194  }
195 
196  data_type get_tmax() const {return 1;}
197  data_type get_t0() const {return 0;}
198 
199  static dimension_type dimension() {return dim__;}
200 
201  void set_control_point(const control_point_type &cp, const index_type &i)
202  {
203  // make sure have valid index
204  if (i>degree())
205  {
206  assert(false);
207  return;
208  }
209 
210  B.row(i)=cp;
211  }
212 
213  control_point_type get_control_point(const index_type &i) const
214  {
215  // make sure have valid index
216  if (i>degree())
217  {
218  assert(false);
219  return B.row(0);
220  }
221 
222  return B.row(i);
223  }
224 
225  void reflect_xy()
226  {
227  B.col(2)=-B.col(2);
228  }
229 
230  void reflect_xz()
231  {
232  B.col(1)=-B.col(1);
233  }
234 
235  void reflect_yz()
236  {
237  B.col(0)=-B.col(0);
238  }
239 
240  void reflect(const point_type &normal)
241  {
242  point_type n(normal);
243 
244  n.normalize();
245  B=B-2*(B*n.transpose())*n;
246  }
247 
248  void reflect(const point_type &normal, const data_type &d)
249  {
250  point_type n(normal);
251 
252  n.normalize();
253  B=B-2*(B*n.transpose()-d*Eigen::Matrix<data_type, Eigen::Dynamic, 1>::Ones(degree()+1, 1))*n;
254  }
255 
256  void reverse()
257  {
258  index_type i, n(degree());
259  control_point_matrix_type B_new(n+1, dim__);
260 
261  for (i=0; i<=n; ++i)
262  {
263  B_new.row(n-i)=B.row(i);
264  }
265 
266  // set the new control points
267  B=B_new;
268  }
269 
270  void get_bounding_box(bounding_box_type &bb) const
271  {
272  index_type i, deg(degree());
273 
274  bb.clear();
275  for (i=0; i<=deg; ++i)
276  {
277  bb.add(B.row(i));
278  }
279  }
280 
281  void rotate(const rotation_matrix_type &rmat)
282  {
283  B*=rmat.transpose();
284  }
285 
286  void rotate(const rotation_matrix_type &rmat, const point_type &rorig)
287  {
288  translate(-rorig);
289  rotate(rmat);
290  translate(rorig);
291  }
292 
293  void translate(const point_type &trans)
294  {
295  index_type i, deg(degree());
296  for (i=0; i<=deg; ++i)
297  {
298  B.row(i)+=trans;
299  }
300  }
301 
302  void scale(const data_type &s)
303  {
304  index_type i, deg(degree());
305  for (i=0; i<=deg; ++i)
306  {
307  B.row(i)*=s;
308  }
309  }
310 
311  bool open() const {return !closed();}
312  bool closed() const
313  {
314  tolerance_type tol;
315 
316  for (index_type i=0; i<dim__; ++i)
317  {
318  if (!tol.approximately_equal(B(0, i), B(degree(), i)))
319  return false;
320  }
321  return true;
322  }
323 
324  point_type f(const data_type &t) const
325  {
326  // check to make sure have valid curve
327  assert(degree()>=0);
328 
329  // check to make sure given valid parametric value
330  assert((t>=0) && (t<=1));
331 
332  // short circuit if degree not high enough
333  if (degree()==0)
334  {
335  return B.row(0);
336  }
337 
338  point_type rtn;
340 
341  return rtn;
342  }
343 
344  point_type fp(const data_type &t) const
345  {
346  // check to make sure have valid curve
347  assert(degree()>=0);
348 
349  // check to make sure given valid parametric value
350  assert((t>=0) && (t<=1));
351 
352  point_type rtn;
353 
354  // short circuit if degree not high enough
355  if (degree()<1)
356  {
357  rtn.setZero();
358  return rtn;
359  }
360 
361  control_point_matrix_type B_p(degree()+1-1, dim__);
362 
365 
366  return rtn;
367  }
368 
369  point_type fpp(const data_type &t) const
370  {
371  // check to make sure have valid curve
372  assert(degree()>=0);
373 
374  // check to make sure given valid parametric value
375  assert((t>=0) && (t<=1));
376 
377  point_type rtn;
378 
379  // short circuit if degree not high enough
380  if (degree()<2)
381  {
382  rtn.setZero();
383  return rtn;
384  }
385 
386  control_point_matrix_type B_pp(degree()+1-2, dim__);
387 
389  eli::geom::utility::de_casteljau(rtn, B_pp, t);
390 
391  return rtn;
392  }
393 
394  point_type fppp(const data_type &t) const
395  {
396  // check to make sure have valid curve
397  assert(degree()>=0);
398 
399  // check to make sure given valid parametric value
400  assert((t>=0) && (t<=1));
401 
402  point_type rtn;
403 
404  // short circuit if degree not high enough
405  if (degree()<3)
406  {
407  rtn.setZero();
408  return rtn;
409  }
410 
411  control_point_matrix_type B_ppp(degree()+1-3, dim__);
412 
414  eli::geom::utility::de_casteljau(rtn, B_ppp, t);
415 
416  return rtn;
417  }
418 
419  point_type tangent(const data_type &t) const
420  {
421  point_type tgt(fp(t));
422 
423  tgt.normalize();
424  return tgt;
425  }
426 
427  void frenet_serret_frame(point_type &t, point_type &n, point_type &b, const data_type &t0)
428  {
429  t=tangent(t0);
430  b=fp(t0).cross(fpp(t0));
431  n=-t.cross(b)/b.norm();
432  b.normalize();
433  n.normalize();
434  }
435 
437  {
438  // create vector of new control points
439  control_point_matrix_type B_new(degree()+2, dim__);
440 
441  // build the new control points
443 
444  // set the new control points
445  B=B_new;
446  }
447 
448  void degree_promote_to(const index_type target_degree)
449  {
450  // create vector of new control points
451  control_point_matrix_type B_new(target_degree+1, dim__);
452 
453  // build the new control points
455 
456  // set the new control points
457  B=B_new;
458  }
459 
461  {
462  // check if can demote
463  int ncon(0);
464  switch(continuity_degree)
465  {
467  ncon=0;
468  break;
470  ncon=2;
471  break;
473  ncon=4;
474  break;
476  ncon=6;
477  break;
478  default:
479  ncon=-1;
480  }
481  if (ncon<0)
482  return false;
483  if (ncon>degree()-2)
484  return false;
485 
486  // demote control points and set them
487  control_point_matrix_type B_new(degree(), dim__);
489  B=B_new;
490 
491  return true;
492  }
493 
495  {
496  // allocate control points and set them
497  control_point_matrix_type B_new(4, dim__);
499  B=B_new;
500  }
501 
502  void split(bezier<data_type, dim__> &bc_l, bezier<data_type, dim__> &bc_r, const data_type &t0) const
503  {
504  if ( (t0>1) || (t0<0) )
505  {
506  assert(false);
507  return;
508  }
509 
510  control_point_matrix_type bl(degree()+1, dim__), br(degree()+1, dim__);
511  index_type n(degree());
512 
513  // resize the curves
514  bc_l.resize(n);
515  bc_r.resize(n);
516 
518 
519  // set the control points
520  for (index_type i=0; i<=n; ++i)
521  {
522  bc_l.set_control_point(bl.row(i), i);
523  bc_r.set_control_point(br.row(i), i);
524  }
525  }
526 
527  void fit(const fit_container_type &fcon, const index_type &deg_in)
528  {
529  std::vector<data_type> t;
530  fit_only(t, fcon, deg_in);
531  }
532 
533  void fit(std::vector<data_type> &t, const fit_container_type &fcon, const index_type &deg_in)
534  {
535  fit_with_error(t, fcon, deg_in);
536  }
537 
538  data_type fit_with_error(const fit_container_type &fcon, const index_type &deg_in)
539  {
540  std::vector<data_type> t;
541  return fit_with_error(t, fcon, deg_in);
542  }
543 
544  data_type fit_with_error(std::vector<data_type> &t, const fit_container_type &fcon, const index_type &deg_in)
545  {
546  fit_only(t, fcon, deg_in);
547  return est_fit_error(t, fcon);
548  }
549 
550  void interpolate(const fit_container_type &fcon)
551  {
552  std::vector<data_type> t;
553  interpolate(t, fcon);
554  }
555 
556  void interpolate(std::vector<data_type> &t, const fit_container_type &fcon)
557  {
558  size_t i, npts(fcon.number_points()), n, nclosed, ai;
559  std::vector<point_type, Eigen::aligned_allocator<point_type> > pts(npts);
560 
561  // get the points from the container
562  fcon.get_points(pts.begin());
563 
564  // account for the closed constraints
565  switch(fcon.get_end_flag())
566  {
568  {
569  nclosed=3;
570  break;
571  }
573  {
574  nclosed=2;
575  break;
576  }
578  {
579  nclosed=1;
580  break;
581  }
582  default:
583  {
584  nclosed=0;
585  break;
586  }
587  }
588 
589  // calculate the t corresponding to input points via approximate arc-length
590  determine_t(t, pts, fcon.closed());
591 
592  // determine the actual degree of curve
593  n=fcon.number_constraints(false)+nclosed+npts-1;
594 
595  // build the fit terms from points
596  mat_type A, x;
597  row_pts_type b;
598  internal::build_fit_Ab(A, b, t, pts, n, dim__);
599 
600  size_t nconpts=fcon.number_constraint_points();
601  std::vector<typename fit_container_type::index_type> indexes(nconpts);
602 
603  // cycle through all of the constraints
604  fcon.get_constraint_indexes(indexes.begin());
605  ai=pts.size();
606  for (size_t i=0; i<nconpts; ++i)
607  {
608  point_type pt;
610  typename fit_container_type::error_code ec;
611 
612  ec=fcon.get_constraint(indexes[i], ci);
614  {
615  assert(false);
616  }
617  else
618  {
619  // set the C1 constraint
621  {
622  col_type Tp;
623  mat_type N;
625 
626  eli::geom::utility::bezier_T_p(Tp, t[indexes[i]], n);
627  A.row(ai)=Tp.transpose()*N;
628  b.row(ai)=ci.get_fp();
629  ++ai;
630  }
631 
632  // set the C2 constraint
634  {
635  col_type Tpp;
636  mat_type N;
638 
639  eli::geom::utility::bezier_T_pp(Tpp, t[indexes[i]], n);
640  A.row(ai)=Tpp.transpose()*N;
641  b.row(ai)=ci.get_fpp();
642  ++ai;
643  }
644  }
645  }
646 
647  // add the closed constraint if needed
648  switch(fcon.get_end_flag())
649  {
651  {
652  A(ai,2)=1.0;
653  A(ai,1)=-2.0;
654  A(ai,0)=1.0;
655  A(ai,n)=-1.0;
656  A(ai,n-1)=2.0;
657  A(ai,n-2)=-1.0;
658  b.row(ai).setZero();
659  ++ai;
660  }
662  {
663  A(ai,1)=1.0;
664  A(ai,0)=-1.0;
665  A(ai,n)=-1.0;
666  A(ai,n-1)=1.0;
667  ++ai;
668  }
670  {
671  A(ai,0)=1.0;
672  A(ai,n)=-1.0;
673  ++ai;
674  }
675  default:
676  {
677  // no need to do anything
678  assert(ai==n+1);
679  }
680  }
681 
682  // solve for the control points
683  x=A.lu().solve(b);
684 
685  // extract the control points and set them
686  control_point_matrix_type ctrl(n+1, dim__);
687  for (i=0; i<n+1; ++i)
688  ctrl.row(i)=x.row(i);
689 
690  // ensure that the last control point and first are the same
691  if (fcon.closed())
692  ctrl.row(n)=ctrl.row(0);
693 
694  B=ctrl;
695  }
696 
697  private:
698  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_point_matrix_type;
699  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> row_pts_type;
700  typedef Eigen::Matrix<data_type, Eigen::Dynamic, 1> col_type;
701  typedef Eigen::Matrix<data_type, 1, Eigen::Dynamic> row_type;
702  typedef Eigen::Matrix<data_type, Eigen::Dynamic, Eigen::Dynamic> mat_type;
703 
704  private:
705  control_point_matrix_type B;
707  private:
708  void determine_t(std::vector<data_type> &t, const std::vector<point_type, Eigen::aligned_allocator<point_type> > &pts, bool closed) const
709  {
710  index_type i, npts(pts.size());
711  data_type len, small_dist(std::numeric_limits<data_type>::epsilon());
712 
713  // calculate the corresponding t values for each point via approximate arc-length
714  t.resize(npts);
715  t[0]=0.0;
716  for (i=1; i<npts; ++i)
717  {
718  data_type temp;
719 
720  temp=geom::point::distance(pts[i-1], pts[i]);
721 
722  // in case two adjacent points are given too close together
723  if (temp<small_dist)
724  temp=small_dist;
725 
726  t[i]=t[i-1]+temp;
727  }
728  len=t[npts-1];
729  if (closed)
730  {
731  len+=geom::point::distance(pts[npts-1], pts[0]);
732  }
733 
734  for (i=0; i<npts; ++i)
735  t[i]/=len;
736  }
737 
738  // This method is private because the t value returned in the first parameter is not the result of a nearest
739  // neighbor fit as would be expected from the fit_with_error methods.
740  void fit_only(std::vector<data_type> &t, const fit_container_type &fcon, const index_type &deg_in)
741  {
742  size_t i, npts(fcon.number_points()), n, nclosed;
743  std::vector<point_type, Eigen::aligned_allocator<point_type> > pts(npts);
744 
745  // get the points from the container
746  fcon.get_points(pts.begin());
747 
748  // account for the closed constraints
749  switch(fcon.get_end_flag())
750  {
752  {
753  nclosed=3;
754  break;
755  }
757  {
758  nclosed=2;
759  break;
760  }
762  {
763  nclosed=1;
764  break;
765  }
766  default:
767  {
768  nclosed=0;
769  break;
770  }
771  }
772 
773  // calculate the t corresponding to input points via approximate arc-length
774  determine_t(t, pts, fcon.closed());
775 
776  // determine the actual degree of curve
777  n=internal::determine_n(deg_in, fcon.number_constraints()+nclosed, npts);
778 
779  // build the fit terms from points
780  mat_type A, x;
781  row_pts_type b;
782  internal::build_fit_Ab(A, b, t, pts, n, dim__);
783 
784  // handle special case of unconstrained optimization problem
785  if ((fcon.number_constraints()==0) && (fcon.open()))
786  {
787  // determine the coefficients
789  }
790  else
791  {
792  // now becomes a constrained least squares problem
793  // the constraints come from closed flag and/or constraint collection
794  size_t bi, ncon=fcon.number_constraints()+nclosed;
795 
796  mat_type B(ncon, n+1);
797  row_pts_type d(ncon, dim__);
798 
799  // construct the system of constraints
800  B.setZero();
801  d.setZero();
802 
803  size_t nconpts=fcon.number_constraint_points();
804  std::vector<typename fit_container_type::index_type> indexes(nconpts);
805 
806  // cycle through all of the constraints
807  fcon.get_constraint_indexes(indexes.begin());
808  bi=0;
809  for (size_t i=0; i<nconpts; ++i)
810  {
811  point_type pt;
813  typename fit_container_type::error_code ec;
814 
815  ec=fcon.get_constraint(indexes[i], ci);
817  {
818  assert(false);
819  }
820  else
821  {
822  // set the C0 constraint
823  col_type T;
824  mat_type N;
825  eli::geom::utility::bezier_T(T, t[indexes[i]], n);
827  B.row(bi)=T.transpose()*N;
828  d.row(bi)=pts[indexes[i]];
829  ++bi;
830 
831  // set the C1 constraint
833  {
834  col_type Tp;
835 
836  eli::geom::utility::bezier_T_p(Tp, t[indexes[i]], n);
837  B.row(bi)=Tp.transpose()*N;
838  d.row(bi)=ci.get_fp();
839  ++bi;
840  }
841 
842  // set the C2 constraint
844  {
845  col_type Tpp;
846 
847  eli::geom::utility::bezier_T_pp(Tpp, t[indexes[i]], n);
848  B.row(bi)=Tpp.transpose()*N;
849  d.row(bi)=ci.get_fpp();
850  ++bi;
851  }
852  }
853  }
854 
855  // add the closed constraint if needed
856  switch(fcon.get_end_flag())
857  {
859  {
860  B(bi,2)=1.0;
861  B(bi,1)=-2.0;
862  B(bi,0)=1.0;
863  B(bi,n)=-1.0;
864  B(bi,n-1)=2.0;
865  B(bi,n-2)=-1.0;
866  ++bi;
867  }
869  {
870  B(bi,1)=1.0;
871  B(bi,0)=-1.0;
872  B(bi,n)=-1.0;
873  B(bi,n-1)=1.0;
874  ++bi;
875  }
877  {
878  B(bi,0)=1.0;
879  B(bi,n)=-1.0;
880  ++bi;
881  }
882  default:
883  {
884  // no need to do anything
885  assert(bi==ncon);
886  }
887  }
888 
889  // determine the coefficients
891  }
892 
893  // extract the control points and set them
894  control_point_matrix_type ctrl(n+1, dim__);
895  for (i=0; i<n+1; ++i)
896  ctrl.row(i)=x.row(i);
897 
898  // ensure that the last control point and first are the same
899  if (fcon.closed())
900  ctrl.row(n)=ctrl.row(0);
901 
902  B=ctrl;
903  }
904 
905  data_type est_fit_error(std::vector<data_type> &t, const fit_container_type &fcon)
906  {
907  size_t i, npts(fcon.number_points());
908  std::vector<point_type, Eigen::aligned_allocator<point_type> > pts(npts);
909 
910  // get the points from the container
911  fcon.get_points(pts.begin());
912 
913  // calculate the error at the point
914  data_type err(0);
915  for (i=0; i<pts.size(); ++i)
916  {
917  err+=eli::geom::intersect::minimum_distance(t[i], *this, pts[i]);
918  }
919 
920  return err;
921  }
922 
923  };
924 
934  }
935  }
936 }
937 #endif
Definition: continuity.hpp:26
void get_bounding_box(bounding_box_type &bb) const
Definition: bezier.hpp:270
void bezier_split_control_points(Eigen::MatrixBase< Derived1 > &cp_lo, Eigen::MatrixBase< Derived1 > &cp_hi, const Eigen::MatrixBase< Derived2 > &cp_in, const typename Derived2::Scalar &t)
Definition: bezier.hpp:291
Definition: math.hpp:20
void clear()
Definition: bounding_box.hpp:106
bool open() const
Definition: bezier.hpp:311
void reflect_xy()
Definition: bezier.hpp:225
bezier< float, 3 > bezier3f
Definition: bezier.hpp:927
void split(bezier< data_type, dim__ > &bc_l, bezier< data_type, dim__ > &bc_r, const data_type &t0) const
Definition: bezier.hpp:502
Eigen::Matrix< data_type, dim__, dim__ > rotation_matrix_type
Definition: bezier.hpp:119
Definition: continuity.hpp:28
curve::piecewise< curve1__, data1__, dim1__, tol1__ >::data_type minimum_distance(typename curve::piecewise< curve1__, data1__, dim1__, tol1__ >::data_type &t, const curve::piecewise< curve1__, data1__, dim1__, tol1__ > &pc, const typename curve::piecewise< curve1__, data1__, dim1__, tol1__ >::point_type &pt)
Derived1__::Scalar distance(const Eigen::MatrixBase< Derived1__ > &p1, const Eigen::MatrixBase< Derived2__ > &p2)
Definition: distance.hpp:33
use_states using_fp() const
Definition: fit_container.hpp:88
void interpolate(std::vector< data_type > &t, const fit_container_type &fcon)
Definition: bezier.hpp:556
void reflect_xz()
Definition: bezier.hpp:230
bezier< float, 2 > bezier2f
Definition: bezier.hpp:926
Definition: bounding_box.hpp:27
bool operator!=(const bezier< data_type, dim__, tolerance_type > &bc) const
Definition: bezier.hpp:148
void least_squares_eqcon(Eigen::MatrixBase< data1__ > &x, const Eigen::MatrixBase< data2__ > &A, const Eigen::MatrixBase< data3__ > &b, const Eigen::MatrixBase< data4__ > &B, const Eigen::MatrixBase< data5__ > &d)
Definition: least_squares.hpp:38
eli::geom::general::bounding_box< data_type, dim__, tolerance_type > bounding_box_type
Definition: bezier.hpp:120
void rotate(const rotation_matrix_type &rmat)
Definition: bezier.hpp:281
bezier()
Definition: bezier.hpp:123
void bezier_control_points_to_cubic(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in)
Definition: bezier.hpp:142
Eigen::Matrix< data_type, 1, dim__ > point_type
Definition: bezier.hpp:114
void bezier_ppp_control_point(Eigen::MatrixBase< Derived1 > &cp_ppp, const Eigen::MatrixBase< Derived2 > &cp)
Definition: bezier.hpp:78
void least_squares_uncon(Eigen::MatrixBase< data1__ > &x, const Eigen::MatrixBase< data2__ > &A, const Eigen::MatrixBase< data3__ > &r)
Definition: least_squares.hpp:25
void degree_promote()
Definition: bezier.hpp:436
void fit_only(std::vector< data_type > &t, const fit_container_type &fcon, const index_type &deg_in)
Definition: bezier.hpp:740
point_type f(const data_type &t) const
Definition: bezier.hpp:324
size_t number_constraint_points() const
Definition: fit_container.hpp:375
Eigen::Matrix< data_type, 1, Eigen::Dynamic > row_type
Definition: bezier.hpp:701
bezier< long double, 3 > bezier3ld
Definition: bezier.hpp:933
bezier & operator=(const bezier< data_type, dim__, tolerance_type > &bc)
Definition: bezier.hpp:128
bool closed() const
Definition: bezier.hpp:312
size_t number_points() const
Definition: fit_container.hpp:374
tol__ tolerance_type
Definition: bezier.hpp:118
void bezier_T_p(Eigen::MatrixBase< Derived > &Tp, const typename Derived::Scalar &t, const typename Derived::Index &n)
Definition: bezier.hpp:356
Definition: fit_container.hpp:34
void bezier_eqp_distance_bound(const Eigen::MatrixBase< Derived1 > &cp_a, const Eigen::MatrixBase< Derived2 > &cp_b, typename Derived1::Scalar &maxd)
Definition: bezier.hpp:256
bezier< double, 1 > bezier1d
Definition: bezier.hpp:928
void degree_promote_to(const index_type target_degree)
Definition: bezier.hpp:448
void build_fit_Ab(Eigen::MatrixBase< Derived1 > &A, Eigen::MatrixBase< Derived2 > &b, std::vector< typename Derived1::Scalar > &t, const std::vector< PointType, Eigen::aligned_allocator< PointType > > &pts, const typename Derived1::Index &n, const size_t &dim)
Definition: bezier.hpp:43
void bezier_T(Eigen::MatrixBase< Derived > &T, const typename Derived::Scalar &t, const typename Derived::Index &n)
Definition: bezier.hpp:341
bool approximately_equal(const bezier< data_type, dim__, tolerance_type > &bc) const
Definition: bezier.hpp:153
data_type fit_with_error(std::vector< data_type > &t, const fit_container_type &fcon, const index_type &deg_in)
Definition: bezier.hpp:544
use_states using_fpp() const
Definition: fit_container.hpp:89
Eigen::Matrix< data_type, Eigen::Dynamic, dim__ > control_point_matrix_type
Definition: bezier.hpp:698
void de_casteljau(Eigen::MatrixBase< Derived1 > &p, const Eigen::MatrixBase< Derived2 > &cp, const typename Derived2::Scalar &t)
Definition: bezier.hpp:27
point_type tangent(const data_type &t) const
Definition: bezier.hpp:419
data_type fit_with_error(const fit_container_type &fcon, const index_type &deg_in)
Definition: bezier.hpp:538
void translate(const point_type &trans)
Definition: bezier.hpp:293
void frenet_serret_frame(point_type &t, point_type &n, point_type &b, const data_type &t0)
Definition: bezier.hpp:427
data_type get_tmax() const
Definition: bezier.hpp:196
geom::curve::fit_container< data_type, index_type, dim__, dim__ > fit_container_type
Definition: bezier.hpp:117
control_point_matrix_type B
Definition: bezier.hpp:705
index_type degree() const
Definition: bezier.hpp:191
void clear()
Definition: bezier.hpp:184
Definition: fit_container.hpp:43
~bezier()
Definition: bezier.hpp:126
error_code get_constraint(const index_type &i, constraint_info &ci) const
Definition: fit_container.hpp:476
point_type fpp(const data_type &t) const
Definition: bezier.hpp:369
control_point_type get_control_point(const index_type &i) const
Definition: bezier.hpp:213
data_type est_fit_error(std::vector< data_type > &t, const fit_container_type &fcon)
Definition: bezier.hpp:905
Definition: continuity.hpp:27
point_type get_fpp() const
Definition: fit_container.hpp:87
const eli::geom::general::continuity & get_end_flag() const
Definition: fit_container.hpp:397
void degree_to_cubic()
Definition: bezier.hpp:494
void bezier_promote_control_points(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in)
Definition: bezier.hpp:93
index_type1 determine_n(const index_type1 &deg_in, const index_type2 &nconstrs, const index_type3 &npts)
Definition: bezier.hpp:78
bool add(const point_type &p)
Definition: bounding_box.hpp:113
data_type get_t0() const
Definition: bezier.hpp:197
void interpolate(const fit_container_type &fcon)
Definition: bezier.hpp:550
void rotate(const rotation_matrix_type &rmat, const point_type &rorig)
Definition: bezier.hpp:286
bezier< long double, 1 > bezier1ld
Definition: bezier.hpp:931
bool operator==(const bezier< data_type, dim__, tolerance_type > &bc) const
Definition: bezier.hpp:137
void reflect(const point_type &normal, const data_type &d)
Definition: bezier.hpp:248
void fit(const fit_container_type &fcon, const index_type &deg_in)
Definition: bezier.hpp:527
void bezier_pp_control_point(Eigen::MatrixBase< Derived1 > &cp_pp, const Eigen::MatrixBase< Derived2 > &cp)
Definition: bezier.hpp:63
void bezier_demote_control_points(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in, int ncon)
Definition: bezier.hpp:179
point_type get_fp() const
Definition: fit_container.hpp:86
point_type control_point_type
Definition: bezier.hpp:115
void resize(const index_type &t_dim)
Definition: bezier.hpp:186
void set_control_point(const control_point_type &cp, const index_type &i)
Definition: bezier.hpp:201
void bezier_p_control_point(Eigen::MatrixBase< Derived1 > &cp_p, const Eigen::MatrixBase< Derived2 > &cp)
Definition: bezier.hpp:48
bezier< float, 1 > bezier1f
Definition: bezier.hpp:925
continuity
Definition: continuity.hpp:24
error_code
Definition: fit_container.hpp:41
bezier(const bezier< data_type, dim__, tolerance_type > &bc)
Definition: bezier.hpp:125
point_type fppp(const data_type &t) const
Definition: bezier.hpp:394
point_type fp(const data_type &t) const
Definition: bezier.hpp:344
bool open() const
Definition: fit_container.hpp:395
bezier< double, 2 > bezier2d
Definition: bezier.hpp:929
bezier< double, 3 > bezier3d
Definition: bezier.hpp:930
void determine_t(std::vector< data_type > &t, const std::vector< point_type, Eigen::aligned_allocator< point_type > > &pts, bool closed) const
Definition: bezier.hpp:708
Definition: fit_container.hpp:50
Eigen::Matrix< data_type, Eigen::Dynamic, 1 > col_type
Definition: bezier.hpp:700
data__ data_type
Definition: bezier.hpp:113
void get_points(it__ itb) const
Definition: fit_container.hpp:451
Eigen::Matrix< data_type, Eigen::Dynamic, dim__ > row_pts_type
Definition: bezier.hpp:699
void bezier_promote_control_points_to(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in)
Definition: bezier.hpp:108
Eigen::Matrix< data_type, Eigen::Dynamic, Eigen::Dynamic > mat_type
Definition: bezier.hpp:702
unsigned short dimension_type
Definition: bezier.hpp:112
Definition: bezier.hpp:109
bezier< long double, 2 > bezier2ld
Definition: bezier.hpp:932
void bezier_N(Eigen::MatrixBase< Derived > &N, const typename Derived::Index &n)
Definition: bezier.hpp:316
void reverse()
Definition: bezier.hpp:256
size_t number_constraints(bool fit=true) const
Definition: fit_container.hpp:377
error_code get_constraint_indexes(it__ itout) const
Definition: fit_container.hpp:688
data_type eqp_distance_bound(const bezier< data_type, dim__, tolerance_type > &bc) const
Definition: bezier.hpp:174
void scale(const data_type &s)
Definition: bezier.hpp:302
void fit(std::vector< data_type > &t, const fit_container_type &fcon, const index_type &deg_in)
Definition: bezier.hpp:533
bool closed() const
Definition: fit_container.hpp:396
bezier(const index_type &n)
Definition: bezier.hpp:124
Definition: continuity.hpp:29
void reflect_yz()
Definition: bezier.hpp:235
void bezier_T_pp(Eigen::MatrixBase< Derived > &Tpp, const typename Derived::Scalar &t, const typename Derived::Index &n)
Definition: bezier.hpp:374
bool degree_demote(const geom::general::continuity &continuity_degree=geom::general::C0)
Definition: bezier.hpp:460
point_type::Index index_type
Definition: bezier.hpp:116
static dimension_type dimension()
Definition: bezier.hpp:199
void reflect(const point_type &normal)
Definition: bezier.hpp:240