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 
131  // NOTHING
132  }
133 
134  /*Operations*/
138  virtual point_t operator()(const time_t t) const {
139  check_conditions();
140  if (Safe & !(T_min_ <= t && t <= T_max_)) {
141  throw std::invalid_argument("can't evaluate bezier curve, time t is out of range"); // TODO
142  }
143  if (size_ == 1) {
144  return mult_T_ * control_points_[0];
145  } else {
146  return evalHorner(t);
147  }
148  }
149 
158  bool isApprox(const bezier_curve_t& other, const Numeric prec = Eigen::NumTraits<Numeric>::dummy_precision()) const {
159  bool equal = ndcurves::isApprox<num_t>(T_min_, other.min()) && ndcurves::isApprox<num_t>(T_max_, other.max()) &&
160  dim_ == other.dim() && degree_ == other.degree() && size_ == other.size_ &&
161  ndcurves::isApprox<Numeric>(mult_T_, other.mult_T_) && bernstein_ == other.bernstein_;
162  if (!equal) return false;
163  for (size_t i = 0; i < size_; ++i) {
164  if (!control_points_.at(i).isApprox(other.control_points_.at(i), prec)) return false;
165  }
166  return true;
167  }
168 
169  virtual bool isApprox(const curve_abc_t* other,
170  const Numeric prec = Eigen::NumTraits<Numeric>::dummy_precision()) const {
171  const bezier_curve_t* other_cast = dynamic_cast<const bezier_curve_t*>(other);
172  if (other_cast)
173  return isApprox(*other_cast, prec);
174  else
175  return false;
176  }
177 
178  virtual bool operator==(const bezier_curve_t& other) const { return isApprox(other); }
179 
180  virtual bool operator!=(const bezier_curve_t& other) const { return !(*this == other); }
181 
186  bezier_curve_t compute_derivate(const std::size_t order) const {
187  check_conditions();
188  if (order == 0) {
189  return *this;
190  }
191  t_point_t derived_wp;
192  for (typename t_point_t::const_iterator pit = control_points_.begin(); pit != control_points_.end() - 1; ++pit) {
193  derived_wp.push_back((num_t)degree_ * (*(pit + 1) - (*pit)));
194  }
195  if (derived_wp.empty()) {
196  derived_wp.push_back(point_t::Zero(dim_));
197  }
198  bezier_curve_t deriv(derived_wp.begin(), derived_wp.end(), T_min_, T_max_, mult_T_ * (1. / (T_max_ - T_min_)));
199  return deriv.compute_derivate(order - 1);
200  }
201 
205  bezier_curve_t* compute_derivate_ptr(const std::size_t order) const {
206  return new bezier_curve_t(compute_derivate(order));
207  }
208 
214  bezier_curve_t compute_primitive(const std::size_t order) const {
215  check_conditions();
216  if (order == 0) {
217  return *this;
218  }
219  num_t new_degree_inv = 1. / ((num_t)(degree_ + 1));
220  t_point_t n_wp;
221  point_t current_sum = point_t::Zero(dim_);
222  // recomputing waypoints q_i from derivative waypoints p_i. q_0 is the given constant.
223  // then q_i = (sum( j = 0 -> j = i-1) p_j) /n+1
224  n_wp.push_back(current_sum);
225  for (typename t_point_t::const_iterator pit = control_points_.begin(); pit != control_points_.end(); ++pit) {
226  current_sum += *pit;
227  n_wp.push_back(current_sum * new_degree_inv);
228  }
229  bezier_curve_t integ(n_wp.begin(), n_wp.end(), T_min_, T_max_, mult_T_ * (T_max_ - T_min_));
230  return integ.compute_primitive(order - 1);
231  }
232 
237  bezier_curve_t elevate(const std::size_t order) const {
238  t_point_t new_waypoints = control_points_, temp_waypoints;
239  for (std::size_t i = 1; i<= order; ++i)
240  {
241  num_t new_degree_inv = 1. / ((num_t)(degree_ + i));
242  temp_waypoints.push_back(*new_waypoints.begin());
243  num_t idx_deg_inv = 0.;
244  for (typename t_point_t::const_iterator pit = new_waypoints.begin()+1; pit != new_waypoints.end(); ++pit) {
245  idx_deg_inv += new_degree_inv;
246  temp_waypoints.push_back(idx_deg_inv * (*(pit-1)) + (1 - idx_deg_inv) * (*pit));
247  }
248  temp_waypoints.push_back(*(new_waypoints.end()-1));
249  new_waypoints = temp_waypoints;
250  temp_waypoints.clear();
251  }
252  return bezier_curve_t (new_waypoints.begin(), new_waypoints.end(), T_min_, T_max_, mult_T_);
253  }
254 
258  void elevate_self(const std::size_t order) {
259  if (order > 0)
260  (*this) = elevate(order);
261  }
262 
270  virtual point_t derivate(const time_t t, const std::size_t order) const { return compute_derivate(order)(t); }
271 
280  point_t evalBernstein(const Numeric t) const {
281  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
282  point_t res = point_t::Zero(dim_);
283  typename t_point_t::const_iterator control_points_it = control_points_.begin();
284  for (typename std::vector<Bern<Numeric> >::const_iterator cit = bernstein_.begin(); cit != bernstein_.end();
285  ++cit, ++control_points_it) {
286  res += cit->operator()(u) * (*control_points_it);
287  }
288  return res * mult_T_;
289  }
290 
303  point_t evalHorner(const Numeric t) const {
304  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
305  typename t_point_t::const_iterator control_points_it = control_points_.begin();
306  Numeric u_op, bc, tn;
307  u_op = 1.0 - u;
308  bc = 1;
309  tn = 1;
310  point_t tmp = (*control_points_it) * u_op;
311  ++control_points_it;
312  for (unsigned int i = 1; i < degree_; i++, ++control_points_it) {
313  tn = tn * u;
314  bc = bc * ((num_t)(degree_ - i + 1)) / i;
315  tmp = (tmp + tn * bc * (*control_points_it)) * u_op;
316  }
317  return (tmp + tn * u * (*control_points_it)) * mult_T_;
318  }
319 
320  const t_point_t& waypoints() const { return control_points_; }
321 
322  const point_t waypointAtIndex(const std::size_t index) const {
323  point_t waypoint;
324  if (index < control_points_.size()) {
325  waypoint = control_points_[index];
326  }
327  return waypoint;
328  }
329 
336  point_t evalDeCasteljau(const Numeric t) const {
337  // normalize time :
338  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
339  t_point_t pts = deCasteljauReduction(waypoints(), u);
340  while (pts.size() > 1) {
341  pts = deCasteljauReduction(pts, u);
342  }
343  return pts[0] * mult_T_;
344  }
345 
346  t_point_t deCasteljauReduction(const Numeric t) const {
347  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
348  return deCasteljauReduction(waypoints(), u);
349  }
350 
360  t_point_t deCasteljauReduction(const t_point_t& pts, const Numeric u) const {
361  if (u < 0 || u > 1) {
362  throw std::out_of_range("In deCasteljau reduction : u is not in [0;1]");
363  }
364  if (pts.size() == 1) {
365  return pts;
366  }
367 
368  t_point_t new_pts;
369  for (cit_point_t cit = pts.begin(); cit != (pts.end() - 1); ++cit) {
370  new_pts.push_back((1 - u) * (*cit) + u * (*(cit + 1)));
371  }
372  return new_pts;
373  }
374 
379  std::pair<bezier_curve_t, bezier_curve_t> split(const Numeric t) const {
380  check_conditions();
381  if (fabs(t - T_max_) < MARGIN) {
382  throw std::runtime_error("can't split curve, interval range is equal to original curve");
383  }
384  t_point_t wps_first(size_), wps_second(size_);
385  const Numeric u = (t - T_min_) / (T_max_ - T_min_);
386  t_point_t casteljau_pts = waypoints();
387  wps_first[0] = casteljau_pts.front();
388  wps_second[degree_] = casteljau_pts.back();
389  size_t id = 1;
390  while (casteljau_pts.size() > 1) {
391  casteljau_pts = deCasteljauReduction(casteljau_pts, u);
392  wps_first[id] = casteljau_pts.front();
393  wps_second[degree_ - id] = casteljau_pts.back();
394  ++id;
395  }
396  bezier_curve_t c_first(wps_first.begin(), wps_first.end(), T_min_, t, mult_T_);
397  bezier_curve_t c_second(wps_second.begin(), wps_second.end(), t, T_max_, mult_T_);
398  return std::make_pair(c_first, c_second);
399  }
400 
406  piecewise_curve_t split(const vector_x_t& times) const {
407  std::vector<bezier_curve_t> curves;
408  bezier_curve_t current = *this;
409  for (int i = 0; i < times.rows(); ++i) {
410  std::pair<bezier_curve_t, bezier_curve_t> pairsplit = current.split(times[i]);
411  curves.push_back(pairsplit.first);
412  current = pairsplit.second;
413  }
414  curves.push_back(current);
415  piecewise_curve_t res;
416  for (typename std::vector<bezier_curve_t>::const_iterator cit = curves.begin(); cit != curves.end(); ++cit) {
417  typename piecewise_curve_t::curve_ptr_t ptr(new bezier_curve_t(*cit));
418  res.add_curve_ptr(ptr);
419  }
420  return res;
421  }
422 
429  bezier_curve_t extract(const Numeric t1, const Numeric t2) {
430  if (t1 < T_min_ || t1 > T_max_ || t2 < T_min_ || t2 > T_max_) {
431  throw std::out_of_range("In Extract curve : times out of bounds");
432  }
433  if (fabs(t1 - T_min_) < MARGIN && fabs(t2 - T_max_) < MARGIN) // t1=T_min and t2=T_max
434  {
435  return bezier_curve_t(waypoints().begin(), waypoints().end(), T_min_, T_max_, mult_T_);
436  }
437  if (fabs(t1 - T_min_) < MARGIN) // t1=T_min
438  {
439  return split(t2).first;
440  }
441  if (fabs(t2 - T_max_) < MARGIN) // t2=T_max
442  {
443  return split(t1).second;
444  }
445  std::pair<bezier_curve_t, bezier_curve_t> c_split = this->split(t1);
446  return c_split.second.split(t2).first;
447  }
448 
456  bezier_curve_t cross(const bezier_curve_t& g) const {
457  //See Farouki and Rajan 1988 Alogirthms for polynomials in Bernstein form and
458  //http://web.mit.edu/hyperbook/Patrikalakis-Maekawa-Cho/node10.html
459  assert_operator_compatible(g);
460  if (dim()!= 3)
461  throw std::invalid_argument("Can't perform cross product on Bezier curves with dimensions != 3 ");
462  int m =(int)(degree());
463  int n =(int)(g.degree());
464  unsigned int mj, n_ij, mn_i;
465  t_point_t new_waypoints;
466  for(int i = 0; i<= m+n; ++i)
467  {
468  bezier_curve_t::point_t current_point = bezier_curve_t::point_t::Zero(dim());
469  for (int j = std::max(0,i-n); j <=std::min(m,i); ++j){
470  mj = bin(m,j);
471  n_ij = bin(n,i-j);
472  mn_i = bin(m+n,i);
473  num_t mul = num_t(mj*n_ij) / num_t(mn_i);
474  current_point += mul*ndcurves::cross(waypointAtIndex(j), g.waypointAtIndex(i-j));
475  }
476  new_waypoints.push_back(current_point);
477  }
478  return bezier_curve_t(new_waypoints.begin(),new_waypoints.end(),min(),max(),mult_T_ * g.mult_T_);
479  }
480 
481 
488  bezier_curve_t cross(const bezier_curve_t::point_t& point) const {
489  //See Farouki and Rajan 1988 Alogirthms for polynomials in Bernstein form and
490  //http://web.mit.edu/hyperbook/Patrikalakis-Maekawa-Cho/node10.html
491  if (dim()!= 3)
492  throw std::invalid_argument("Can't perform cross product on Bezier curves with dimensions != 3 ");
493  t_point_t new_waypoints;
494  for(typename t_point_t::const_iterator cit = waypoints().begin(); cit != waypoints().end(); ++cit){
495  new_waypoints.push_back(ndcurves::cross(*cit, point));
496  }
497  return bezier_curve_t(new_waypoints.begin(),new_waypoints.end(),min(),max(),mult_T_);
498  }
499 
500  bezier_curve_t& operator+=(const bezier_curve_t& other) {
501  assert_operator_compatible(other);
502  bezier_curve_t other_elevated = other * (other.mult_T_ / this->mult_T_); // TODO remove mult_T_ from Bezier
503  if(other.degree() > degree()){
504  elevate_self(other.degree() - degree());
505  }
506  else if(other_elevated.degree() < degree()){
507  other_elevated.elevate_self(degree() - other_elevated.degree());
508  }
509  typename t_point_t::const_iterator otherit = other_elevated.control_points_.begin();
510  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it, ++otherit){
511  (*it)+=(*otherit);
512  }
513  return *this;
514  }
515 
516  bezier_curve_t& operator-=(const bezier_curve_t& other) {
517  assert_operator_compatible(other);
518  bezier_curve_t other_elevated = other * (other.mult_T_ / this->mult_T_);
519  if(other.degree() > degree()){
520  elevate_self(other.degree() - degree());
521  }
522  else if(other_elevated.degree() < degree()){
523  other_elevated.elevate_self(degree() - other_elevated.degree());
524  }
525  typename t_point_t::const_iterator otherit = other_elevated.control_points_.begin();
526  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it, ++otherit){
527  (*it)-=(*otherit);
528  }
529  return *this;
530  }
531 
532  bezier_curve_t& operator+=(const bezier_curve_t::point_t& point) {
533  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
534  (*it)+=point;
535  }
536  return *this;
537  }
538 
539  bezier_curve_t& operator-=(const bezier_curve_t::point_t& point) {
540  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
541  (*it)-=point;
542  }
543  return *this;
544  }
545 
546  bezier_curve_t& operator/=(const double d) {
547  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
548  (*it)/=d;
549  }
550  return *this;
551  }
552 
553  bezier_curve_t& operator*=(const double d) {
554  for (typename t_point_t::iterator it = control_points_.begin(); it!=control_points_.end(); ++it){
555  (*it)*=d;
556  }
557  return *this;
558  }
559 
560 
561  private:
566  template <typename In>
567  t_point_t add_constraints(In PointsBegin, In PointsEnd, const curve_constraints_t& constraints) {
568  t_point_t res;
569  num_t T = T_max_ - T_min_;
570  num_t T_square = T * T;
571  point_t P0, P1, P2, P_n_2, P_n_1, PN;
572  P0 = *PointsBegin;
573  PN = *(PointsEnd - 1);
574  P1 = P0 + constraints.init_vel * T / (num_t)degree_;
575  P_n_1 = PN - constraints.end_vel * T / (num_t)degree_;
576  P2 = constraints.init_acc * T_square / (num_t)(degree_ * (degree_ - 1)) + 2 * P1 - P0;
577  P_n_2 = constraints.end_acc * T_square / (num_t)(degree_ * (degree_ - 1)) + 2 * P_n_1 - PN;
578  res.push_back(P0);
579  res.push_back(P1);
580  res.push_back(P2);
581  for (In it = PointsBegin + 1; it != PointsEnd - 1; ++it) {
582  res.push_back(*it);
583  }
584  res.push_back(P_n_2);
585  res.push_back(P_n_1);
586  res.push_back(PN);
587  return res;
588  }
589 
590  void check_conditions() const {
591  if (control_points_.size() == 0) {
592  throw std::runtime_error(
593  "Error in bezier curve : there is no control points set / did you use empty constructor ?");
594  } else if (dim_ == 0) {
595  throw std::runtime_error(
596  "Error in bezier curve : Dimension of points is zero / did you use empty constructor ?");
597  }
598  }
599 
600 
601  void assert_operator_compatible(const bezier_curve_t& other) const{
602  if ((fabs(min() - other.min()) > MARGIN) || (fabs(max() - other.max()) > MARGIN)){
603  throw std::invalid_argument("Can't perform base operation (+ - ) on two Bezier curves with different time ranges");
604  }
605  }
606 
607  /*Operations*/
608 
609  public:
610  /*Helpers*/
613  std::size_t virtual dim() const { return dim_; };
616  virtual time_t min() const { return T_min_; }
619  virtual time_t max() const { return T_max_; }
622  virtual std::size_t degree() const { return degree_; }
623  /*Helpers*/
624 
625  /* Attributes */
627  std::size_t dim_;
629  /*const*/ time_t T_min_;
631  /*const*/ time_t T_max_;
632  /*const*/ time_t mult_T_;
633  /*const*/ std::size_t size_;
634  /*const*/ std::size_t degree_;
635  /*const*/ std::vector<Bern<Numeric> > bernstein_;
636  /*const*/ t_point_t control_points_;
637  /* Attributes */
638 
639  static bezier_curve_t zero(const std::size_t dim, const time_t T = 1.) {
640  std::vector<point_t> ts;
641  ts.push_back(point_t::Zero(dim));
642  return bezier_curve_t(ts.begin(), ts.end(), 0., T);
643  }
644 
645  // Serialization of the class
646  friend class boost::serialization::access;
647 
648  template <class Archive>
649  void serialize(Archive& ar, const unsigned int version) {
650  if (version) {
651  // Do something depending on version ?
652  }
653  ar& BOOST_SERIALIZATION_BASE_OBJECT_NVP(curve_abc_t);
654  ar& boost::serialization::make_nvp("dim", dim_);
655  ar& boost::serialization::make_nvp("T_min", T_min_);
656  ar& boost::serialization::make_nvp("T_max", T_max_);
657  ar& boost::serialization::make_nvp("mult_T", mult_T_);
658  ar& boost::serialization::make_nvp("size", size_);
659  ar& boost::serialization::make_nvp("degree", degree_);
660  ar& boost::serialization::make_nvp("bernstein", bernstein_);
661  ar& boost::serialization::make_nvp("control_points", control_points_);
662  }
663 }; // End struct bezier_curve
664 
665 template <typename T, typename N, bool S, typename P >
667  bezier_curve<T,N,S,P> res(p1);
668  return res+=p2;
669 }
670 
671 template <typename T, typename N, bool S, typename P >
673  std::vector<typename bezier_curve<T,N,S,P>::point_t> ts;
674  for (std::size_t i = 0; i <= p1.degree(); ++i){
675  ts.push_back(bezier_curve<T,N,S,P>::point_t::Zero(p1.dim()));
676  }
677  bezier_curve<T,N,S,P> res (ts.begin(),ts.end(),p1.min(),p1.max());
678  res-=p1;
679  return res;
680 }
681 
682 template <typename T, typename N, bool S, typename P >
684  bezier_curve<T,N,S,P> res(p1);
685  return res-=p2;
686 }
687 
688 
689 template <typename T, typename N, bool S, typename P >
691  bezier_curve<T,N,S,P> res(p1);
692  return res-=point;
693 }
694 
695 template <typename T, typename N, bool S, typename P >
697  bezier_curve<T,N,S,P> res(-p1);
698  return res+=point;
699 }
700 
701 template <typename T, typename N, bool S, typename P >
703  bezier_curve<T,N,S,P> res(p1);
704  return res+=point;
705 }
706 
707 template <typename T, typename N, bool S, typename P >
709  bezier_curve<T,N,S,P> res(p1);
710  return res+=point;
711 }
712 
713 template <typename T, typename N, bool S, typename P >
715  bezier_curve<T,N,S,P> res(p1);
716  return res/=k;
717 }
718 
719 template <typename T, typename N, bool S, typename P >
721  bezier_curve<T,N,S,P> res(p1);
722  return res*=k;
723 }
724 
725 template <typename T, typename N, bool S, typename P >
727  bezier_curve<T,N,S,P> res(p1);
728  return res*=k;
729 }
730 
731 } // namespace ndcurves
732 
733 DEFINE_CLASS_TEMPLATE_VERSION(SINGLE_ARG(typename Time, typename Numeric, bool Safe, typename Point),
735 
736 #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:270
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:237
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:186
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:158
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:180
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:205
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:672
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:138
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:280
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
Compute the primitive of the curve at order N. Computes the primitive at order N of bezier curve of p...
Definition: bezier_curve.h:214
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:178
double Time
Definition: effector_spline.h:27
~bezier_curve()
Destructor.
Definition: bezier_curve.h:130
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:169
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:258
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:714
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:720
point_t end_vel
Definition: curve_constraint.h:65
Definition: curve_constraint.h:22
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:666