DGtal 1.3.0
Loading...
Searching...
No Matches
FrechetShortcut.h
1
18#pragma once
19
32#if defined(FrechetShortcut_RECURSES)
33#error Recursive header files inclusion detected in FrechetShortcut.h
34#else // defined(FrechetShortcut_RECURSES)
36#define FrechetShortcut_RECURSES
37
38#if !defined FrechetShortcut_h
40#define FrechetShortcut_h
41
43// Inclusions
44#include "DGtal/base/Common.h"
45#include "DGtal/base/ConceptUtils.h"
46#include "DGtal/kernel/PointVector.h"
47#include "DGtal/arithmetic/IntegerComputer.h"
48#include <boost/icl/interval_set.hpp>
49#include <map>
50
52
53#include "DGtal/geometry/curves/SegmentComputerUtils.h"
54
55namespace DGtal
56{
57
59 // class FrechetShortcut
106 template <typename TIterator,typename TInteger = typename IteratorCirculatorTraits<TIterator>::Value::Coordinate>
108 {
109 // ----------------------- Standard services ------------------------------
110 public:
111
112 //entier
113
115 typedef TInteger Integer;
116
117
118
119 //required types
120 typedef TIterator ConstIterator;
123
124 //2D point and 2D vector
127 typedef typename Vector::Coordinate Coordinate;
128
134 {
135 private:
140
141 protected:
142
147 typedef struct occulter_attributes{
148 double angle_min; //
149 double angle_max; //
150 } occulter_attributes;
151
156 typedef std::map <ConstIterator,occulter_attributes > occulter_list;
157
158 public:
160
161
162 public:
163
164 typedef boost::icl::interval_set<double> IntervalSet;
165
169 int myQuad;
170
175 bool myFlag;
176
178
184
189
194
201
206 Backpath(const Backpath & other);
207
208
214 Backpath& operator=(const Backpath & other);
215
216
221
225 void reset();
226
231
236
244
249
254
255
256 }; // End of class Backpath
257
258
263 class Cone{
264
265 public:
269 double myMin;
270
274 double myMax;
275
279 bool myInf; // true if the cone is infinite (the whole plane)
280
282
288 Cone(double a0, double a1);
289
300 Cone(double x, double y, double x0, double y0, double x1, double y1);
301
306 Cone(const Cone& c) = default;
307
313 Cone& operator=(const Cone& c);
314
319 bool isEmpty() const;
320
326
333
339
344 void selfDisplay ( std::ostream & out) const;
345
346 };
347
348
349
350
351 struct Tools
352 {
362 static bool isBetween(double i, double a, double b, double n)
363 {
364 if(a<=b)
365 return (i>=a && i<=b);
366 else
367 return ((i>=a && i<=n) || (i>=0 && i<=b));
368 }
369
370
371 // Code by Tim Voght
372 // http://local.wasp.uwa.edu.au/~pbourke/geometry/2circle/tvoght.c
373
374 /* circle_circle_intersection() *
375 * Determine the points where 2 circles in a common plane intersect.
376 *
377 * int circle_circle_intersection(
378 * // center and radius of 1st circle
379 * double x0, double y0, double r0,
380 * // center and radius of 2nd circle
381 * double x1, double y1, double r1,
382 * // 1st intersection point
383 * double *xi, double *yi,
384 * // 2nd intersection point
385 * double *xi_prime, double *yi_prime)
386 *
387 * This is a public domain work. 3/26/2005 Tim Voght
388 *
389 */
397 static int circle_circle_intersection(double x0, double y0, double r0,
398 double x1, double y1, double r1,
399 double *xi, double *yi,
400 double *xi_prime, double *yi_prime)
401 {
402 double a, dx, dy, d, h, rx, ry;
403 double x2, y2;
404
405 /* dx and dy are the vertical and horizontal distances between
406 * the circle centers.
407 */
408 dx = x1 - x0;
409 dy = y1 - y0;
410
411 /* Determine the straight-line distance between the centers. */
412 //d = sqrt((dy*dy) + (dx*dx));
413 d = hypot(dx,dy); // Suggested by Keith Briggs
414
415 /* Check for solvability. */
416 if (d > (r0 + r1))
417 {
418 /* no solution. circles do not intersect. */
419 std::cerr << "Warning : the two circles do not intersect -> should never happen" << std::endl;
420 return 0; //I. Sivignon 05/2010 should never happen for our specific use.
421 }
422 if (d < fabs(r0 - r1))
423 {
424 /* no solution. one circle is contained in the other */
425 return 0;
426 }
427
428 /* 'point 2' is the point where the line through the circle
429 * intersection points crosses the line between the circle
430 * centers.
431 */
432
433 /* Determine the distance from point 0 to point 2. */
434 a = ((r0*r0) - (r1*r1) + (d*d)) / (2.0 * d) ;
435
436 /* Determine the coordinates of point 2. */
437 x2 = x0 + (dx * a/d);
438 y2 = y0 + (dy * a/d);
439
440 /* Determine the distance from point 2 to either of the
441 * intersection points.
442 */
443 h = sqrt((r0*r0) - (a*a));
444
445 /* Now determine the offsets of the intersection points from
446 * point 2.
447 */
448 rx = -dy * (h/d);
449 ry = dx * (h/d);
450
451 /* Determine the absolute intersection points. */
452 *xi = x2 + rx;
453 *xi_prime = x2 - rx;
454 *yi = y2 + ry;
455 *yi_prime = y2 - ry;
456
457 return 1;
458 }
459
460
461
481 static int circleTangentPoints(double x, double y, double x1, double y1, double r1, double *xi, double *yi,
482 double *xi_prime, double *yi_prime)
483 {
484 double x0 = (x+x1)/2;
485 double y0 = (y+y1)/2;
486 double r0 = sqrt((x-x1)*(x-x1) + (y-y1)*(y-y1))/2;
487
488 int res =
489 circle_circle_intersection(x0,y0,r0,x1,y1,r1,xi,yi,xi_prime,yi_prime);
490
491 return res;
492
493 }
494
495
505 static double computeAngle(double x0, double y0, double x1, double y1)
506 {
507 double x = x1-x0;
508 double y = y1-y0;
509
510 if(x!=0)
511 {
512 double alpha = y/x;
513
514 if(x>0 && y>=0)
515 return atan(alpha);
516 else
517 if(x>0 && y<0)
518 return atan(alpha)+2*M_PI;
519 else
520 if(x<0)
521 return atan(alpha)+M_PI;
522 }
523 else
524 {
525 if(y>0)
526 return M_PI_2;
527 else
528 return 3*M_PI_2;
529 }
530 return -1;
531 }
532
533
534
540 static double angleVectVect(Vector u, Vector v)
541 {
543 return acos((double)ic.dotProduct(u,v)/(u.norm()*v.norm()));
544 }
545
551 static int computeChainCode(Point p, Point q)
552 {
553 int d;
554 Coordinate x2 = q[0];
555 Coordinate y2 = q[1];
556 Coordinate x1 = p[0];
557 Coordinate y1 = p[1];
558
559 if(x2-x1==0)
560 if(y2-y1==1)
561 d=2;
562 else
563 d=6;
564 else
565 if(x2-x1==1)
566 if(y2-y1==0)
567 d=0;
568 else
569 if(y2-y1==1)
570 d=1;
571 else
572 d=7;
573 else
574 if(y2-y1==0)
575 d=4;
576 else
577 if(y2-y1==1)
578 d=3;
579 else
580 d=5;
581 return d;
582 }
583
589 static int computeOctant(Point p, Point q)
590 {
591 int d = 0;
592 Coordinate x = q[0]-p[0];
593 Coordinate y = q[1]-p[1];
594
595 if(x>=0)
596 {
597 if(y>=0)
598 {
599 if(x>y)
600 d=0; // 0 <= y < x
601 else
602 if(x!=0)
603 d=1; // 0 <= x <= y
604 }
605 else
606 {
607 if(x>=abs(y))
608 d=7; // 0 < abs(y) <= x
609 else
610 d=6; // 0 <= x < abs(y)
611 }
612 }
613 if(x<=0)
614 {
615 if(y>0)
616 if(abs(x)>=y)
617 d=3; //
618 else
619 d=2;
620 else
621 {
622 if(abs(x)>abs(y))
623 d=4;
624 else
625 if(x!=0)
626 d=5;
627 }
628 }
629 return d;
630
631 }
632
641 static int rot(int d, int quad)
642 {
643 return (d-quad+8)%8;
644 }
645
651 static Vector chainCode2Vect(int d)
652 {
653 Vector v;
654 switch(d){
655
656 case 0:{
657 v[0] = 1;
658 v[1] = 0;
659 break;
660 }
661 case 1:{
662 v[0] = 1;
663 v[1] = 1;
664 break;
665 }
666 case 2:{
667 v[0] = 0;
668 v[1] = 1;
669 break;
670 }
671 case 3:{
672 v[0] = -1;
673 v[1] = 1;
674 break;
675 }
676 case 4:{
677 v[0] = -1;
678 v[1] = 0;
679 break;
680 }
681 case 5:{
682 v[0] = -1;
683 v[1] = -1;
684 break;
685 }
686 case 6:{
687 v[0] = 0;
688 v[1] = -1;
689 break;
690 }
691 case 7:{
692 v[0] = 1;
693 v[1] = -1;
694 break;
695 }
696
697 }
698
699 return v;
700 }
701
702
703 };
704
705
711
716 FrechetShortcut(double error);
717
723 void init(const ConstIterator& it);
724
726
731 FrechetShortcut ( const Self & other );
732
738 FrechetShortcut & operator= ( const Self & other );
739
744
751 bool operator==( const Self & other ) const;
752
759 bool operator!=( const Self & other ) const;
760
761
762
767
768 // ----------------------- Interface --------------------------------------
769public:
770
777
784
785
786 // ---------------------------- Accessors ----------------------------------
787
788
798
799
800
801
802 public:
803
807 std::string className() const;
808
809
814 bool isValid() const;
815
816 // ------------------------- Protected Datas ------------------------------
817 protected:
818
822 double myError;
823
828 std::vector <Backpath> myBackpath;
829
834
843
844
845
846 // ------------------------- Private Datas --------------------------------
847 private:
848
849
850
851 // ------------------------- Hidden services ------------------------------
852
853 public:
854
862
869
875
880
884 void resetCone();
885
891
898
904
905
906
907
912 void selfDisplay ( std::ostream & out ) const ;
913
914
915
916 // ------------------------- Internals ------------------------------------
917 private:
918
919 }; // end of class FrechetShortcut
920
921
922 // Utils
923
924
931 template <typename TIterator,typename TInteger>
932 std::ostream&
933 operator<< ( std::ostream & out, const FrechetShortcut<TIterator,TInteger> & object );
934
935
936
937
938
939
940} // namespace DGtal
941
942
944// Includes inline functions.
945#include "DGtal/geometry/curves/FrechetShortcut.ih"
946// //
948
949#endif // !defined FrechetShortcut_h
950
951#undef FrechetShortcut_RECURSES
952#endif // else defined(FrechetShortcut_RECURSES)
const FrechetShortcut< ConstIterator, Integer > * myS
boost::icl::interval_set< double > IntervalSet
void updateBackPathFirstQuad(int d, const ConstIterator &)
std::map< ConstIterator, occulter_attributes > occulter_list
Backpath & operator=(const Backpath &other)
Backpath(const Backpath &other)
Backpath(const FrechetShortcut< ConstIterator, Integer > *s, int q)
Cone & operator=(const Cone &c)
void selfDisplay(std::ostream &out) const
Cone(const Cone &c)=default
Cone(double a0, double a1)
Cone(double x, double y, double x0, double y0, double x1, double y1)
Aim: On-line computation Computation of the longest shortcut according to the Fr├ęchet distance for a ...
std::string className() const
FrechetShortcut & operator=(const Self &other)
FrechetShortcut< std::reverse_iterator< ConstIterator >, Integer > Reverse
ConstIterator begin() const
IteratorCirculatorTraits< ConstIterator >::Value Point
bool operator!=(const Self &other) const
Vector::Coordinate Coordinate
Reverse getReverse() const
IteratorCirculatorTraits< ConstIterator >::Value Vector
FrechetShortcut(const Self &other)
void selfDisplay(std::ostream &out) const
void init(const ConstIterator &it)
FrechetShortcut(double error)
FrechetShortcut getSelf()
BOOST_CONCEPT_ASSERT((concepts::CInteger< TInteger >))
bool operator==(const Self &other) const
ConstIterator end() const
std::vector< Backpath > myBackpath
FrechetShortcut< ConstIterator, Integer > Self
Aim: This class gathers several types and methods to make computation with integers.
Integer dotProduct(const Vector2I &u, const Vector2I &v) const
DGtal is the top-level namespace which contains all DGtal functions and types.
std::ostream & operator<<(std::ostream &out, const ClosedIntegerHalfPlane< TSpace > &object)
static int circle_circle_intersection(double x0, double y0, double r0, double x1, double y1, double r1, double *xi, double *yi, double *xi_prime, double *yi_prime)
static double angleVectVect(Vector u, Vector v)
static Vector chainCode2Vect(int d)
static int circleTangentPoints(double x, double y, double x1, double y1, double r1, double *xi, double *yi, double *xi_prime, double *yi_prime)
static bool isBetween(double i, double a, double b, double n)
static double computeAngle(double x0, double y0, double x1, double y1)
static int rot(int d, int quad)
static int computeChainCode(Point p, Point q)
static int computeOctant(Point p, Point q)
Aim: Concept checking for Integer Numbers. More precisely, this concept is a refinement of both CEucl...
Definition: CInteger.h:88