DGtal 1.3.0
Loading...
Searching...
No Matches
FrechetShortcut.ih
1/**
2 * This program is free software: you can redistribute it and/or modify
3 * it under the terms of the GNU Lesser General Public License as
4 * published by the Free Software Foundation, either version 3 of the
5 * License, or (at your option) any later version.
6 *
7 * This program is distributed in the hope that it will be useful,
8 * but WITHOUT ANY WARRANTY; without even the implied warranty of
9 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
10 * GNU General Public License for more details.
11 *
12 * You should have received a copy of the GNU General Public License
13 * along with this program. If not, see <http://www.gnu.org/licenses/>.
14 *
15 **/
16
17/**
18 * @file FrechetShortcut.ih
19 * @author Isabelle Sivignon (\c isabelle.sivignon@gipsa-lab.grenoble-inp.fr )
20 * gipsa-lab Grenoble Images Parole Signal Automatique (CNRS, UMR 5216), CNRS, France
21 *
22 * @date 2012/02/24
23 *
24 * Implementation of inline methods defined in FrechetShortcut.h
25 *
26 * This file is part of the DGtal library.
27 */
28
29///////////////////////////////////////////////////////////////////////////////
30// IMPLEMENTATION of inline methods.
31///////////////////////////////////////////////////////////////////////////////
32
33//////////////////////////////////////////////////////////////////////////////
34#include <cstdlib>
35//////////////////////////////////////////////////////////////////////////////
36
37
38///////////////////////////////////////////////////////////////
39// Class backpath
40////////////////////////////////////////////////////////////////
41
42#define PRECISION 0.00001
43
44
45//creation of a backPath
46// Default constructor
47template <typename TIterator, typename TInteger>
48inline
49DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::Backpath()
50{
51 myQuad = 0;
52 myFlag = false;
53}
54
55
56//creation of a backPath
57template <typename TIterator, typename TInteger>
58inline
59DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::Backpath(const FrechetShortcut<TIterator,TInteger> *s,int q): myS(s),myQuad(q),myFlag(false)
60{
61}
62
63template <typename TIterator, typename TInteger>
64inline
65DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::Backpath(const Backpath & other): myS(other.myS),myQuad(other.myQuad),myFlag(other.myFlag),myOcculters(other.myOcculters),myForbiddenIntervals(other.myForbiddenIntervals)
66{
67}
68
69
70template <typename TIterator, typename TInteger>
71inline
72typename DGtal::FrechetShortcut<TIterator,TInteger>::Backpath& DGtal::FrechetShortcut<TIterator,TInteger>::Backpath::operator=(const Backpath & other)
73{
74 myS = other.myS;
75 myQuad = other.myQuad;
76 myFlag = other.myFlag;
77 myOcculters = other.myOcculters;
78 myForbiddenIntervals = IntervalSet(other.myForbiddenIntervals);
79 return *this;
80}
81
82
83template <typename TIterator, typename TInteger>
84inline
85void DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::reset()
86{
87 myFlag = false;
88 myOcculters.clear();
89 myForbiddenIntervals.clear();
90}
91
92
93//destruction of a backPath
94template <typename TIterator, typename TInteger>
95inline
96DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::~Backpath()
97{ }
98
99template <typename TIterator, typename TInteger>
100inline
101void DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::updateBackPathFirstQuad(int d, const ConstIterator& it)
102{
103
104 myIt = it;
105
106 switch(d)
107 {
108 case 0:
109 case 1:
110 case 2:
111 case 7:
112 {
113 addPositivePoint();
114 break;
115 }
116 case 3:
117 case 4:
118 case 5:
119 case 6:
120 {
121 addNegativePoint();
122 break;
123 }
124 }
125}
126
127// update the list of active occulters
128template <typename TIterator, typename TInteger>
129inline
130void DGtal::FrechetShortcut<TIterator,TInteger>::Backpath::updateOcculters()
131{
132
133 // The potential new occulter is the last-but-one point
134 Point p = Point(*(myIt-1));
135
136
137 Point pi,v;
138 Vector dir = Tools::chainCode2Vect(myQuad);
139 Vector dir_ortho = Tools::chainCode2Vect((myQuad+6)%8);
140
141
142 Point u1,u2;
143 u1 = dir;
144 u2 = Tools::chainCode2Vect((myQuad+1)%8);;
145
146 double angle_min=0;
147 double angle_max=M_PI_4;
148 bool occ = false;
149
150 IntegerComputer<TInteger> ic;
151
152 if(myOcculters.size()==0)
153 {
154 occ =true;
155 angle_min=0;
156 angle_max=M_PI_4;
157 }
158 else
159 {
160 typename occulter_list::iterator iter, next;
161 bool ok = true;
162 iter = myOcculters.begin();
163 while(iter!=myOcculters.end() && ok)
164 {
165 pi = Point(*(iter->first));
166 v = p-pi;
167
168 next = iter;
169 next++;
170 // pi is after p for all directions -> p is not an occulter
171 if(ic.dotProduct(v,u1) < 0 && ic.dotProduct(v,u2) <0)
172 {
173 ok = false;
174 occ = false;
175 }
176 else
177 // p is after pi for all directions -> pi is not an occulter
178 // anymore, p is a new occulter.
179 if(ic.dotProduct(v,u1) > 0 && ic.dotProduct(v,u2) > 0)
180 {
181 myOcculters.erase(iter);
182 occ = true;
183 angle_min = 0;
184 angle_max = M_PI_4;
185 }
186 else
187 // p is after pi on [0,alpha], before pi on [alpha,pi/4]
188 if(ic.dotProduct(v,u1) > 0 && ic.dotProduct(v,u2) <= 0)
189 {
190 double alpha = Tools::angleVectVect(v,dir_ortho);
191
192 if(alpha >= iter->second.angle_min && alpha <=
193 iter->second.angle_max)
194 {
195 // p is a new occulter
196 occ = true;
197 angle_min = 0;
198 angle_max = alpha;
199 // pi's angle_min is updated
200 iter->second.angle_min = alpha;
201 }
202 else
203 if(alpha > iter->second.angle_max)
204 {
205 //pi is not an occulter anymore
206 myOcculters.erase(iter);
207 occ=true;
208 angle_min = 0;
209 angle_max = M_PI_4;
210 }
211 // if alpha < iter->second.angle_min, pi does not
212 // change, p may be an occulter -> do nothing
213 }
214 else // scalar_product(v,u1) < 0 && scalar_product(v,u2) > 0
215 // p is after pi on [alpha,pi/4], before pi on [0,alpha]
216 {
217 double alpha = Tools::angleVectVect(v,dir_ortho);
218 alpha = M_PI - alpha;
219
220 if(alpha >= iter->second.angle_min && alpha <=
221 iter->second.angle_max)
222 {
223 occ = true;
224 angle_min = alpha;
225 angle_max = M_PI_4;
226 // pi's angle_max is updated
227 iter->second.angle_max = alpha;
228 }
229 else
230 if(alpha < iter->second.angle_min)
231 {
232 //pi is not an occulter anymore
233 myOcculters.erase(iter);
234 occ=true;
235 angle_min = 0;
236 angle_max = M_PI_4;
237 }
238 // if(alpha > iter->second.angle_max), pi does not
239 // change, p may be an occulter -> do nothing
240
241 }
242 iter = next;
243 }
244 }
245
246 if(occ)
247 {
248 occulter_attributes new_occ;
249 new_occ.angle_min = angle_min;
250 new_occ.angle_max = angle_max;
251 myOcculters.insert(myOcculters.end(),std::pair<const ConstIterator,occulter_attributes>(myIt-1,new_occ));
252
253 }
254
255
256}
257
258
259
260// update the set of intervals
261template <typename TIterator, typename TInteger>
262inline
263void DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::updateIntervals()
264{
265 Point p = Point(*myIt);
266
267 Point pi,v;
268 Vector dir,dir1;
269
270 IntegerComputer<TInteger> ic;
271
272 for(typename occulter_list::iterator iter =
273 myOcculters.begin(); iter!=myOcculters.end() ;++iter)
274 {
275 pi = Point(*(iter->first));
276
277 v = p-pi;
278
279 dir = Tools::chainCode2Vect(myQuad);
280 dir1 = Tools::chainCode2Vect((myQuad+1)%8);
281
282 if(ic.dotProduct(v,dir)<0 || ic.dotProduct(v,dir1)<0)
283 {
284 if(v.norm()>=myS->myError/sqrt(2.0F))
285 {
286 if(ic.crossProduct(dir,v)<=0)
287 {
288 v[0] = -v[0];
289 v[1] = -v[1];
290 }
291 double angle_v = Tools::angleVectVect(v,dir);
292
293 double tmp = acos((double) myS->myError/(sqrt(2.0F)*v.norm()));
294 double angle1 = -tmp+angle_v;
295 double angle2 = tmp+angle_v;
296 if(angle1 < 0)
297 angle1 = 0;
298 if(angle2 > M_PI_4)
299 angle2 = M_PI_4;
300
301 // Define a new interval of forbidden angles and insert it in the list.
302 boost::icl::interval<double>::type s = boost::icl::interval<double>::closed(angle1,angle2);
303 myForbiddenIntervals.insert(s);
304
305 }
306 }
307 }
308
309}
310
311
312// update the length of the longest backpath on a curve part
313template <typename TIterator, typename TInteger>
314inline
315void DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::addPositivePoint()
316{
317 // if we were on a monotone backpath, the point is an end of backpath
318 // otherwise, do nothing
319 if(myFlag)
320 {
321 myFlag=false;
322 }
323}
324
325
326
327/************************************************************/
328
329
330template <typename TIterator, typename TInteger>
331inline
332void DGtal::FrechetShortcut<TIterator, TInteger>::Backpath::addNegativePoint()
333{
334
335 // if we were on a monotone backpath, do nothing, the backpath
336 // continues
337 // otherwise it is the beggining of a new monotone backpath,
338 // possibly a locally maximal occulting point
339
340 //trace.info() << "add negative point" << std::endl;
341
342 if(!myFlag)
343 {
344 myFlag=true;
345 updateOcculters();
346 updateIntervals();
347 }
348 else
349 {
350 updateIntervals();
351 }
352
353
354}
355///////////////////////////////////////////////////////////////////////////
356// End of class backpath
357////////////////////////////////////////////////////////////////////////////
358
359//////////////////////////////////////////////////////////////////////////
360// Class cone
361////////////////////////////////////////////////////////////////////////////
362
363
364//creation of a cone
365
366template <typename TIterator, typename TInteger>
367inline
368DGtal::FrechetShortcut<TIterator,TInteger>::Cone::Cone()
369{
370 myInf = true;
371 myMin = 0;
372 myMax = 0;
373}
374
375
376
377template <typename TIterator, typename TInteger>
378inline
379DGtal::FrechetShortcut<TIterator,TInteger>::Cone::Cone(double angle0, double angle1)
380{
381
382 // angle0 and angle1 are ordered in direct orientation such that the
383 // angle made by the two directions is lower than PI.
384
385
386 // case angle0-angle1 = PI -> infinite cone
387 if(fabs(fabs(angle0-angle1)-M_PI) < PRECISION)
388 {
389 // the orientation is supposed to be ok (depends on the points involved)
390 myMin = angle0;
391 myMax = angle1;
392 }
393 else
394 if(fabs(angle0-angle1)<M_PI)
395 {
396 if(angle0-angle1>0)
397 {
398 myMin = angle1;
399 myMax = angle0;
400 }
401 else
402 {
403 myMin = angle0;
404 myMax = angle1;
405 }
406 }
407 else
408 {
409 // the cone includes the direction of angle=0
410 if(angle0>angle1)
411 {
412 myMin = angle0;
413 myMax = angle1;
414 }
415 else
416 {
417 myMin = angle1;
418 myMax = angle0;
419 }
420 }
421 myInf = false;
422}
423
424template <typename TIterator, typename TInteger>
425inline
426DGtal::FrechetShortcut<TIterator,TInteger>::Cone::Cone(double x, double y, double x0, double y0, double x1, double y1)
427{
428 double angle0 = Tools::computeAngle(x, y, x0, y0);
429 double angle1 = Tools::computeAngle(x, y, x1, y1);
430
431 assert(angle0 != -1 && angle1 != -1);
432
433 *this = Cone(angle0,angle1);
434 myInf = false;
435}
436
437template <typename TIterator, typename TInteger>
438inline
439bool DGtal::FrechetShortcut<TIterator,TInteger>::Cone::isEmpty() const
440{
441 if(myInf)
442 return false;
443 else
444 if(myMin==myMax)
445 return true;
446 else
447 return false;
448}
449
450template <typename TIterator, typename TInteger>
451inline
452typename DGtal::FrechetShortcut<TIterator,TInteger>::Cone& DGtal::FrechetShortcut<TIterator,TInteger>::Cone::operator=(const Cone& c)
453{
454 myMin =c.myMin;
455 myMax=c.myMax;
456 myInf = c.myInf;
457 return *this;
458}
459
460// // Computes the symmetrical cone
461template <typename TIterator, typename TInteger>
462inline
463typename DGtal::FrechetShortcut<TIterator,TInteger>::Cone DGtal::FrechetShortcut<TIterator,TInteger>::Cone::symmetricalCone()
464{
465 Cone cnew(myMin+M_PI,myMax+M_PI);
466 return cnew;
467}
468
469// Computes the intersection between the self Cone and another one.
470template <typename TIterator, typename TInteger>
471inline
472void DGtal::FrechetShortcut<TIterator,TInteger>::Cone::intersectCones(Cone c)
473{
474 Cone res;
475
476 // computes the intersection between the self cone and one half of another cone
477 res = intersectConesSimple(c);
478
479 // if they are disjoint, try the intersection with the other half
480 if(res.isEmpty())
481 {
482 Cone sym = c.symmetricalCone();
483 res = intersectConesSimple(sym);
484 }
485
486 *this = res;
487}
488
489
490//intersection of the self cone with another cone: considers only one half of the cone
491template <typename TIterator, typename TInteger>
492inline
493typename DGtal::FrechetShortcut<TIterator,TInteger>::Cone DGtal::FrechetShortcut<TIterator,TInteger>::Cone::intersectConesSimple(Cone c)
494{
495 Cone res;
496
497 // if the cone is infinite, the new cone is c
498 if(myInf)
499 {
500 res = c;
501 res.myInf = false;
502 }
503 else
504 // the directions of the new cone are not included in the old one
505 if(!Tools::isBetween(c.myMin, myMin, myMax, 2*M_PI) && !Tools::isBetween(c.myMax, myMin,
506 myMax,
507 2*M_PI))
508 {
509 // first possibility: the cones are disjoint
510 if(!Tools::isBetween(myMin, c.myMin, c.myMax, 2*M_PI) && !Tools::isBetween(myMax, c.myMin,
511 c.myMax, 2*M_PI))
512 res = Cone(0,0);
513 else
514 // or the new cone includes the old one, nothing changes, the cone remains the same.
515 res = *this;
516 }
517 else
518 // the old cone is "cut" by the new one
519 if(Tools::isBetween(c.myMin, myMin, myMax, 2*M_PI))
520 if(Tools::isBetween(c.myMax, myMin, myMax, 2*M_PI))
521 res = c;
522 else
523 res = Cone(c.myMin, myMax);
524 else
525 res = Cone(myMin,c.myMax);
526
527 return res;
528}
529
530template <typename TIterator, typename TInteger>
531inline
532void DGtal::FrechetShortcut<TIterator,TInteger>::Cone::selfDisplay ( std::ostream & out) const
533{
534 out << "[Cone]" << std::endl;
535 if(myInf)
536 out << "Infinite" << std::endl;
537 else
538 out << "[Cone min = " << myMin << " max = " << myMax << "]" << std::endl;
539 out << "[End Cone]" << std::endl;
540}
541
542
543/////////////////////////////////////////////////////////////////////////////
544/// FrechetShortcut class
545/////////////////////////////////////////////////////////////////////////////
546
547///////////////////////////////////////////////////////////////////////////////
548// Implementation of inline methods //
549
550template <typename TIterator, typename TInteger>
551inline
552DGtal::FrechetShortcut<TIterator,TInteger>::FrechetShortcut()
553{
554 myError = 0;
555 myCone = Cone();
556
557 for(int i=0;i<8;i++)
558 {
559 //backpath b(i,0);
560 Backpath b(this,i);
561 myBackpath.push_back(b);
562 }
563
564}
565
566
567template <typename TIterator, typename TInteger>
568inline
569DGtal::FrechetShortcut<TIterator,TInteger>::FrechetShortcut(double error)
570{
571 myError = error;
572 myCone = Cone();
573
574 for(int i=0;i<8;i++)
575 {
576 //backpath b(i,error);
577 Backpath b(this,i);
578 myBackpath.push_back(b);
579 }
580}
581
582
583
584
585template <typename TIterator, typename TInteger>
586inline
587void DGtal::FrechetShortcut<TIterator,TInteger>::init(const ConstIterator& it)
588{
589 myBegin = it;
590 myEnd = it;
591 resetCone();
592 resetBackpath();
593}
594
595
596template <typename TIterator, typename TInteger>
597inline
598DGtal::FrechetShortcut<TIterator,TInteger> DGtal::FrechetShortcut<TIterator,TInteger>::getSelf()
599{
600 FrechetShortcut<TIterator,TInteger> other = FrechetShortcut(myError);
601 return other;
602}
603
604
605template <typename TIterator, typename TInteger>
606inline
607DGtal::FrechetShortcut<TIterator,TInteger>::FrechetShortcut (const FrechetShortcut<TIterator,TInteger> & other ) : myError(other.myError), myBackpath(other.myBackpath), myCone(other.myCone), myBegin(other.myBegin), myEnd(other.myEnd){
608 resetBackpath();
609 resetCone();
610
611}
612
613template <typename TIterator, typename TInteger>
614inline
615DGtal::FrechetShortcut<TIterator,TInteger> & DGtal::FrechetShortcut<TIterator,TInteger>::operator=(const FrechetShortcut<TIterator,TInteger> & other)
616{
617
618 if(this != &other)
619 {
620 myError = other.myError;
621 myBackpath = other.myBackpath;
622 myCone = other.myCone;
623 myBegin = other.myBegin;
624 myEnd = other.myEnd;
625 }
626 return *this;
627}
628
629template <typename TIterator, typename TInteger>
630inline
631DGtal::FrechetShortcut<std::reverse_iterator<TIterator>,TInteger>
632DGtal::FrechetShortcut<TIterator,TInteger>
633::getReverse() const
634{
635 return Reverse(myError);
636}
637
638template <typename TIterator, typename TInteger>
639inline
640bool
641DGtal::FrechetShortcut<TIterator,TInteger>::operator==(
642 const DGtal::FrechetShortcut<TIterator,TInteger>& other) const {
643 return ((myBegin == other.myBegin) && (myEnd == other.myEnd) && (myError == other.myError));
644}
645
646
647
648template <typename TIterator, typename TInteger>
649inline
650bool
651DGtal::FrechetShortcut<TIterator,TInteger>::operator!=(
652 const DGtal::FrechetShortcut<TIterator,TInteger>& other) const {
653 return (!(*this == other));
654}
655
656
657template <typename TIterator, typename TInteger>
658inline
659bool
660DGtal::FrechetShortcut<TIterator,TInteger>::extendFront()
661{
662 bool flag = (updateWidth() && updateBackpath());
663
664 if(flag)
665 ++myEnd;
666
667 return flag;
668}
669
670
671template <typename TIterator, typename TInteger>
672inline
673bool
674DGtal::FrechetShortcut<TIterator,TInteger>::isExtendableFront()
675{
676
677 return (testUpdateWidth() && testUpdateBackpath());
678
679}
680
681template <typename TIterator, typename TInteger>
682inline
683typename DGtal::FrechetShortcut<TIterator,TInteger>::Cone
684DGtal::FrechetShortcut<TIterator,TInteger>::computeNewCone()
685{
686 double x0, y0,x1,y1;
687
688 Point firstP = Point(*myBegin);
689 Point newP = Point(*(myEnd+1));
690
691
692 Cone newCone=myCone;
693
694 // compute the tangent points defined by the first point and the
695 // circle C(newP,error)
696 bool intersect = Tools::circleTangentPoints(firstP[0],firstP[1], newP[0], newP[1], myError/(sqrt(2.0F)), &x0, &y0,
697 &x1, &y1);
698
699 if(intersect)
700 {
701 // define a cone according to the new tangent points
702 Cone c;
703 // case where there is one single tangent point
704 if(fabs(x0-x1) < PRECISION && fabs(y0-y1) < PRECISION)
705 {
706 double angle = Tools::computeAngle(firstP[0],firstP[1],newP[0],newP[1]);
707 assert(angle != -1);
708 double angle0 = angle - M_PI_2;
709 if(angle0<0)
710 angle0 = angle0+2*M_PI;
711 double angle1 = angle + M_PI_2;
712 if(angle1>2*M_PI)
713 angle1 = angle1-2*M_PI;
714 c = Cone(angle0,angle1);
715 }
716 else
717 c = Cone(firstP[0],firstP[1],x0,y0,x1,y1);
718
719 newCone.intersectCones(c);
720 }
721
722 return newCone;
723
724}
725
726// Test if the new direction belongs to the new cone, but does not
727// modify myCone
728template <typename TIterator, typename TInteger>
729inline
730bool DGtal::FrechetShortcut<TIterator,TInteger>::testUpdateWidth()
731{
732 Cone c = computeNewCone();
733
734 Point firstP = Point(*myBegin);
735 Point newP = Point(*(myEnd+1));
736
737 if(!(c.isEmpty()))
738 if(c.myInf)
739 return true;
740 else
741 {
742 double angle = Tools::computeAngle(firstP[0], firstP[1], newP[0], newP[1]);
743 assert(angle != -1);
744 return Tools::isBetween(angle,c.myMin,c.myMax,2*M_PI);
745 }
746 else
747 return false;
748
749}
750
751
752template <typename TIterator, typename TInteger>
753inline
754bool DGtal::FrechetShortcut<TIterator,TInteger>::testUpdateBackpath()
755{
756 // Save the current value of the backpath
757 std::vector <typename DGtal::FrechetShortcut<TIterator,TInteger>::Backpath> BackpathSave;
758
759 for(unsigned int i=0;i<8;i++)
760 {
761 Backpath b(myBackpath[i]);
762 BackpathSave.push_back(b);
763 }
764
765 // Check whether the next point could be added or not with respect to the backpath
766 bool flag = updateBackpath();
767
768 // Copy back the values of backpath before the test.
769 for(unsigned int i=0;i<8;i++)
770 myBackpath[i] = Backpath(BackpathSave[i]);
771
772 return flag;
773
774}
775
776
777// Same as testUpdateWidth() but myCone is modified.
778template <typename TIterator, typename TInteger>
779inline
780bool DGtal::FrechetShortcut<TIterator,TInteger>::updateWidth()
781{
782 Cone c = computeNewCone();
783
784 myCone = c;
785
786 Point firstP = Point(*myBegin);
787 Point newP = Point(*(myEnd+1));
788
789 bool flag = true;
790
791 if(!(c.isEmpty()))
792 if(c.myInf)
793 flag = true;
794 else
795 {
796 double angle = Tools::computeAngle(firstP[0], firstP[1], newP[0],
797 newP[1]);
798 assert(angle != -1);
799 flag = Tools::isBetween(angle,c.myMin,c.myMax,2*M_PI);
800 }
801 else
802 flag = false;
803
804 return flag;
805}
806
807template <typename TIterator, typename TInteger>
808inline
809bool DGtal::FrechetShortcut<TIterator,TInteger>::updateBackpath()
810{
811 Point prevP = Point(*myEnd);
812 Point P = Point(*(myEnd+1));
813
814 int d = Tools::computeChainCode(prevP,P);
815
816 for(unsigned int j=0;j<8;j++)
817 myBackpath[j].updateBackPathFirstQuad(Tools::rot(d,j),myEnd+1);
818
819
820 return isBackpathOk();
821
822}
823
824template <typename TIterator, typename TInteger>
825inline
826bool DGtal::FrechetShortcut<TIterator,TInteger>::isBackpathOk()
827{
828 // compute the quadrant of the direction of P(i,j)
829
830 Point firstP = Point(*myBegin);
831 Point P = Point(*(myEnd+1));
832
833 int q = Tools::computeOctant(firstP,P);
834
835
836 // to handle non simple curves (a point is visited twice)
837 if(firstP==P)
838 return true;
839
840 // compute the direction vector pipj
841 Point v;
842 v[0] = P[0]-firstP[0];
843 v[1] = P[1]-firstP[1];
844
845 // compute the angle between the direction vector and the elementary
846 // direction (defined by the quadrant)
847 Point dir_elem = Tools::chainCode2Vect(q);
848
849 double angle = Tools::angleVectVect(v,dir_elem);
850
851 boost::icl::interval_set<double> intervals = myBackpath[q].myForbiddenIntervals;
852
853 if(boost::icl::contains(intervals,angle))
854 return false;
855
856 return true;
857
858}
859
860
861template <typename TIterator, typename TInteger>
862inline
863void DGtal::FrechetShortcut<TIterator,TInteger>::resetBackpath()
864{
865 for(unsigned int i=0;i<8;i++)
866 {
867 myBackpath[i].reset();
868 }
869}
870
871template <typename TIterator, typename TInteger>
872inline
873void DGtal::FrechetShortcut<TIterator,TInteger>::resetCone()
874{
875 myCone.myMin = 0;
876 myCone.myMax = 0;
877 myCone.myInf = true;
878}
879
880
881template <typename TIterator, typename TInteger>
882inline
883TIterator
884DGtal::FrechetShortcut<TIterator,TInteger>::begin() const {
885 return myBegin;
886}
887
888template <typename TIterator, typename TInteger>
889inline
890TIterator
891DGtal::FrechetShortcut<TIterator,TInteger>::end() const {
892 ConstIterator i(myEnd); ++i;
893 return i;
894}
895
896
897
898template <typename TIterator, typename TInteger>
899inline
900std::string
901DGtal::FrechetShortcut<TIterator,TInteger>::className() const
902{
903 return "FrechetShortcut";
904}
905
906
907
908/**
909 * Writes/Displays the object on an output stream.
910 * @param out the output stream where the object is written.
911 */
912
913template <typename TIterator, typename TInteger>
914inline
915void
916DGtal::FrechetShortcut<TIterator,TInteger>::selfDisplay ( std::ostream & out) const
917{
918
919 out << "[FrechetShortcut]" << std::endl;
920 out << "(Begin, End)=";
921 out << "("<< Point(*myBegin) << ", " << Point(*myEnd) << ")\n";
922 out << "[End FrechetShortcut]" << std::endl;
923
924}
925
926// Implementation of inline functions //
927
928template <typename TIterator, typename TInteger>
929inline
930std::ostream&
931DGtal::operator<< ( std::ostream & out,
932 const DGtal::FrechetShortcut<TIterator,TInteger> & object )
933{
934 object.selfDisplay( out );
935 return out;
936}
937
938
939
940
941// //
942///////////////////////////////////////////////////////////////////////////////
943
944