bezier_curve.h
Go to the documentation of this file.
1 
9 #ifndef _CLASS_BEZIERCURVE
10 #define _CLASS_BEZIERCURVE
11 
12 #include "curve_abc.h"
13 #include "cross_implementation.h"
14 #include "bernstein.h"
15 #include "curve_constraint.h"
16 #include "piecewise_curve.h"
17 
18 #include "MathDefs.h"
19 
20 #include <vector>
21 #include <stdexcept>
22 
23 #include <iostream>
24 
25 namespace ndcurves {
31 template <typename Time = double, typename Numeric = Time, bool Safe = false,
32  typename Point = Eigen::Matrix<Numeric, Eigen::Dynamic, 1> >
33 struct bezier_curve : public curve_abc<Time, Numeric, Safe, Point> {
34  typedef Point point_t;
35  typedef Eigen::Matrix<Numeric, Eigen::Dynamic, 1> vector_x_t;
36  typedef Eigen::Ref<const vector_x_t> vector_x_ref_t;
37  typedef Time time_t;
38  typedef Numeric num_t;
40  typedef std::vector<point_t, Eigen::aligned_allocator<point_t> > t_point_t;
41  typedef typename t_point_t::const_iterator cit_point_t;
43  typedef boost::shared_ptr<bezier_curve_t> bezier_curve_ptr_t;
47 
48  /* Constructors - destructors */
49  public:
52  bezier_curve() : dim_(0), T_min_(0), T_max_(0) {}
53 
62  template <typename In>
63  bezier_curve(In PointsBegin, In PointsEnd, const time_t T_min = 0., const time_t T_max = 1.,
64  const time_t mult_T = 1.)
65  : dim_(PointsBegin->size()),
66  T_min_(T_min),
67  T_max_(T_max),
68  mult_T_(mult_T),
69  size_(std::distance(PointsBegin, PointsEnd)),
70  degree_(size_ - 1),
71  bernstein_(ndcurves::makeBernstein<num_t>((unsigned int)degree_))
72  {
73  if (bernstein_.size() != size_) {
74  throw std::invalid_argument("Invalid size of polynomial");
75  }
76  In it(PointsBegin);
77  if (Safe && (size_ < 1 || T_max_ <= T_min_)) {
78  throw std::invalid_argument("can't create bezier min bound is higher than max bound");
79  }
80  for (; it != PointsEnd; ++it) {
81  if(Safe && static_cast<size_t>(it->size()) != dim_)
82  throw std::invalid_argument("All the control points must have the same dimension.");
83  control_points_.push_back(*it);
84  }
85  }
86 
97  template <typename In>
98  bezier_curve(In PointsBegin, In PointsEnd, const curve_constraints_t& constraints, const time_t T_min = 0.,
99  const time_t T_max = 1., const time_t mult_T = 1.)
100  : dim_(PointsBegin->size()),
101  T_min_(T_min),
102  T_max_(T_max),
103  mult_T_(mult_T),
104  size_(std::distance(PointsBegin, PointsEnd) + 4),
105  degree_(size_ - 1),
106  bernstein_(ndcurves::makeBernstein<num_t>((unsigned int)degree_))
107  {
108  if (Safe && (size_ < 1 || T_max_ <= T_min_)) {
109  throw std::invalid_argument("can't create bezier min bound is higher than max bound");
110  }
111  t_point_t updatedList = add_constraints<In>(PointsBegin, PointsEnd, constraints);
112  for (cit_point_t cit = updatedList.begin(); cit != updatedList.end(); ++cit) {
113  if(Safe && static_cast<size_t>(cit->size()) != dim_)
114  throw std::invalid_argument("All the control points must have the same dimension.");
115  control_points_.push_back(*cit);
116  }
117  }
118 
120  : dim_(other.dim_),
121  T_min_(other.T_min_),
122  T_max_(other.T_max_),
123  mult_T_(other.mult_T_),
124  size_(other.size_),
125  degree_(other.degree_),
126  bernstein_(other.bernstein_),
127  control_points_(other.control_points_) {}
128 
130  virtual ~bezier_curve() {}
131 
132  /*Operations*/
136  virtual point_t operator()(const time_t t) const {
137  check_conditions();
138  if (Safe & !(T_min_ <= t && t <= T_max_)) {
139  throw std::invalid_argument("can't evaluate bezier curve, time t is out of range"); // TODO
140  }
141  if (size_ == 1) {
142  return mult_T_ * control_points_[0];
143  } else {
144  return evalHorner(t);
145  }
146  }
147 
156  bool isApprox(const bezier_curve_t& other, const Numeric prec = Eigen::NumTraits<Numeric>::dummy_precision()) const {
157  bool equal = ndcurves::isApprox<num_t>(T_min_, other.min()) && ndcurves::isApprox<num_t>(T_max_, other.max()) &&
158  dim_ == other.dim() && degree_ == other.degree() && size_ == other.size_ &&
159  ndcurves::isApprox<Numeric>(mult_T_, other.mult_T_) && bernstein_ == other.bernstein_;
160  if (!equal) return false;
161  for (size_t i = 0; i < size_; ++i) {
162  if (!control_points_.at(i).isApprox(other.control_points_.at(i), prec)) return false;
163  }
164  return true;
165  }
166 
167  virtual bool isApprox(const curve_abc_t* other,
168  const Numeric prec = Eigen::NumTraits<Numeric>::dummy_precision()) const {
169  const bezier_curve_t* other_cast = dynamic_cast<const bezier_curve_t*>(other);
170  if (other_cast)
171  return isApprox(*other_cast, prec);
172  else
173  return false;
174  }
175 
176  virtual bool operator==(const bezier_curve_t& other) const { return isApprox(other); }
177 
178  virtual bool operator!=(const bezier_curve_t& other) const { return !(*this == other); }
179 
184  bezier_curve_t compute_derivate(const std::size_t order) const {
185  check_conditions();
186  if (order == 0) {
187  return *this;
188  }
189  t_point_t derived_wp;
190  for (typename t_point_t::const_iterator pit = control_points_.begin(); pit != control_points_.end() - 1; ++pit) {
191  derived_wp.push_back((num_t)degree_ * (*(pit + 1) - (*pit)));
192  }
193  if (derived_wp.empty()) {
194  derived_wp.push_back(point_t::Zero(dim_));
195  }
196  bezier_curve_t deriv(derived_wp.begin(), derived_wp.end(), T_min_, T_max_, mult_T_ * (1. / (T_max_ - T_min_)));
197  return deriv.compute_derivate(order - 1);
198  }
199 
203  bezier_curve_t* compute_derivate_ptr(const std::size_t order) const {
204  return new bezier_curve_t(compute_derivate(order));
205  }
206 
213  bezier_curve_t compute_primitive(const std::size_t order, const point_t& init) const {
214  check_conditions();
215  if (order == 0) {
216  return *this;
217  }
218  num_t new_degree_inv = 1. / ((num_t)(degree_ + 1));
219  t_point_t n_wp;
220  point_t current_sum (init);
221  n_wp.push_back(current_sum);
222  for (typename t_point_t::const_iterator pit = control_points_.begin(); pit != control_points_.end(); ++pit) {
223  current_sum += *pit;
224  n_wp.push_back(current_sum * new_degree_inv);
225  }
226  bezier_curve_t integ(n_wp.begin(), n_wp.end(), T_min_, T_max_, mult_T_ * (T_max_ - T_min_));
227  return integ.compute_primitive(order - 1);
228  }
229 
230  bezier_curve_t compute_primitive(const std::size_t order) const {
231  return compute_primitive(order, point_t::Zero(dim_));
232  }
233 
234  bezier_curve_t* compute_primitive_ptr(const std::size_t order, const point_t& init) const {
235  return new bezier_curve_t(compute_primitive(order, init));
236  }
237 
242  bezier_curve_t elevate(const std::size_t order) const {
243  t_point_t new_waypoints = control_points_, temp_waypoints;
244  for (std::size_t i = 1; i<= order; ++i)
245  {
246  num_t new_degree_inv = 1. / ((num_t)(degree_ + i));
247  temp_waypoints.push_back(*new_waypoints.begin());
248  num_t idx_deg_inv = 0.;
249  for (typename t_point_t::const_iterator pit = new_waypoints.begin()+1; pit != new_waypoints.end(); ++pit) {
250  idx_deg_inv += new_degree_inv;
251  temp_waypoints.push_back(idx_deg_inv * (*(pit-1)) + (1 - idx_deg_inv) * (*pit));
252  }
253  temp_waypoints.push_back(*(new_waypoints.end()-1));
254  new_waypoints = temp_waypoints;
255  temp_waypoints.clear();
256  }
257  return bezier_curve_t (new_waypoints.begin(), new_waypoints.end(), T_min_, T_max_, mult_T_);
258  }
259 
263  void elevate_self(const std::size_t order) {
264  if (order > 0)
265  (*this) = elevate(order);
266  }
267 
275  virtual point_t derivate(const time_t t, const std::size_t order) const { return compute_derivate(order)(t); }
276 
285  point_t evalBernstein(const Numeric t) const {
286  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
287  point_t res = point_t::Zero(dim_);
288  typename t_point_t::const_iterator control_points_it = control_points_.begin();
289  for (typename std::vector<Bern<Numeric> >::const_iterator cit = bernstein_.begin(); cit != bernstein_.end();
290  ++cit, ++control_points_it) {
291  res += cit->operator()(u) * (*control_points_it);
292  }
293  return res * mult_T_;
294  }
295 
308  point_t evalHorner(const Numeric t) const {
309  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
310  typename t_point_t::const_iterator control_points_it = control_points_.begin();
311  Numeric u_op, bc, tn;
312  u_op = 1.0 - u;
313  bc = 1;
314  tn = 1;
315  point_t tmp = (*control_points_it) * u_op;
316  ++control_points_it;
317  for (unsigned int i = 1; i < degree_; i++, ++control_points_it) {
318  tn = tn * u;
319  bc = bc * ((num_t)(degree_ - i + 1)) / i;
320  tmp = (tmp + tn * bc * (*control_points_it)) * u_op;
321  }
322  return (tmp + tn * u * (*control_points_it)) * mult_T_;
323  }
324 
325  const t_point_t& waypoints() const { return control_points_; }
326 
327  const point_t waypointAtIndex(const std::size_t index) const {
328  point_t waypoint;
329  if (index < control_points_.size()) {
330  waypoint = control_points_[index];
331  }
332  return waypoint;
333  }
334 
341  point_t evalDeCasteljau(const Numeric t) const {
342  // normalize time :
343  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
344  t_point_t pts = deCasteljauReduction(waypoints(), u);
345  while (pts.size() > 1) {
346  pts = deCasteljauReduction(pts, u);
347  }
348  return pts[0] * mult_T_;
349  }
350 
351  t_point_t deCasteljauReduction(const Numeric t) const {
352  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
353  return deCasteljauReduction(waypoints(), u);
354  }
355 
365  t_point_t deCasteljauReduction(const t_point_t& pts, const Numeric u) const {
366  if (u < 0 || u > 1) {
367  throw std::out_of_range("In deCasteljau reduction : u is not in [0;1]");
368  }
369  if (pts.size() == 1) {
370  return pts;
371  }
372 
373  t_point_t new_pts;
374  for (cit_point_t cit = pts.begin(); cit != (pts.end() - 1); ++cit) {
375  new_pts.push_back((1 - u) * (*cit) + u * (*(cit + 1)));
376  }
377  return new_pts;
378  }
379 
384  std::pair<bezier_curve_t, bezier_curve_t> split(const Numeric t) const {
385  check_conditions();
386  if (fabs(t - T_max_) < MARGIN) {
387  throw std::runtime_error("can't split curve, interval range is equal to original curve");
388  }
389  t_point_t wps_first(size_), wps_second(size_);
390  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
391  t_point_t casteljau_pts = waypoints();
392  wps_first[0] = casteljau_pts.front();
393  wps_second[degree_] = casteljau_pts.back();
394  size_t id = 1;
395  while (casteljau_pts.size() > 1) {
396  casteljau_pts = deCasteljauReduction(casteljau_pts, u);
397  wps_first[id] = casteljau_pts.front();
398  wps_second[degree_ - id] = casteljau_pts.back();
399  ++id;
400  }
401  bezier_curve_t c_first(wps_first.begin(), wps_first.end(), T_min_, t, mult_T_);
402  bezier_curve_t c_second(wps_second.begin(), wps_second.end(), t, T_max_, mult_T_);
403  return std::make_pair(c_first, c_second);
404  }
405 
411  piecewise_curve_t split(const vector_x_t& times) const {
412  std::vector<bezier_curve_t> curves;
413  bezier_curve_t current = *this;
414  for (int i = 0; i < times.rows(); ++i) {
415  std::pair<bezier_curve_t, bezier_curve_t> pairsplit = current.split(times[i]);
416  curves.push_back(pairsplit.first);
417  current = pairsplit.second;
418  }
419  curves.push_back(current);
420  piecewise_curve_t res;
421  for (typename std::vector<bezier_curve_t>::const_iterator cit = curves.begin(); cit != curves.end(); ++cit) {
422  typename piecewise_curve_t::curve_ptr_t ptr(new bezier_curve_t(*cit));
423  res.add_curve_ptr(ptr);
424  }
425  return res;
426  }
427 
434  bezier_curve_t extract(const Numeric t1, const Numeric t2) {
435  if (t1 < T_min_ || t1 > T_max_ || t2 < T_min_ || t2 > T_max_) {
436  throw std::out_of_range("In Extract curve : times out of bounds");
437  }
438  if (fabs(t1 - T_min_) < MARGIN && fabs(t2 - T_max_) < MARGIN) // t1=T_min and t2=T_max
439  {
440  return bezier_curve_t(waypoints().begin(), waypoints().end(), T_min_, T_max_, mult_T_);
441  }
442  if (fabs(t1 - T_min_) < MARGIN) // t1=T_min
443  {
444  return split(t2).first;
445  }
446  if (fabs(t2 - T_max_) < MARGIN) // t2=T_max
447  {
448  return split(t1).second;
449  }
450  std::pair<bezier_curve_t, bezier_curve_t> c_split = this->split(t1);
451  return c_split.second.split(t2).first;
452  }
453 
461  bezier_curve_t cross(const bezier_curve_t& g) const {
462  //See Farouki and Rajan 1988 Alogirthms for polynomials in Bernstein form and
463  //http://web.mit.edu/hyperbook/Patrikalakis-Maekawa-Cho/node10.html
464  assert_operator_compatible(g);
465  if (dim()!= 3)
466  throw std::invalid_argument("Can't perform cross product on Bezier curves with dimensions != 3 ");
467  int m =(int)(degree());
468  int n =(int)(g.degree());
469  unsigned int mj, n_ij, mn_i;
470  t_point_t new_waypoints;
471  for(int i = 0; i<= m+n; ++i)
472  {
473  bezier_curve_t::point_t current_point = bezier_curve_t::point_t::Zero(dim());
474  for (int j = std::max(0,i-n); j <=std::min(m,i); ++j){
475  mj = bin(m,j);
476  n_ij = bin(n,i-j);
477  mn_i = bin(m+n,i);
478  num_t mul = num_t(mj*n_ij) / num_t(mn_i);
479  current_point += mul*ndcurves::cross(waypointAtIndex(j), g.waypointAtIndex(i-j));
480  }
481  new_waypoints.push_back(current_point);
482  }
483  return bezier_curve_t(new_waypoints.begin(),new_waypoints.end(),min(),max(),mult_T_ * g.mult_T_);
484  }
485 
486 
493  bezier_curve_t cross(const bezier_curve_t::point_t& point) const {
494  //See Farouki and Rajan 1988 Alogirthms for polynomials in Bernstein form and
495  //http://web.mit.edu/hyperbook/Patrikalakis-Maekawa-Cho/node10.html
496  if (dim()!= 3)
497  throw std::invalid_argument("Can't perform cross product on Bezier curves with dimensions != 3 ");
498  t_point_t new_waypoints;
499  for(typename t_point_t::const_iterator cit = waypoints().begin(); cit != waypoints().end(); ++cit){
500  new_waypoints.push_back(ndcurves::cross(*cit, point));
501  }
502  return bezier_curve_t(new_waypoints.begin(),new_waypoints.end(),min(),max(),mult_T_);
503  }
504 
505  bezier_curve_t& operator+=(const bezier_curve_t& other) {
506  assert_operator_compatible(other);
507  bezier_curve_t other_elevated = other * (other.mult_T_ / this->mult_T_); // TODO remove mult_T_ from Bezier
508  if(other.degree() > degree()){
509  elevate_self(other.degree() - degree());
510  }
511  else if(other_elevated.degree() < degree()){
512  other_elevated.elevate_self(degree() - other_elevated.degree());
513  }
514  typename t_point_t::const_iterator otherit = other_elevated.control_points_.begin();
515  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it, ++otherit){
516  (*it)+=(*otherit);
517  }
518  return *this;
519  }
520 
521  bezier_curve_t& operator-=(const bezier_curve_t& other) {
522  assert_operator_compatible(other);
523  bezier_curve_t other_elevated = other * (other.mult_T_ / this->mult_T_);
524  if(other.degree() > degree()){
525  elevate_self(other.degree() - degree());
526  }
527  else if(other_elevated.degree() < degree()){
528  other_elevated.elevate_self(degree() - other_elevated.degree());
529  }
530  typename t_point_t::const_iterator otherit = other_elevated.control_points_.begin();
531  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it, ++otherit){
532  (*it)-=(*otherit);
533  }
534  return *this;
535  }
536 
537  bezier_curve_t& operator+=(const bezier_curve_t::point_t& point) {
538  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
539  (*it)+=point;
540  }
541  return *this;
542  }
543 
544  bezier_curve_t& operator-=(const bezier_curve_t::point_t& point) {
545  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
546  (*it)-=point;
547  }
548  return *this;
549  }
550 
551  bezier_curve_t& operator/=(const double d) {
552  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
553  (*it)/=d;
554  }
555  return *this;
556  }
557 
558  bezier_curve_t& operator*=(const double d) {
559  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
560  (*it)*=d;
561  }
562  return *this;
563  }
564 
565 
566  private:
571  template <typename In>
572  t_point_t add_constraints(In PointsBegin, In PointsEnd, const curve_constraints_t& constraints) {
573  t_point_t res;
574  num_t T = T_max_ - T_min_;
575  num_t T_square = T * T;
576  point_t P0, P1, P2, P_n_2, P_n_1, PN;
577  P0 = *PointsBegin;
578  PN = *(PointsEnd - 1);
579  P1 = P0 + constraints.init_vel * T / (num_t)degree_;
580  P_n_1 = PN - constraints.end_vel * T / (num_t)degree_;
581  P2 = constraints.init_acc * T_square / (num_t)(degree_ * (degree_ - 1)) + 2 * P1 - P0;
582  P_n_2 = constraints.end_acc * T_square / (num_t)(degree_ * (degree_ - 1)) + 2 * P_n_1 - PN;
583  res.push_back(P0);
584  res.push_back(P1);
585  res.push_back(P2);
586  for (In it = PointsBegin + 1; it != PointsEnd - 1; ++it) {
587  res.push_back(*it);
588  }
589  res.push_back(P_n_2);
590  res.push_back(P_n_1);
591  res.push_back(PN);
592  return res;
593  }
594 
595  void check_conditions() const {
596  if (control_points_.size() == 0) {
597  throw std::runtime_error(
598  "Error in bezier curve : there is no control points set / did you use empty constructor ?");
599  } else if (dim_ == 0) {
600  throw std::runtime_error(
601  "Error in bezier curve : Dimension of points is zero / did you use empty constructor ?");
602  }
603  }
604 
605 
606  void assert_operator_compatible(const bezier_curve_t& other) const{
607  if ((fabs(min() - other.min()) > MARGIN) || (fabs(max() - other.max()) > MARGIN)){
608  throw std::invalid_argument("Can't perform base operation (+ - ) on two Bezier curves with different time ranges");
609  }
610  }
611 
612  /*Operations*/
613 
614  public:
615  /*Helpers*/
618  std::size_t virtual dim() const { return dim_; };
621  virtual time_t min() const { return T_min_; }
624  virtual time_t max() const { return T_max_; }
627  virtual std::size_t degree() const { return degree_; }
628  /*Helpers*/
629 
630  /* Attributes */
632  std::size_t dim_;
634  /*const*/ time_t T_min_;
636  /*const*/ time_t T_max_;
637  /*const*/ time_t mult_T_;
638  /*const*/ std::size_t size_;
639  /*const*/ std::size_t degree_;
640  /*const*/ std::vector<Bern<Numeric> > bernstein_;
641  /*const*/ t_point_t control_points_;
642  /* Attributes */
643 
644  static bezier_curve_t zero(const std::size_t dim, const time_t T = 1.) {
645  std::vector<point_t> ts;
646  ts.push_back(point_t::Zero(dim));
647  return bezier_curve_t(ts.begin(), ts.end(), 0., T);
648  }
649 
650  // Serialization of the class
651  friend class boost::serialization::access;
652 
653  template <class Archive>
654  void serialize(Archive& ar, const unsigned int version) {
655  if (version) {
656  // Do something depending on version ?
657  }
658  ar& BOOST_SERIALIZATION_BASE_OBJECT_NVP(curve_abc_t);
659  ar& boost::serialization::make_nvp("dim", dim_);
660  ar& boost::serialization::make_nvp("T_min", T_min_);
661  ar& boost::serialization::make_nvp("T_max", T_max_);
662  ar& boost::serialization::make_nvp("mult_T", mult_T_);
663  ar& boost::serialization::make_nvp("size", size_);
664  ar& boost::serialization::make_nvp("degree", degree_);
665  ar& boost::serialization::make_nvp("bernstein", bernstein_);
666  ar& boost::serialization::make_nvp("control_points", control_points_);
667  }
668 }; // End struct bezier_curve
669 
670 template <typename T, typename N, bool S, typename P >
672  bezier_curve<T,N,S,P> res(p1);
673  return res+=p2;
674 }
675 
676 template <typename T, typename N, bool S, typename P >
678  std::vector<typename bezier_curve<T,N,S,P>::point_t> ts;
679  for (std::size_t i = 0; i <= p1.degree(); ++i){
680  ts.push_back(bezier_curve<T,N,S,P>::point_t::Zero(p1.dim()));
681  }
682  bezier_curve<T,N,S,P> res (ts.begin(),ts.end(),p1.min(),p1.max());
683  res-=p1;
684  return res;
685 }
686 
687 template <typename T, typename N, bool S, typename P >
689  bezier_curve<T,N,S,P> res(p1);
690  return res-=p2;
691 }
692 
693 
694 template <typename T, typename N, bool S, typename P >
696  bezier_curve<T,N,S,P> res(p1);
697  return res-=point;
698 }
699 
700 template <typename T, typename N, bool S, typename P >
702  bezier_curve<T,N,S,P> res(-p1);
703  return res+=point;
704 }
705 
706 template <typename T, typename N, bool S, typename P >
708  bezier_curve<T,N,S,P> res(p1);
709  return res+=point;
710 }
711 
712 template <typename T, typename N, bool S, typename P >
714  bezier_curve<T,N,S,P> res(p1);
715  return res+=point;
716 }
717 
718 template <typename T, typename N, bool S, typename P >
720  bezier_curve<T,N,S,P> res(p1);
721  return res/=k;
722 }
723 
724 template <typename T, typename N, bool S, typename P >
726  bezier_curve<T,N,S,P> res(p1);
727  return res*=k;
728 }
729 
730 template <typename T, typename N, bool S, typename P >
732  bezier_curve<T,N,S,P> res(p1);
733  return res*=k;
734 }
735 
736 } // namespace ndcurves
737 
738 DEFINE_CLASS_TEMPLATE_VERSION(SINGLE_ARG(typename Time, typename Numeric, bool Safe, typename Point),
740 
741 #endif //_CLASS_BEZIERCURVE
Time time_t
Definition: bezier_curve.h:37
Definition: bernstein.h:20
virtual point_t derivate(const time_t t, const std::size_t order) const
Evaluate the derivative order N of curve at time t. If derivative is to be evaluated several times...
Definition: bezier_curve.h:275
bezier_curve_t elevate(const std::size_t order) const
Computes a Bezier curve of order degrees higher than the current curve, but strictly equivalent...
Definition: bezier_curve.h:242
bezier_curve_t compute_derivate(const std::size_t order) const
Compute the derived curve at order N. Computes the derivative order N, of bezier curve of parametric...
Definition: bezier_curve.h:184
curve_abc_t::curve_ptr_t curve_ptr_t
Definition: bezier_curve.h:46
boost::shared_ptr< curve_t > curve_ptr_t
Definition: piecewise_curve.h:40
bool isApprox(const bezier_curve_t &other, const Numeric prec=Eigen::NumTraits< Numeric >::dummy_precision()) const
isApprox check if other and *this are approximately equals. Only two curves of the same class can be ...
Definition: bezier_curve.h:156
bezier_curve_t compute_primitive(const std::size_t order, const point_t &init) const
Compute the primitive of the curve at order N. Computes the primitive at order N of bezier curve of p...
Definition: bezier_curve.h:213
Eigen::Vector3d cross(const Eigen::VectorXd &a, const Eigen::VectorXd &b)
Definition: cross_implementation.h:15
virtual bool operator!=(const bezier_curve_t &other) const
Definition: bezier_curve.h:178
bezier_curve()
Empty constructor. Curve obtained this way can not perform other class functions. ...
Definition: bezier_curve.h:52
bezier_curve_t * compute_derivate_ptr(const std::size_t order) const
Compute the derived curve at order N.
Definition: bezier_curve.h:203
curve_constraints< point_t > curve_constraints_t
Definition: bezier_curve.h:39
t_point_t::const_iterator cit_point_t
Definition: bezier_curve.h:41
virtual time_t max() const =0
Get the maximum time for which the curve is defined.
class allowing to create a cubic hermite spline of any dimension.
interface for a Curve of arbitrary dimension.
void add_curve_ptr(const curve_ptr_t &cf)
Add a new curve to piecewise curve, which should be defined in where is equal to of the actual pie...
Definition: piecewise_curve.h:153
bezier_curve< T, N, S, P > operator-(const bezier_curve< T, N, S, P > &p1)
Definition: bezier_curve.h:677
bezier_curve(In PointsBegin, In PointsEnd, const curve_constraints_t &constraints, const time_t T_min=0., const time_t T_max=1., const time_t mult_T=1.)
Constructor with constraints. This constructor will add 4 points (2 after the first one...
Definition: bezier_curve.h:98
Eigen::Matrix< Numeric, Eigen::Dynamic, 1 > vector_x_t
Definition: bezier_curve.h:35
virtual point_t operator()(const time_t t) const
Evaluation of the bezier curve at time t.
Definition: bezier_curve.h:136
bezier_curve(In PointsBegin, In PointsEnd, const time_t T_min=0., const time_t T_max=1., const time_t mult_T=1.)
Constructor. Given the first and last point of a control points set, create the bezier curve...
Definition: bezier_curve.h:63
brief Computes all Bernstein polynomes for a certain degree std::vector< Bern< Numeric > > makeBernstein(const unsigned int n)
Definition: bernstein.h:88
Numeric num_t
Definition: bezier_curve.h:38
boost::shared_ptr< curve_t > curve_ptr_t
Definition: curve_abc.h:41
point_t evalBernstein(const Numeric t) const
Evaluate all Bernstein polynomes for a certain degree. A bezier curve with N control points is repres...
Definition: bezier_curve.h:285
std::vector< point_t, Eigen::aligned_allocator< point_t > > t_point_t
Definition: bezier_curve.h:40
boost::shared_ptr< bezier_curve_t > bezier_curve_ptr_t
Definition: bezier_curve.h:43
void serialize(Archive &ar, const unsigned int version)
Definition: curve_abc.h:139
bezier_curve_t compute_primitive(const std::size_t order) const
Definition: bezier_curve.h:230
virtual std::size_t degree() const =0
Get the degree of the curve.
Definition: fwd.h:49
virtual bool operator==(const bezier_curve_t &other) const
Definition: bezier_curve.h:176
double Time
Definition: effector_spline.h:27
std::vector< bezier_curve< Numeric, Numeric, true, linear_variable< Numeric > > > split(const problem_definition< Point, Numeric > &pDef, problem_data< Point, Numeric > &pData)
Definition: details.h:208
Eigen::Ref< const vector_x_t > vector_x_ref_t
Definition: bezier_curve.h:36
virtual bool isApprox(const curve_abc_t *other, const Numeric prec=Eigen::NumTraits< Numeric >::dummy_precision()) const
Definition: bezier_curve.h:167
bezier_curve< Time, Numeric, Safe, Point > bezier_curve_t
Definition: bezier_curve.h:42
void elevate_self(const std::size_t order)
Elevate the Bezier curve of order degrees higher than the current curve, but strictly equivalent...
Definition: bezier_curve.h:263
virtual time_t min() const =0
Get the minimum time for which the curve is defined.
point_t init_acc
Definition: curve_constraint.h:63
point_t end_acc
Definition: curve_constraint.h:66
virtual std::size_t dim() const =0
Get dimension of curve.
Eigen::Matrix< Numeric, Eigen::Dynamic, 1 > Point
Definition: effector_spline.h:28
Definition: bezier_curve.h:33
Point point_t
Definition: bezier_curve.h:34
struct to define constraints on start / end velocities and acceleration on a curve ...
piecewise_curve< Time, Numeric, Safe, point_t, point_t, bezier_curve_t > piecewise_curve_t
Definition: bezier_curve.h:44
double Numeric
Definition: effector_spline.h:26
point_t init_vel
Definition: curve_constraint.h:62
class allowing to create a piecewise curve.
bezier_curve< T, N, S, P > operator/(const bezier_curve< T, N, S, P > &p1, const double k)
Definition: bezier_curve.h:719
bezier_curve(const bezier_curve &other)
Definition: bezier_curve.h:119
Definition: fwd.h:34
bezier_curve< T, N, S, P > operator*(const bezier_curve< T, N, S, P > &p1, const double k)
Definition: bezier_curve.h:725
point_t end_vel
Definition: curve_constraint.h:65
Definition: curve_constraint.h:22
virtual ~bezier_curve()
Destructor.
Definition: bezier_curve.h:130
Represents a curve of dimension Dim. If value of parameter Safe is false, no verification is made on ...
Definition: curve_abc.h:34
curve_abc< Time, Numeric, Safe, point_t > curve_abc_t
Definition: bezier_curve.h:45
bezier_curve< T, N, S, P > operator+(const bezier_curve< T, N, S, P > &p1, const bezier_curve< T, N, S, P > &p2)
Definition: bezier_curve.h:671
bezier_curve_t * compute_primitive_ptr(const std::size_t order, const point_t &init) const
Definition: bezier_curve.h:234