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_surface_bezier_hpp
14 #define eli_geom_surface_bezier_hpp
15 
16 #include <cstddef> // nullptr
17 
18 #include <vector>
19 
20 #include "eli/code_eli.hpp"
21 
22 #include "eli/util/tolerance.hpp"
23 
30 
31 namespace eli
32 {
33  namespace geom
34  {
35  namespace surface
36  {
37  template<typename data__, unsigned short dim__, typename tol__=eli::util::tolerance<data__> >
38  class bezier
39  {
40  public:
41  typedef unsigned short dimension_type;
42  typedef data__ data_type;
43  typedef Eigen::Matrix<data_type, 1, dim__> point_type;
44  typedef point_type control_point_type;
45  typedef typename control_point_type::Index index_type;
46  typedef tol__ tolerance_type;
47  typedef Eigen::Matrix<data_type, dim__, dim__> rotation_matrix_type;
50 
51  private:
52  typedef Eigen::Map<Eigen::Matrix<data_type, Eigen::Dynamic, dim__>,
53  Eigen::Unaligned,
54  Eigen::Stride<1, dim__> > control_point_matrix_type;
55  typedef Eigen::Map<Eigen::Matrix<data_type, Eigen::Dynamic, dim__>,
56  Eigen::Unaligned,
57  Eigen::Stride<1, Eigen::Dynamic> > v_dir_control_point_matrix_type;
58  typedef std::vector<data_type> control_point_container;
59  typedef std::vector<control_point_matrix_type> u_control_point_matrix_container;
60  typedef std::vector<v_dir_control_point_matrix_type> v_control_point_matrix_container;
61 
62  private:
63  control_point_container point_data;
64  u_control_point_matrix_container B_u;
65  v_control_point_matrix_container B_v;
67  public:
68  bezier() : point_data(3, 0)
69  {
70  // set the B_u and B_v maps
71  set_Bs(1, 1);
72  }
73  bezier(const index_type &u_dim, const index_type &v_dim) : point_data(dim__*(u_dim+1)*(v_dim+1))
74  {
75  // set the B_u and B_v maps
76  set_Bs(u_dim, v_dim);
77  }
78  bezier(const bezier<data_type, dim__, tol__> &bs) : point_data(bs.point_data)
79  {
80  // set the B_u and B_v maps
81  set_Bs(bs.degree_u(), bs.degree_v());
82  }
83  ~bezier() {}
84 
86  {
87  if (this!=&bs)
88  {
89  point_data=bs.point_data;
90  set_Bs(bs.degree_u(), bs.degree_v());
91  }
92 
93  return (*this);
94  }
95 
97  {
98  if (this==&bs)
99  return true;
100 
101  if (point_data!=bs.point_data)
102  return false;
103 
104  if (B_u.size()!=bs.B_u.size())
105  return false;
106 
107  if (B_v.size()!=bs.B_v.size())
108  return false;
109 
110  return true;
111  }
112 
113  bool abouteq(const bezier<data_type, dim__, tol__> &bs, const data_type &ttol2 ) const
114  {
115  if (this==&bs)
116  return true;
117 
118  if (B_u.size()!=bs.B_u.size())
119  return false;
120 
121  if (B_v.size()!=bs.B_v.size())
122  return false;
123 
124  index_type i, j, degu(degree_u()), degv(degree_v());
125  for (j=0; j<=degv; ++j)
126  {
127  for (i=0; i<=degu; ++i)
128  {
129  if ( eli::geom::point::distance2( B_u[j].row(i), bs.B_u[j].row(i) ) > ttol2 )
130  {
131  return false;
132  }
133  }
134  }
135 
136  return true;
137  }
138 
140  {
141  return !operator==(bs);
142  }
143 
144  static dimension_type dimension() {return dim__;}
145 
146  index_type degree_u() const {return static_cast<index_type>(B_v.size())-1;}
147  index_type degree_v() const {return static_cast<index_type>(B_u.size())-1;}
148 
149  void get_parameter_min(data_type &umin, data_type &vmin) const
150  {
151  umin=0;
152  vmin=0;
153  }
154 
155  void get_parameter_max(data_type &umax, data_type &vmax) const
156  {
157  umax=1;
158  vmax=1;
159  }
160 
161  data_type get_umin() const {return static_cast<data_type>(0);}
162  data_type get_vmin() const {return static_cast<data_type>(0);}
163  data_type get_umax() const {return static_cast<data_type>(1);}
164  data_type get_vmax() const {return static_cast<data_type>(1);}
165 
166  void resize(const index_type &u_dim, const index_type &v_dim)
167  {
168  // allocate the control points
169  point_data.resize(dim__*(u_dim+1)*(v_dim+1));
170 
171  // set the B_u and B_v maps
172  set_Bs(u_dim, v_dim);
173  }
174 
175  point_type get_control_point(const index_type &i, const index_type &j) const
176  {
177  // make sure have valid indexes
178  if ( (i>degree_u()) || (j>degree_v()) )
179  {
180  assert(false);
181  return B_u[0].row(0);
182  }
183 
184  // make sure that the two maps point to same items
185  assert(B_u[j].row(i)==B_v[i].row(j));
186 
187  return B_u[j].row(i);
188  }
189 
190  void get_bounding_box(bounding_box_type &bb) const
191  {
192  index_type i, j, degu(degree_u()), degv(degree_v());
193 
194  bb.clear();
195  for (i=0; i<=degu; ++i)
196  {
197  for (j=0; j<=degv; ++j)
198  {
199  bb.add(B_u[j].row(i));
200  }
201  }
202  }
203 
204  void rotate(const rotation_matrix_type &rmat)
205  {
206  index_type j, degv(degree_v());
207  for (j=0; j<=degv; ++j)
208  {
209  B_u[j]*=rmat.transpose();
210  }
211  }
212 
213  void rotate(const rotation_matrix_type &rmat, const point_type &rorig)
214  {
215  translate(-rorig);
216  rotate(rmat);
217  translate(rorig);
218  }
219 
220  void translate(const point_type &trans)
221  {
222  index_type i, j, degu(degree_u()), degv(degree_v());
223  for (j=0; j<=degv; ++j)
224  {
225  for (i=0; i<=degu; ++i)
226  {
227  B_u[j].row(i)+=trans;
228  }
229  }
230  }
231 
232  void scale(const data_type &s)
233  {
234  index_type i, j, degu(degree_u()), degv(degree_v());
235  for (j=0; j<=degv; ++j)
236  {
237  for (i=0; i<=degu; ++i)
238  {
239  B_u[j].row(i)*=s;
240  }
241  }
242  }
243 
244  bool open_u() const {return !closed_u();}
245  bool closed_u() const
246  {
247  curve_type bc0, bc1;
248 
249  get_uconst_curve(bc0, 0);
250  get_uconst_curve(bc1, 1);
251 
252  return eli::geom::curve::equivalent_curves(bc0, bc1);
253  }
254  bool open_v() const {return !closed_v();}
255  bool closed_v() const
256  {
257  curve_type bc0, bc1;
258 
259  get_vconst_curve(bc0, 0);
260  get_vconst_curve(bc1, 1);
261 
262  return eli::geom::curve::equivalent_curves(bc0, bc1);
263  }
264 
265  void set_control_point(const point_type &cp, const index_type &i, const index_type &j)
266  {
267  // make sure have valid indexes
268  if ( (i>degree_u()) || (j>degree_v()) )
269  {
270  assert(false);
271  return;
272  }
273 
274  B_u[j].row(i) = cp;
275  }
276 
277  void reverse_u()
278  {
279  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_col_type;
280  typedef std::vector<control_col_type, Eigen::aligned_allocator<control_col_type> > control_col_collection_type;
281 
282  index_type i, n(degree_u()), m(degree_v());
283  control_col_collection_type current_col(n+1, control_col_type(m+1, dim__));
284 
285  // copy the control cols
286  for (i=0; i<=n; ++i)
287  current_col[i]=B_v[i];
288 
289  // set the new control points
290  for (i=0; i<=n; ++i)
291  {
292  B_v[i]=current_col[n-i];
293  }
294  }
295 
296  void reverse_v()
297  {
298  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
299  typedef std::vector<control_row_type, Eigen::aligned_allocator<control_row_type> > control_row_collection_type;
300 
301  index_type i, n(degree_u()), m(degree_v());
302  control_row_collection_type current_row(m+1, control_row_type(n+1, dim__));
303 
304  // copy the control rows
305  for (i=0; i<=m; ++i)
306  current_row[i]=B_u[i];
307 
308  // set the new control points
309  for (i=0; i<=m; ++i)
310  {
311  B_u[i]=current_row[m-i];
312  }
313  }
314 
315  void swap_uv()
316  {
317  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
318  typedef std::vector<control_row_type, Eigen::aligned_allocator<control_row_type> > control_row_collection_type;
319 
320  index_type i, j, n(degree_u()), m(degree_v());
321  control_row_collection_type current_row(m+1, control_row_type(n+1, dim__));
322 
323  // copy the control rows
324  for (i=0; i<=m; ++i)
325  current_row[i]=B_u[i];
326 
327  // resize current surface
328  resize(m, n);
329 
330  // set the new control points
331  for (i=0; i<=m; ++i)
332  {
333  for (j=0; j<=n; ++j)
334  {
335  set_control_point(current_row[i].row(j), i, j);
336  }
337  }
338  }
339 
340  void get_uconst_curve(curve_type &bc, const data_type &u) const
341  {
342  index_type j, m(degree_v());
343 
344  // check to make sure have valid curve
345  assert(m>=0);
346 
347  // check to make sure given valid parametric value
348  assert((u>=0) && (u<=1));
349 
350  // build curve
351  point_type cp;
352 
353  bc.resize(m);
354  for (j=0; j<=m; ++j)
355  {
356  eli::geom::utility::de_casteljau(cp, B_u[j], u);
357  bc.set_control_point(cp, j);
358  }
359  }
360 
361  void get_vconst_curve(curve_type &bc, const data_type &v) const
362  {
363  index_type i, n(degree_u());
364 
365  // check to make sure have valid curve
366  assert(n>=0);
367 
368  // check to make sure given valid parametric value
369  assert((v>=0) && (v<=1));
370 
371  // build curve
372  point_type cp;
373 
374  bc.resize(n);
375  for (i=0; i<=n; ++i)
376  {
377  eli::geom::utility::de_casteljau(cp, B_v[i], v);
378  bc.set_control_point(cp, i);
379  }
380  }
381 
382  point_type f(const data_type &u, const data_type &v) const
383  {
384  point_type ans, tmp;
385  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp;
386  index_type i, n(degree_u()), m(degree_v());
387 
388  // check to make sure have valid curve
389  assert(n>=0);
390  assert(m>=0);
391 
392  // check to make sure given valid parametric value
393  assert((u>=0) && (u<=1));
394  assert((v>=0) && (v<=1));
395 
396  if (n<=m)
397  {
398  temp_cp.resize(m+1, dim__);
399  // build the temporary control points
400  for (i=0; i<=m; ++i)
401  {
402  eli::geom::utility::de_casteljau(tmp, B_u[i], u);
403  temp_cp.row(i)=tmp;
404  }
405  eli::geom::utility::de_casteljau(ans, temp_cp, v);
406  }
407  else
408  {
409  temp_cp.resize(n+1, dim__);
410  // build the temporary control points
411  for (i=0; i<=n; ++i)
412  {
413  eli::geom::utility::de_casteljau(tmp, B_v[i], v);
414  temp_cp.row(i)=tmp;
415  }
416  eli::geom::utility::de_casteljau(ans, temp_cp, u);
417  }
418 
419  return ans;
420  }
421 
422  point_type f_u(const data_type &u, const data_type &v) const
423  {
424  point_type ans, tmp;
425  index_type i, n(degree_u()), m(degree_v());
426 
427  // check to make sure have valid curve
428  assert(n>=0);
429  assert(m>=0);
430 
431  // check to make sure given valid parametric value
432  assert((u>=0) && (u<=1));
433  assert((v>=0) && (v<=1));
434 
435  if (this->degree_u()<1)
436  {
437  ans.setZero();
438  return ans;
439  }
440 
441  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_up(n+1-1, dim__);
442 
443  if ((n-1)<=m)
444  {
445  temp_cp.resize(m+1, dim__);
446  // build the temporary control points
447  for (i=0; i<=m; ++i)
448  {
449  // create the control points for first derivative curve
451  eli::geom::utility::de_casteljau(tmp, B_up, u);
452  temp_cp.row(i)=tmp;
453  }
454  eli::geom::utility::de_casteljau(ans, temp_cp, v);
455  }
456  else
457  {
458  temp_cp.resize(n+1, dim__);
459  // build the temporary control points
460  for (i=0; i<=n; ++i)
461  {
462  eli::geom::utility::de_casteljau(tmp, B_v[i], v);
463  temp_cp.row(i)=tmp;
464  }
465  // create the control points for first derivative curve
467  eli::geom::utility::de_casteljau(ans, B_up, u);
468  }
469 
470  return ans;
471  }
472 
473  point_type f_v(const data_type &u, const data_type &v) const
474  {
475  point_type ans, tmp;
476  index_type i, n(degree_u()), m(degree_v());
477 
478  // check to make sure have valid curve
479  assert(n>=0);
480  assert(m>=0);
481 
482  // check to make sure given valid parametric value
483  assert((u>=0) && (u<=1));
484  assert((v>=0) && (v<=1));
485 
486  if (this->degree_v()<1)
487  {
488  ans.setZero();
489  return ans;
490  }
491 
492  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_vp(m+1-1, dim__);
493 
494  if (n<=(m-1))
495  {
496  temp_cp.resize(m+1, dim__);
497  // build the temporary control points
498  for (i=0; i<=m; ++i)
499  {
500  eli::geom::utility::de_casteljau(tmp, B_u[i], u);
501  temp_cp.row(i)=tmp;
502  }
503  // create the control points for first derivative curve
505  eli::geom::utility::de_casteljau(ans, B_vp, v);
506  }
507  else
508  {
509  temp_cp.resize(n+1, dim__);
510  // build the temporary control points
511  for (i=0; i<=n; ++i)
512  {
513  // create the control points for first derivative curve
515  eli::geom::utility::de_casteljau(tmp, B_vp, v);
516  temp_cp.row(i)=tmp;
517  }
518  eli::geom::utility::de_casteljau(ans, temp_cp, u);
519  }
520 
521  return ans;
522  }
523 
524  point_type f_uu(const data_type &u, const data_type &v) const
525  {
526  point_type ans, tmp;
527  index_type i, n(degree_u()), m(degree_v());
528 
529  // check to make sure have valid curve
530  assert(n>=0);
531  assert(m>=0);
532 
533  // check to make sure given valid parametric value
534  assert((u>=0) && (u<=1));
535  assert((v>=0) && (v<=1));
536 
537  if (this->degree_u()<2)
538  {
539  ans.setZero();
540  return ans;
541  }
542 
543  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_upp(n+1-2, dim__);
544 
545  if ((n-2)<=m)
546  {
547  temp_cp.resize(m+1, dim__);
548  // build the temporary control points
549  for (i=0; i<=m; ++i)
550  {
551  // create the control points for second derivative curve
553  eli::geom::utility::de_casteljau(tmp, B_upp, u);
554  temp_cp.row(i)=tmp;
555  }
556  eli::geom::utility::de_casteljau(ans, temp_cp, v);
557  }
558  else
559  {
560  temp_cp.resize(n+1, dim__);
561  // build the temporary control points
562  for (i=0; i<=n; ++i)
563  {
564  eli::geom::utility::de_casteljau(tmp, B_v[i], v);
565  temp_cp.row(i)=tmp;
566  }
567 
568  // create the control points for second derivative curve
570  eli::geom::utility::de_casteljau(ans, B_upp, u);
571  }
572 
573  return ans;
574  }
575 
576  point_type f_uv(const data_type &u, const data_type &v) const
577  {
578  point_type ans, tmp;
579  index_type i, n(degree_u()), m(degree_v());
580 
581  // check to make sure have valid curve
582  assert(n>=0);
583  assert(m>=0);
584 
585  // check to make sure given valid parametric value
586  assert((u>=0) && (u<=1));
587  assert((v>=0) && (v<=1));
588 
589  if ( (this->degree_u()<1) || (this->degree_v()<1) )
590  {
591  ans.setZero();
592  return ans;
593  }
594 
595  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_up(n+1-1, dim__), B_vp(m+1-1, dim__);
596 
597  if (n<=m)
598  {
599  temp_cp.resize(m+1, dim__);
600  // build the temporary control points
601  for (i=0; i<=m; ++i)
602  {
603  // create the control points for first u-derivative curve
605  eli::geom::utility::de_casteljau(tmp, B_up, u);
606  temp_cp.row(i)=tmp;
607  }
608  // create the control points for first v-derivative curve
610  eli::geom::utility::de_casteljau(ans, B_vp, v);
611  }
612  else
613  {
614  temp_cp.resize(n+1, dim__);
615  // build the temporary control points
616  for (i=0; i<=n; ++i)
617  {
618  // create the control points for first v-derivative curve
620  eli::geom::utility::de_casteljau(tmp, B_vp, v);
621  temp_cp.row(i)=tmp;
622  }
623  // create the control points for first u-derivative curve
625  eli::geom::utility::de_casteljau(ans, B_up, u);
626  }
627 
628  return ans;
629  }
630 
631  point_type f_vv(const data_type &u, const data_type &v) const
632  {
633  point_type ans, tmp;
634  index_type i, n(degree_u()), m(degree_v());
635 
636  // check to make sure have valid curve
637  assert(n>=0);
638  assert(m>=0);
639 
640  // check to make sure given valid parametric value
641  assert((u>=0) && (u<=1));
642  assert((v>=0) && (v<=1));
643 
644  if (this->degree_v()<2)
645  {
646  ans.setZero();
647  return ans;
648  }
649 
650  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_vpp(m+1-2, dim__);
651 
652  if (n<=(m-2))
653  {
654  temp_cp.resize(m+1, dim__);
655  // build the temporary control points
656  for (i=0; i<=m; ++i)
657  {
658  eli::geom::utility::de_casteljau(tmp, B_u[i], u);
659  temp_cp.row(i)=tmp;
660  }
661  // create the control points for second derivative curve
663  eli::geom::utility::de_casteljau(ans, B_vpp, v);
664  }
665  else
666  {
667  temp_cp.resize(n+1, dim__);
668  // build the temporary control points
669  for (i=0; i<=n; ++i)
670  {
671  // create the control points for second derivative curve
673  eli::geom::utility::de_casteljau(tmp, B_vpp, v);
674  temp_cp.row(i)=tmp;
675  }
676  eli::geom::utility::de_casteljau(ans, temp_cp, u);
677  }
678 
679  return ans;
680  }
681 
682  point_type f_uuu(const data_type &u, const data_type &v) const
683  {
684 
685  point_type ans, tmp;
686  index_type i, n(degree_u()), m(degree_v());
687 
688  // check to make sure have valid curve
689  assert(n>=0);
690  assert(m>=0);
691 
692  // check to make sure given valid parametric value
693  assert((u>=0) && (u<=1));
694  assert((v>=0) && (v<=1));
695 
696  if (this->degree_u()<3)
697  {
698  ans.setZero();
699  return ans;
700  }
701 
702  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_uppp(n+1-3, dim__);
703 
704  if ((n-3)<=m)
705  {
706  temp_cp.resize(m+1, dim__);
707  // build the temporary control points
708  for (i=0; i<=m; ++i)
709  {
710  // create the control points for third derivative curve
712  eli::geom::utility::de_casteljau(tmp, B_uppp, u);
713  temp_cp.row(i)=tmp;
714  }
715  eli::geom::utility::de_casteljau(ans, temp_cp, v);
716  }
717  else
718  {
719  temp_cp.resize(n+1, dim__);
720  // build the temporary control points
721  for (i=0; i<=n; ++i)
722  {
723  eli::geom::utility::de_casteljau(tmp, B_v[i], v);
724  temp_cp.row(i)=tmp;
725  }
726 
727  // create the control points for third derivative curve
729  eli::geom::utility::de_casteljau(ans, B_uppp, u);
730  }
731 
732  return ans;
733  }
734 
735  point_type f_uuv(const data_type &u, const data_type &v) const
736  {
737  point_type ans, tmp;
738  index_type i, n(degree_u()), m(degree_v());
739 
740  // check to make sure have valid curve
741  assert(n>=0);
742  assert(m>=0);
743 
744  // check to make sure given valid parametric value
745  assert((u>=0) && (u<=1));
746  assert((v>=0) && (v<=1));
747 
748  if ( (this->degree_u()<2) || (this->degree_v()<1) )
749  {
750  ans.setZero();
751  return ans;
752  }
753 
754  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_upp(n+1-2, dim__), B_vp(m+1-1, dim__);
755 
756  if ((n-1)<=m)
757  {
758  temp_cp.resize(m+1, dim__);
759  // build the temporary control points
760  for (i=0; i<=m; ++i)
761  {
762  // create the control points for second u-derivative curve
764  eli::geom::utility::de_casteljau(tmp, B_upp, u);
765  temp_cp.row(i)=tmp;
766  }
767  // create the control points for first v-derivative curve
769  eli::geom::utility::de_casteljau(ans, B_vp, v);
770  }
771  else
772  {
773  temp_cp.resize(n+1, dim__);
774  // build the temporary control points
775  for (i=0; i<=n; ++i)
776  {
777  // create the control points for first v-derivative curve
779  eli::geom::utility::de_casteljau(tmp, B_vp, v);
780  temp_cp.row(i)=tmp;
781  }
782  // create the control points for second u-derivative curve
784  eli::geom::utility::de_casteljau(ans, B_upp, u);
785  }
786 
787  return ans;
788  }
789 
790  point_type f_uvv(const data_type &u, const data_type &v) const
791  {
792  point_type ans, tmp;
793  index_type i, n(degree_u()), m(degree_v());
794 
795  // check to make sure have valid curve
796  assert(n>=0);
797  assert(m>=0);
798 
799  // check to make sure given valid parametric value
800  assert((u>=0) && (u<=1));
801  assert((v>=0) && (v<=1));
802 
803  if ( (this->degree_u()<1) || (this->degree_v()<2) )
804  {
805  ans.setZero();
806  return ans;
807  }
808 
809  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_up(n+1-1, dim__), B_vpp(m+1-2, dim__);
810 
811  if (n<=(m-1))
812  {
813  temp_cp.resize(m+1, dim__);
814  // build the temporary control points
815  for (i=0; i<=m; ++i)
816  {
817  // create the control points for first u-derivative curve
819  eli::geom::utility::de_casteljau(tmp, B_up, u);
820  temp_cp.row(i)=tmp;
821  }
822  // create the control points for second v-derivative curve
824  eli::geom::utility::de_casteljau(ans, B_vpp, v);
825  }
826  else
827  {
828  temp_cp.resize(n+1, dim__);
829  // build the temporary control points
830  for (i=0; i<=n; ++i)
831  {
832  // create the control points for second v-derivative curve
834  eli::geom::utility::de_casteljau(tmp, B_vpp, v);
835  temp_cp.row(i)=tmp;
836  }
837  // create the control points for first u-derivative curve
839  eli::geom::utility::de_casteljau(ans, B_up, u);
840  }
841 
842  return ans;
843  }
844 
845  point_type f_vvv(const data_type &u, const data_type &v) const
846  {
847  point_type ans, tmp;
848  index_type i, n(degree_u()), m(degree_v());
849 
850  // check to make sure have valid curve
851  assert(n>=0);
852  assert(m>=0);
853 
854  // check to make sure given valid parametric value
855  assert((u>=0) && (u<=1));
856  assert((v>=0) && (v<=1));
857 
858  if (this->degree_v()<3)
859  {
860  ans.setZero();
861  return ans;
862  }
863 
864  Eigen::Matrix<data_type, Eigen::Dynamic, dim__> temp_cp, B_vppp(m+1-3, dim__);
865 
866  if (n<=(m-2))
867  {
868  temp_cp.resize(m+1, dim__);
869  // build the temporary control points
870  for (i=0; i<=m; ++i)
871  {
872  eli::geom::utility::de_casteljau(tmp, B_u[i], u);
873  temp_cp.row(i)=tmp;
874  }
875  // create the control points for third derivative curve
877  eli::geom::utility::de_casteljau(ans, B_vppp, v);
878  }
879  else
880  {
881  temp_cp.resize(n+1, dim__);
882  // build the temporary control points
883  for (i=0; i<=n; ++i)
884  {
885  // create the control points for third derivative curve
887  eli::geom::utility::de_casteljau(tmp, B_vppp, v);
888  temp_cp.row(i)=tmp;
889  }
890  eli::geom::utility::de_casteljau(ans, temp_cp, u);
891  }
892 
893  return ans;
894  }
895 
896  point_type normal(const data_type &u, const data_type &v) const
897  {
898  point_type n=f_u(u, v).cross(f_v(u, v));
899  data_type nlen(n.norm());
900  tolerance_type tol;
901 
902  // If have degenerate surface try higher order terms.
903  // Since want direction and don't care about magnitude,
904  // can use du=dv=1 in Taylor serier expansion:
905  // N(du, dv)=N+N_u*du+N_v*dv+1/2*(N_uu*du^2+2*N_uv*du*dv+N_vv*dv^2)+H.O.T.
906  // see "Bezier Normal Vector Surface and Its Applications" by Yamaguchi
907  if (tol.approximately_equal(nlen, 0))
908  {
909  point_type S_u, S_v, S_uu, S_uv, S_vv, N_u, N_v;
910  data_type du(1), dv(1);
911 
912  // calculate Taylor series second term
913  S_u=f_u(u, v);
914  S_v=f_v(u, v);
915  S_uu=f_uu(u, v);
916  S_uv=f_uv(u, v);
917  S_vv=f_vv(u, v);
918  N_u=S_uu.cross(S_v)+S_u.cross(S_uv);
919  N_v=S_uv.cross(S_v)+S_u.cross(S_vv);
920  n=N_u*du+N_v*dv;
921  nlen=n.norm();
922 
923  // if still have zero normal then calculate Taylor series third term
924  if (tol.approximately_equal(nlen, 0))
925  {
926  point_type S_uuu, S_uuv, S_uvv, S_vvv, N_uu, N_uv, N_vv;
927 
928  S_uuu=f_uuu(u, v);
929  S_uuv=f_uuv(u, v);
930  S_uvv=f_uvv(u, v);
931  S_vvv=f_vvv(u, v);
932  N_uu=S_uuu.cross(S_v)+2*S_uu.cross(S_uv)+S_u.cross(S_uuv);
933  N_uv=S_uuv.cross(S_v)+S_uu.cross(S_vv)+S_uv.cross(S_uv)+S_u.cross(S_uvv);
934  N_vv=S_uvv.cross(S_v)+2*S_uv.cross(S_vv)+S_u.cross(S_vvv);
935  n=0.5*(N_uu*du*du+2*N_uv*du*dv+N_vv*dv*dv);
936  nlen=n.norm();
937 
938  // if still get zero normal then give up
939  if (tol.approximately_equal(nlen, 0))
940  {
941  n.setZero();
942  nlen=1;
943  }
944  }
945  }
946 
947  n/=nlen;
948  return n;
949  }
950 
951  void promote_u()
952  {
953  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
954  typedef std::vector<control_row_type, Eigen::aligned_allocator<control_row_type> > control_row_collection_type;
955 
956  index_type i, n(degree_u()), m(degree_v());
957  control_row_collection_type current_row(m+1, control_row_type(n+1, dim__));
958 
959  // copy the control rows
960  for (i=0; i<=m; ++i)
961  current_row[i]=B_u[i];
962 
963  // resize current surface
964  resize(n+1, m);
965 
966  // set the new control points
967  control_row_type tmp_cp(n+2, dim__);
968  for (i=0; i<=m; ++i)
969  {
971  B_u[i]=tmp_cp;
972  }
973  }
974 
975  void promote_u_to(index_type target_degree)
976  {
977  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
978  typedef std::vector<control_row_type, Eigen::aligned_allocator<control_row_type> > control_row_collection_type;
979 
980  index_type i, n(degree_u()), m(degree_v());
981  control_row_collection_type current_row(m+1, control_row_type(n+1, dim__));
982 
983  // copy the control rows
984  for (i=0; i<=m; ++i)
985  {
986  current_row[i]=B_u[i];
987  }
988 
989  // resize current surface
990  resize(target_degree, m);
991 
992  // set the new control points
993  control_row_type tmp_cp(target_degree+1, dim__);
994  for (i=0; i<=m; ++i)
995  {
997  B_u[i]=tmp_cp;
998  }
999  }
1000 
1001  void promote_v()
1002  {
1003  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_col_type;
1004  typedef std::vector<control_col_type, Eigen::aligned_allocator<control_col_type> > control_col_collection_type;
1005 
1006  index_type i, n(degree_u()), m(degree_v());
1007  control_col_collection_type current_col(n+1, control_col_type(m+1, dim__));
1008 
1009  // copy the control cols
1010  for (i=0; i<=n; ++i)
1011  {
1012  current_col[i]=B_v[i];
1013  }
1014 
1015  // resize current surface
1016  resize(n, m+1);
1017 
1018  // set the new control points
1019  control_col_type tmp_cp(m+2, dim__);
1020  for (i=0; i<=n; ++i)
1021  {
1022  eli::geom::utility::bezier_promote_control_points(tmp_cp, current_col[i]);
1023  B_v[i]=tmp_cp;
1024  }
1025  }
1026 
1027  void promote_v_to(index_type target_degree)
1028  {
1029  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_col_type;
1030  typedef std::vector<control_col_type, Eigen::aligned_allocator<control_col_type> > control_col_collection_type;
1031 
1032  index_type i, n(degree_u()), m(degree_v());
1033  control_col_collection_type current_col(n+1, control_col_type(m+1, dim__));
1034 
1035  // copy the control cols
1036  for (i=0; i<=n; ++i)
1037  {
1038  current_col[i]=B_v[i];
1039  }
1040 
1041  // resize current surface
1042  resize(n, target_degree);
1043 
1044  // set the new control points
1045  control_col_type tmp_cp(target_degree+1, dim__);
1046  for (i=0; i<=n; ++i)
1047  {
1049  B_v[i]=tmp_cp;
1050  }
1051  }
1052 
1053  bool demote_u(const geom::general::continuity &u_continuity_degree=geom::general::C0)
1054  {
1055  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
1056  typedef std::vector<control_row_type, Eigen::aligned_allocator<control_row_type> > control_row_collection_type;
1057 
1058  index_type i, n(degree_u()), m(degree_v());
1059  control_row_collection_type current_row(m+1, control_row_type(n+1, dim__));
1060 
1061  // check if can demote
1062  int ncon(0);
1063  switch(u_continuity_degree)
1064  {
1066  ncon=0;
1067  break;
1068  case(eli::geom::general::C0):
1069  ncon=2;
1070  break;
1071  case(eli::geom::general::C1):
1072  ncon=4;
1073  break;
1074  case(eli::geom::general::C2):
1075  ncon=6;
1076  break;
1077  default:
1078  ncon=-1;
1079  }
1080  if (ncon<0)
1081  return false;
1082  if (ncon>n-2)
1083  return false;
1084 
1085  // copy the control rows
1086  for (i=0; i<=m; ++i)
1087  current_row[i]=B_u[i];
1088 
1089  // resize current surface
1090  resize(n-1, m);
1091 
1092  // set the new control points
1093  control_row_type tmp_cp(n, dim__);
1094  for (i=0; i<=m; ++i)
1095  {
1096  eli::geom::utility::bezier_demote_control_points(tmp_cp, current_row[i], ncon);
1097  B_u[i]=tmp_cp;
1098  }
1099 
1100  return true;
1101  }
1102 
1103  bool demote_v(const geom::general::continuity &v_continuity_degree=geom::general::C0)
1104  {
1105  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_col_type;
1106  typedef std::vector<control_col_type, Eigen::aligned_allocator<control_col_type> > control_col_collection_type;
1107 
1108  index_type i, n(degree_u()), m(degree_v());
1109  control_col_collection_type current_col(n+1, control_col_type(m+1, dim__));
1110 
1111  // check if can demote
1112  int ncon(0);
1113  switch(v_continuity_degree)
1114  {
1116  ncon=0;
1117  break;
1118  case(eli::geom::general::C0):
1119  ncon=2;
1120  break;
1121  case(eli::geom::general::C1):
1122  ncon=4;
1123  break;
1124  case(eli::geom::general::C2):
1125  ncon=6;
1126  break;
1127  default:
1128  ncon=-1;
1129  }
1130  if (ncon<0)
1131  return false;
1132  if (ncon>m-2)
1133  return false;
1134 
1135  // copy the control rows
1136  for (i=0; i<=n; ++i)
1137  current_col[i]=B_v[i];
1138 
1139  // resize current surface
1140  resize(n, m-1);
1141 
1142  // set the new control points
1143  control_col_type tmp_cp(m, dim__);
1144  for (i=0; i<=n; ++i)
1145  {
1146  eli::geom::utility::bezier_demote_control_points(tmp_cp, current_col[i], ncon);
1147  B_v[i]=tmp_cp;
1148  }
1149 
1150  return true;
1151  }
1152 
1153  void to_cubic_u()
1154  {
1155  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
1156  typedef std::vector<control_row_type, Eigen::aligned_allocator<control_row_type> > control_row_collection_type;
1157 
1158  index_type i, n(degree_u()), m(degree_v());
1159  control_row_collection_type current_row(m+1, control_row_type(n+1, dim__));
1160 
1161  // copy the control rows
1162  for (i=0; i<=m; ++i)
1163  {
1164  current_row[i]=B_u[i];
1165  }
1166 
1167  // resize current surface
1168  resize(3, m);
1169 
1170  // set the new control points
1171  control_row_type tmp_cp(4, dim__);
1172  for (i=0; i<=m; ++i)
1173  {
1175  B_u[i]=tmp_cp;
1176  }
1177  }
1178 
1179  void to_cubic_v()
1180  {
1181  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_col_type;
1182  typedef std::vector<control_col_type, Eigen::aligned_allocator<control_col_type> > control_col_collection_type;
1183 
1184  index_type i, n(degree_u()), m(degree_v());
1185  control_col_collection_type current_col(n+1, control_col_type(m+1, dim__));
1186 
1187  // copy the control rows
1188  for (i=0; i<=n; ++i)
1189  {
1190  current_col[i]=B_v[i];
1191  }
1192 
1193  // resize current surface
1194  resize(n, 3);
1195 
1196  // set the new control points
1197  control_col_type tmp_cp(4, dim__);
1198  for (i=0; i<=n; ++i)
1199  {
1201  B_v[i]=tmp_cp;
1202  }
1203  }
1204 
1205  void split_u(bezier<data_type, dim__, tol__> &bs_lo, bezier<data_type, dim__, tol__> &bs_hi, const data_type &u0) const
1206  {
1207  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_row_type;
1208 
1209  index_type i, j, n(degree_u()), m(degree_v());
1210  control_row_type cp_lo(n+1, dim__), cp_hi(n+1, dim__);
1211 
1212  // make sure have valid index
1213  assert((u0>=0) && (u0<=1));
1214 
1215  // resize the surfaces
1216  bs_lo.resize(n, m);
1217  bs_hi.resize(n, m);
1218 
1219  // cycle through each row and split each it
1220  for (j=0; j<=m; ++j)
1221  {
1222  eli::geom::utility::bezier_split_control_points(cp_lo, cp_hi, B_u[j], u0);
1223  for (i=0; i<=n; ++i)
1224  {
1225  bs_lo.set_control_point(cp_lo.row(i), i, j);
1226  bs_hi.set_control_point(cp_hi.row(i), i, j);
1227  }
1228  }
1229  }
1230 
1231  void split_v(bezier<data_type, dim__, tol__> &bs_lo, bezier<data_type, dim__, tol__> &bs_hi, const data_type &v0) const
1232  {
1233  typedef Eigen::Matrix<data_type, Eigen::Dynamic, dim__> control_col_type;
1234 
1235  index_type i, j, n(degree_u()), m(degree_v());
1236  control_col_type cp_lo(m+1, dim__), cp_hi(m+1, dim__);
1237 
1238  // make sure have valid index
1239  assert((v0>=0) && (v0<=1));
1240 
1241  // resize the surfaces
1242  bs_lo.resize(n, m);
1243  bs_hi.resize(n, m);
1244 
1245  // cycle through each col and split each it
1246  for (i=0; i<=n; ++i)
1247  {
1248  eli::geom::utility::bezier_split_control_points(cp_lo, cp_hi, B_v[i], v0);
1249  for (j=0; j<=m; ++j)
1250  {
1251  bs_lo.set_control_point(cp_lo.row(j), i, j);
1252  bs_hi.set_control_point(cp_hi.row(j), i, j);
1253  }
1254  }
1255  }
1256 
1258  {
1259  typedef bezier<data_type, dim__, tol__> surf_type;
1260 
1261  // Make working copies of surfaces.
1262  surf_type bsa(*this);
1263  surf_type bsb(bs);
1264 
1265  // Find maximum common order.
1266  index_type n(bsa.degree_u()), m(bsa.degree_v());
1267  if(bsb.degree_u() > n)
1268  {
1269  n = bsb.degree_u();
1270  }
1271 
1272  if(bsb.degree_v() > m)
1273  {
1274  m = bsb.degree_v();
1275  }
1276 
1277  // Promote both to max common order.
1278  bsa.promote_u_to(n);
1279  bsa.promote_v_to(m);
1280 
1281  bsb.promote_u_to(n);
1282  bsb.promote_v_to(m);
1283 
1284  // Find maximum distance between control points.
1285  index_type i, j;
1286  data_type d, maxd(0);
1287  for (i=0; i<=n; ++i)
1288  {
1289  for (j=0; j<=m; ++j)
1290  {
1291  d = (bsa.get_control_point(i, j) - bsb.get_control_point(i, j)).norm();
1292  if(d > maxd)
1293  {
1294  maxd=d;
1295  }
1296  }
1297  }
1298 
1299  return maxd;
1300  }
1301 
1302  private:
1303  void set_Bs(index_type n, index_type m)
1304  {
1305  // allocate vectors of control point maps
1306  B_u.resize(m+1, control_point_matrix_type(nullptr, m+1, dim__, Eigen::Stride<1, dim__>()));
1307  for (index_type j=0; j<=m; ++j)
1308  {
1309 #ifdef ELI_NO_VECTOR_DATA
1310  new (&(B_u.at(0))+j) control_point_matrix_type(&(point_data.at(0))+j*(n+1)*dim__, n+1, dim__, Eigen::Stride<1, dim__>());
1311 #else
1312  new (B_u.data()+j) control_point_matrix_type(point_data.data()+j*(n+1)*dim__, n+1, dim__, Eigen::Stride<1, dim__>());
1313 #endif
1314  }
1315 
1316  B_v.resize(n+1, v_dir_control_point_matrix_type(nullptr, n+1, dim__, Eigen::Stride<1, Eigen::Dynamic>(1, (n+1)*dim__)));
1317  for (index_type i=0; i<=n; ++i)
1318  {
1319 #ifdef ELI_NO_VECTOR_DATA
1320  new (&(B_v.at(0))+i) v_dir_control_point_matrix_type(&(point_data.at(0))+i*dim__, m+1, dim__, Eigen::Stride<1, Eigen::Dynamic>(1, (n+1)*dim__));
1321 #else
1322  new (B_v.data()+i) v_dir_control_point_matrix_type(point_data.data()+i*dim__, m+1, dim__, Eigen::Stride<1, Eigen::Dynamic>(1, (n+1)*dim__));
1323 #endif
1324  }
1325  }
1326  };
1327 
1333  }
1334  }
1335 }
1336 
1337 #endif
std::vector< data_type > control_point_container
Definition: bezier.hpp:58
void set_control_point(const point_type &cp, const index_type &i, const index_type &j)
Definition: bezier.hpp:265
point_type f_uuv(const data_type &u, const data_type &v) const
Definition: bezier.hpp:735
Eigen::Map< Eigen::Matrix< data_type, Eigen::Dynamic, dim__ >, Eigen::Unaligned, Eigen::Stride< 1, dim__ > > control_point_matrix_type
Definition: bezier.hpp:54
point_type f_uuu(const data_type &u, const data_type &v) const
Definition: bezier.hpp:682
bool operator==(const bezier< data_type, dim__, tol__ > &bs) const
Definition: bezier.hpp:96
Definition: continuity.hpp:26
point_type f_vvv(const data_type &u, const data_type &v) const
Definition: bezier.hpp:845
Eigen::Map< Eigen::Matrix< data_type, Eigen::Dynamic, dim__ >, Eigen::Unaligned, Eigen::Stride< 1, Eigen::Dynamic > > v_dir_control_point_matrix_type
Definition: bezier.hpp:57
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
bezier< double, 2 > bezier2d
Definition: bezier.hpp:1329
Definition: continuity.hpp:28
data_type get_vmax() const
Definition: bezier.hpp:164
Definition: bounding_box.hpp:27
void promote_v_to(index_type target_degree)
Definition: bezier.hpp:1027
point_type f_uu(const data_type &u, const data_type &v) const
Definition: bezier.hpp:524
bool closed_u() const
Definition: bezier.hpp:245
void bezier_control_points_to_cubic(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in)
Definition: bezier.hpp:142
bool operator!=(const bezier< data_type, dim__, tol__ > &bs) const
Definition: bezier.hpp:139
void bezier_ppp_control_point(Eigen::MatrixBase< Derived1 > &cp_ppp, const Eigen::MatrixBase< Derived2 > &cp)
Definition: bezier.hpp:78
std::vector< control_point_matrix_type > u_control_point_matrix_container
Definition: bezier.hpp:59
data_type eqp_distance_bound(const bezier< data_type, dim__, tol__ > &bs) const
Definition: bezier.hpp:1257
static dimension_type dimension()
Definition: bezier.hpp:144
point_type normal(const data_type &u, const data_type &v) const
Definition: bezier.hpp:896
Definition: bezier.hpp:38
bool demote_u(const geom::general::continuity &u_continuity_degree=geom::general::C0)
Definition: bezier.hpp:1053
void to_cubic_u()
Definition: bezier.hpp:1153
bezier< double, 3 > bezier3d
Definition: bezier.hpp:1330
void get_bounding_box(bounding_box_type &bb) const
Definition: bezier.hpp:190
~bezier()
Definition: bezier.hpp:83
bezier(const index_type &u_dim, const index_type &v_dim)
Definition: bezier.hpp:73
void resize(const index_type &u_dim, const index_type &v_dim)
Definition: bezier.hpp:166
control_point_container point_data
Definition: bezier.hpp:63
Eigen::Matrix< data_type, dim__, dim__ > rotation_matrix_type
Definition: bezier.hpp:47
point_type f_uvv(const data_type &u, const data_type &v) const
Definition: bezier.hpp:790
void scale(const data_type &s)
Definition: bezier.hpp:232
void de_casteljau(Eigen::MatrixBase< Derived1 > &p, const Eigen::MatrixBase< Derived2 > &cp, const typename Derived2::Scalar &t)
Definition: bezier.hpp:27
bezier< long double, 3 > bezier3ld
Definition: bezier.hpp:1332
index_type degree_u() const
Definition: bezier.hpp:146
bool closed_v() const
Definition: bezier.hpp:255
void get_uconst_curve(curve_type &bc, const data_type &u) const
Definition: bezier.hpp:340
bezier(const bezier< data_type, dim__, tol__ > &bs)
Definition: bezier.hpp:78
void reverse_v()
Definition: bezier.hpp:296
Definition: continuity.hpp:27
point_type f_u(const data_type &u, const data_type &v) const
Definition: bezier.hpp:422
bool open_v() const
Definition: bezier.hpp:254
bezier()
Definition: bezier.hpp:68
eli::geom::curve::bezier< data_type, dim__, tolerance_type > curve_type
Definition: bezier.hpp:49
data_type get_umax() const
Definition: bezier.hpp:163
void bezier_promote_control_points(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in)
Definition: bezier.hpp:93
point_type get_control_point(const index_type &i, const index_type &j) const
Definition: bezier.hpp:175
void rotate(const rotation_matrix_type &rmat)
Definition: bezier.hpp:204
bool add(const point_type &p)
Definition: bounding_box.hpp:113
data_type get_vmin() const
Definition: bezier.hpp:162
unsigned short dimension_type
Definition: bezier.hpp:41
Eigen::Matrix< data_type, 1, dim__ > point_type
Definition: bezier.hpp:43
bool abouteq(const bezier< data_type, dim__, tol__ > &bs, const data_type &ttol2) const
Definition: bezier.hpp:113
eli::geom::general::bounding_box< data_type, dim__, tolerance_type > bounding_box_type
Definition: bezier.hpp:48
void rotate(const rotation_matrix_type &rmat, const point_type &rorig)
Definition: bezier.hpp:213
bezier & operator=(const bezier< data_type, dim__, tol__ > &bs)
Definition: bezier.hpp:85
data__ data_type
Definition: bezier.hpp:42
void bezier_pp_control_point(Eigen::MatrixBase< Derived1 > &cp_pp, const Eigen::MatrixBase< Derived2 > &cp)
Definition: bezier.hpp:63
point_type f(const data_type &u, const data_type &v) const
Definition: bezier.hpp:382
bool open_u() const
Definition: bezier.hpp:244
void swap_uv()
Definition: bezier.hpp:315
tol__ tolerance_type
Definition: bezier.hpp:46
void reverse_u()
Definition: bezier.hpp:277
void split_v(bezier< data_type, dim__, tol__ > &bs_lo, bezier< data_type, dim__, tol__ > &bs_hi, const data_type &v0) const
Definition: bezier.hpp:1231
void bezier_demote_control_points(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in, int ncon)
Definition: bezier.hpp:179
control_point_type::Index index_type
Definition: bezier.hpp:45
void get_vconst_curve(curve_type &bc, const data_type &v) const
Definition: bezier.hpp:361
void resize(const index_type &t_dim)
Definition: bezier.hpp:186
bezier< long double, 2 > bezier2ld
Definition: bezier.hpp:1331
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
void to_cubic_v()
Definition: bezier.hpp:1179
void split_u(bezier< data_type, dim__, tol__ > &bs_lo, bezier< data_type, dim__, tol__ > &bs_hi, const data_type &u0) const
Definition: bezier.hpp:1205
void promote_u_to(index_type target_degree)
Definition: bezier.hpp:975
bool demote_v(const geom::general::continuity &v_continuity_degree=geom::general::C0)
Definition: bezier.hpp:1103
void set_Bs(index_type n, index_type m)
Definition: bezier.hpp:1303
continuity
Definition: continuity.hpp:24
bool equivalent_curves(const bezier< data__, dim__, tol__ > &c0, const bezier< data__, dim__, tol__ > &c1)
Definition: equivalent_curves.hpp:27
std::vector< v_dir_control_point_matrix_type > v_control_point_matrix_container
Definition: bezier.hpp:60
void promote_u()
Definition: bezier.hpp:951
void get_parameter_max(data_type &umax, data_type &vmax) const
Definition: bezier.hpp:155
Derived1__::Scalar distance2(const Eigen::MatrixBase< Derived1__ > &p1, const Eigen::MatrixBase< Derived2__ > &p2)
Definition: distance.hpp:27
void translate(const point_type &trans)
Definition: bezier.hpp:220
index_type degree_v() const
Definition: bezier.hpp:147
void bezier_promote_control_points_to(Eigen::MatrixBase< Derived1 > &cp_out, const Eigen::MatrixBase< Derived2 > &cp_in)
Definition: bezier.hpp:108
point_type f_vv(const data_type &u, const data_type &v) const
Definition: bezier.hpp:631
Definition: bezier.hpp:109
void get_parameter_min(data_type &umin, data_type &vmin) const
Definition: bezier.hpp:149
point_type f_v(const data_type &u, const data_type &v) const
Definition: bezier.hpp:473
v_control_point_matrix_container B_v
Definition: bezier.hpp:65
point_type control_point_type
Definition: bezier.hpp:44
Definition: continuity.hpp:29
u_control_point_matrix_container B_u
Definition: bezier.hpp:64
void promote_v()
Definition: bezier.hpp:1001
bezier< float, 3 > bezier3f
Definition: bezier.hpp:1328
data_type get_umin() const
Definition: bezier.hpp:161
point_type f_uv(const data_type &u, const data_type &v) const
Definition: bezier.hpp:576