DGtal  0.6.devel
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
COBANaivePlane.ih
1 
30 
31 #include <cstdlib>
33 
35 // IMPLEMENTATION of inline methods.
37 
39 // ----------------------- Standard services ------------------------------
40 
41 //-----------------------------------------------------------------------------
42 template <typename TSpace, typename TInternalInteger>
43 inline
44 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
45 ~COBANaivePlane()
46 { // Nothing to do.
47 }
48 //-----------------------------------------------------------------------------
49 template <typename TSpace, typename TInternalInteger>
50 inline
51 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
52 COBANaivePlane()
53  : myG( NumberTraits<TInternalInteger>::ZERO )
54 { // Object is invalid
55 }
56 //-----------------------------------------------------------------------------
57 template <typename TSpace, typename TInternalInteger>
58 inline
59 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
60 COBANaivePlane( const COBANaivePlane & other )
61  : myAxis( other.myAxis ),
62  myG( other.myG ),
63  myWidth( other.myWidth ),
64  myPointSet( other.myPointSet ),
65  myState( other.myState ),
66  myCst1( other.myCst1 ),
67  myCst2( other.myCst2 )
68 {
69 }
70 //-----------------------------------------------------------------------------
71 template <typename TSpace, typename TInternalInteger>
72 inline
73 DGtal::COBANaivePlane<TSpace, TInternalInteger> &
74 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
75 operator=( const COBANaivePlane & other )
76 {
77  if ( this != &other )
78  {
79  myAxis = other.myAxis;
80  myG = other.myG;
81  myWidth = other.myWidth;
82  myPointSet = other.myPointSet;
83  myState = other.myState;
84  myCst1 = other.myCst1;
85  myCst2 = other.myCst2;
86  }
87  return *this;
88 }
89 //-----------------------------------------------------------------------------
90 template <typename TSpace, typename TInternalInteger>
91 inline
92 typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::MyIntegerComputer &
93 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
94 ic() const
95 {
96  return myState.cip.ic();
97 }
98 //-----------------------------------------------------------------------------
99 template <typename TSpace, typename TInternalInteger>
100 inline
101 void
102 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
103 clear()
104 {
105  myPointSet.clear();
106  myState.cip.clear();
107  // initialize the search space as a square.
108  myState.cip.pushBack( InternalPoint2( -myG, -myG ) );
109  myState.cip.pushBack( InternalPoint2( myG, -myG ) );
110  myState.cip.pushBack( InternalPoint2( myG, myG ) );
111  myState.cip.pushBack( InternalPoint2( -myG, myG ) );
112  computeCentroidAndNormal( myState );
113 }
114 //-----------------------------------------------------------------------------
115 template <typename TSpace, typename TInternalInteger>
116 inline
117 void
118 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
119 init( Dimension axis, InternalInteger diameter,
120  InternalInteger widthNumerator,
121  InternalInteger widthDenominator )
122 {
123  myAxis = axis;
124  myWidth[ 0 ] = widthNumerator;
125  myWidth[ 1 ] = widthDenominator;
126  // initialize the grid step.
127  myG = 2*diameter; myG *= diameter; myG *= diameter;
128  // Initializes some constants
129  // _cst1 = ( (int) ceil( get_si( myG ) * myWidth ) + 1 );
130  // _cst2 = ( (int) floor( get_si( myG ) * myWidth ) - 1 );
131  myCst1 = ( ( myG * myWidth[ 0 ] - 1 ) / myWidth[ 1 ] ) + 2;
132  myCst2 = ( ( myG * myWidth[ 0 ] ) / myWidth[ 1 ] ) - 1;
133  clear();
134 }
135 //-----------------------------------------------------------------------------
136 template <typename TSpace, typename TInternalInteger>
137 inline
138 typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::Size
139 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
140 size() const
141 {
142  return myPointSet.size();
143 }
144 //-----------------------------------------------------------------------------
145 template <typename TSpace, typename TInternalInteger>
146 inline
147 bool
148 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
149 empty() const
150 {
151  return myPointSet.empty();
152 }
153 //-----------------------------------------------------------------------------
154 template <typename TSpace, typename TInternalInteger>
155 inline
156 typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::Size
157 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
158 max_size() const
159 {
160  return myPointSet.max_size();
161 }
162 //-----------------------------------------------------------------------------
163 template <typename TSpace, typename TInternalInteger>
164 inline
165 typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::Size
166 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
167 maxSize() const
168 {
169  return max_size();
170 }
171 //-----------------------------------------------------------------------------
172 template <typename TSpace, typename TInternalInteger>
173 inline
174 typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::Size
175 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
176 complexity() const
177 {
178  return myState.cip.size();
179 }
180 //-----------------------------------------------------------------------------
181 template <typename TSpace, typename TInternalInteger>
182 inline
183 bool
184 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
185 operator()( const Point & p ) const
186 {
187  ic().getDotProduct( _v, myState.N, p );
188  return ( _v >= myState.min ) && ( _v <= myState.max );
189 }
190 //-----------------------------------------------------------------------------
191 template <typename TSpace, typename TInternalInteger>
192 inline
193 bool
194 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
195 extendAsIs( const Point & p )
196 {
197  ASSERT( isValid() && ! empty() );
198  bool ok = this->operator()( p );
199  if ( ok ) myPointSet.insert( p );
200  return ok;
201 }
202 
203 //-----------------------------------------------------------------------------
204 template <typename TSpace, typename TInternalInteger>
205 bool
206 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
207 extend( const Point & p )
208 {
209  ASSERT( isValid() );
210  // Checks if first point.
211  if ( empty() )
212  {
213  myPointSet.insert( p );
214  ic().getDotProduct( myState.max, myState.N, p );
215  myState.min = myState.max;
216  myState.ptMax = myState.ptMin = p;
217  return true;
218  }
219 
220  // Check first if p is already a point of the plane.
221  if ( myPointSet.find( p ) != myPointSet.end() ) // already in set
222  return true;
223  // Check if p lies within the current bounds of the plane.
224  _state.N = myState.N;
225  _state.min = myState.min;
226  _state.max = myState.max;
227  _state.ptMin = myState.ptMin;
228  _state.ptMax = myState.ptMax;
229  bool changed = updateMinMax( _state, &p, (&p)+1 );
230  // Check if point is already within bounds.
231  if ( ! changed )
232  {
233  myPointSet.insert( p );
234  return true;
235  }
236  // Check if width is still ok
237  if ( checkPlaneWidth( _state ) )
238  {
239  myState.min = _state.min;
240  myState.max = _state.max;
241  myState.ptMin = _state.ptMin;
242  myState.ptMax = _state.ptMax;
243  myPointSet.insert( p );
244  return true;
245  }
246  // We have to find a new normal. First, update gradient.
247  computeGradient( _grad, _state );
248 
249  // Checks if we can change the normal so as to find another digital plane.
250  if( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
251  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
252  {
253  // Unable to update solution.
254  return false;
255  }
256 
257  // There is a gradient, tries to optimize.
258  _state.cip = myState.cip;
259  doubleCut( _grad, _state );
260 
261  // While at least 1 point left on the search space
262  while ( ! _state.cip.empty() )
263  {
264  computeCentroidAndNormal( _state );
265  // Calls oracle
266  computeMinMax( _state, myPointSet.begin(), myPointSet.end() );
267  updateMinMax( _state, &p, (&p)+1 );
268  // Check if width is now ok
269  if ( checkPlaneWidth( _state ) )
270  { // Found a plane.
271  myState.min = _state.min;
272  myState.max = _state.max;
273  myState.ptMin = _state.ptMin;
274  myState.ptMax = _state.ptMax;
275  myState.cip.swap( _state.cip );
276  myState.centroid = _state.centroid;
277  myState.N = _state.N;
278  myPointSet.insert( p );
279  return true;
280  }
281 
282  // We have to find a new normal. First, update gradient.
283  computeGradient( _grad, _state );
284 
285  // Checks if we can change the normal so as to find another digital plane.
286  if ( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
287  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
288  {
289  // Unable to update solution. Removes point from set.
290  break;
291  }
292 
293  // There is a gradient, tries to optimize.
294  doubleCut( _grad, _state );
295  }
296  // was unable to find a correct plane.
297  return false;
298 }
299 //-----------------------------------------------------------------------------
300 template <typename TSpace, typename TInternalInteger>
301 bool
302 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
303 isExtendable( const Point & p ) const
304 {
305  ASSERT( isValid() );
306  // Checks if first point.
307  if ( empty() ) return true;
308 
309  // Check first if p is already a point of the plane.
310  if ( myPointSet.find( p ) != myPointSet.end() ) // already in set
311  return true;
312  // Check if p lies within the current bounds of the plane.
313  _state.N = myState.N;
314  _state.min = myState.min;
315  _state.max = myState.max;
316  _state.ptMin = myState.ptMin;
317  _state.ptMax = myState.ptMax;
318  bool changed = updateMinMax( _state, (&p), (&p)+1 );
319  // Check if point is already within bounds.
320  if ( ! changed ) return true;
321  // Check if width is still ok
322  if ( checkPlaneWidth( _state ) )
323  return true;
324  // We have to find a new normal. First, update gradient.
325  computeGradient( _grad, _state );
326 
327  // Checks if we can change the normal so as to find another digital plane.
328  if( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
329  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
330  // Unable to update solution.
331  return false;
332 
333  // There is a gradient, tries to optimize.
334  _state.cip = myState.cip;
335  doubleCut( _grad, _state );
336 
337  // While at least 1 point left on the search space
338  while ( ! _state.cip.empty() )
339  {
340  computeCentroidAndNormal( _state );
341  // Calls oracle
342  computeMinMax( _state, myPointSet.begin(), myPointSet.end() );
343  updateMinMax( _state, (&p), (&p)+1 );
344  // Check if width is now ok
345  if ( checkPlaneWidth( _state ) )
346  // Found a plane.
347  return true;
348 
349  // We have to find a new normal. First, update gradient.
350  computeGradient( _grad, _state );
351 
352  // Checks if we can change the normal so as to find another digital plane.
353  if ( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
354  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
355  // Unable to update solution.
356  return false;
357 
358  // There is a gradient, tries to optimize.
359  doubleCut( _grad, _state );
360  }
361  // was unable to find a correct plane.
362  return false;
363 }
364 //-----------------------------------------------------------------------------
365 template <typename TSpace, typename TInternalInteger>
366 template <typename TInputIterator>
367 bool
368 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
369 extend( TInputIterator it, TInputIterator itE )
370 {
371  BOOST_CONCEPT_ASSERT(( boost::InputIterator<TInputIterator> ));
372 
373  ASSERT( isValid() );
374 
375  // Check if points lies within the current bounds of the plane.
376  bool changed;
377  _state.N = myState.N;
378  if ( empty() )
379  {
380  changed = true;
381  computeMinMax( _state, it, itE );
382  }
383  else
384  {
385  _state.min = myState.min;
386  _state.max = myState.max;
387  _state.ptMin = myState.ptMin;
388  _state.ptMax = myState.ptMax;
389  changed = updateMinMax( _state, it, itE );
390  }
391  // Check if points are already within bounds.
392  if ( ! changed )
393  { // All points are within bounds. Put them in pointset.
394  for ( TInputIterator tmpIt = it; tmpIt != itE; ++tmpIt )
395  myPointSet.insert( *tmpIt );
396  return true;
397  }
398  // Check if width is still ok
399  if ( checkPlaneWidth( _state ) )
400  {
401  myState.min = _state.min;
402  myState.max = _state.max;
403  myState.ptMin = _state.ptMin;
404  myState.ptMax = _state.ptMax;
405  for ( TInputIterator tmpIt = it; tmpIt != itE; ++tmpIt )
406  myPointSet.insert( *tmpIt );
407  return true;
408  }
409  // We have to find a new normal. First, update gradient.
410  computeGradient( _grad, _state );
411 
412  // Checks if we can change the normal so as to find another digital plane.
413  if( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
414  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
415  {
416  // Unable to update solution.
417  return false;
418  }
419 
420  // There is a gradient, tries to optimize.
421  _state.cip = myState.cip;
422  doubleCut( _grad, _state );
423 
424  // While at least 1 point left on the search space
425  while ( ! _state.cip.empty() )
426  {
427  computeCentroidAndNormal( _state );
428  // Calls oracle
429  computeMinMax( _state, myPointSet.begin(), myPointSet.end() );
430  updateMinMax( _state, it, itE );
431  // Check if width is now ok
432  if ( checkPlaneWidth( _state ) )
433  { // Found a plane.
434  myState.min = _state.min;
435  myState.max = _state.max;
436  myState.ptMin = _state.ptMin;
437  myState.ptMax = _state.ptMax;
438  myState.cip.swap( _state.cip );
439  myState.centroid = _state.centroid;
440  myState.N = _state.N;
441  for ( TInputIterator tmpIt = it; tmpIt != itE; ++tmpIt )
442  myPointSet.insert( *tmpIt );
443  return true;
444  }
445 
446  // We have to find a new normal. First, update gradient.
447  computeGradient( _grad, _state );
448 
449  // Checks if we can change the normal so as to find another digital plane.
450  if ( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
451  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
452  {
453  // Unable to update solution. Removes point from set.
454  break;
455  }
456 
457  // There is a gradient, tries to optimize.
458  doubleCut( _grad, _state );
459  }
460  // was unable to find a correct plane.
461  return false;
462 }
463 //-----------------------------------------------------------------------------
464 template <typename TSpace, typename TInternalInteger>
465 template <typename TInputIterator>
466 bool
467 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
468 isExtendable( TInputIterator it, TInputIterator itE ) const
469 {
470  BOOST_CONCEPT_ASSERT(( boost::InputIterator<TInputIterator> ));
471 
472  ASSERT( isValid() );
473  // Check if points lies within the current bounds of the plane.
474  bool changed;
475  _state.N = myState.N;
476  if ( empty() )
477  {
478  changed = true;
479  computeMinMax( _state, it, itE );
480  }
481  else
482  {
483  _state.N = myState.N;
484  _state.min = myState.min;
485  _state.max = myState.max;
486  _state.ptMin = myState.ptMin;
487  _state.ptMax = myState.ptMax;
488  changed = updateMinMax( _state, it, itE );
489  }
490 
491  // Check if point is already within bounds.
492  if ( ! changed ) return true;
493  // Check if width is still ok
494  if ( checkPlaneWidth( _state ) )
495  return true;
496  // We have to find a new normal. First, update gradient.
497  computeGradient( _grad, _state );
498 
499  // Checks if we can change the normal so as to find another digital plane.
500  if( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
501  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
502  // Unable to update solution.
503  return false;
504 
505  // There is a gradient, tries to optimize.
506  _state.cip = myState.cip;
507  doubleCut( _grad, _state );
508 
509  // While at least 1 point left on the search space
510  while ( ! _state.cip.empty() )
511  {
512  computeCentroidAndNormal( _state );
513  // Calls oracle
514  computeMinMax( _state, myPointSet.begin(), myPointSet.end() );
515  updateMinMax( _state, it, itE );
516  // Check if width is now ok
517  if ( checkPlaneWidth( _state ) )
518  // Found a plane.
519  return true;
520 
521  // We have to find a new normal. First, update gradient.
522  computeGradient( _grad, _state );
523 
524  // Checks if we can change the normal so as to find another digital plane.
525  if ( ( ( _grad[ 0 ] == NumberTraits<InternalInteger>::ZERO )
526  && ( _grad[ 1 ] == NumberTraits<InternalInteger>::ZERO ) ) )
527  // Unable to update solution.
528  return false;
529 
530  // There is a gradient, tries to optimize.
531  doubleCut( _grad, _state );
532  }
533  // was unable to find a correct plane.
534  return false;
535 }
536 
537 
538 //-----------------------------------------------------------------------------
539 template <typename TSpace, typename TInternalInteger>
540 template <typename Vector3D>
541 inline
542 void
543 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
544 getNormal( Vector3D & normal ) const
545 {
546  switch( myAxis ) {
547  case 0 :
548  normal[0] = 1.0;
549  normal[1] = NumberTraits<InternalInteger>::castToDouble( myState.N[ 1 ]) / NumberTraits<InternalInteger>::castToDouble( myState.N[ 0 ] );
550  normal[2] = NumberTraits<InternalInteger>::castToDouble( myState.N[ 2 ]) / NumberTraits<InternalInteger>::castToDouble( myState.N[ 0 ] );
551  break;
552  case 1 :
553  normal[0] = NumberTraits<InternalInteger>::castToDouble( myState.N[ 0 ]) / NumberTraits<InternalInteger>::castToDouble( myState.N[ 1 ] );
554  normal[1] = 1.0;
555  normal[2] = NumberTraits<InternalInteger>::castToDouble( myState.N[ 2 ]) / NumberTraits<InternalInteger>::castToDouble( myState.N[ 1 ] );
556  break;
557  case 2 :
558  normal[0] = NumberTraits<InternalInteger>::castToDouble( myState.N[ 0 ]) / NumberTraits<InternalInteger>::castToDouble( myState.N[ 2 ] );
559  normal[1] = NumberTraits<InternalInteger>::castToDouble( myState.N[ 1 ]) / NumberTraits<InternalInteger>::castToDouble( myState.N[ 2 ] );
560  normal[2] = 1.0;
561  break;
562  }
563 }
564 //-----------------------------------------------------------------------------
565 //-----------------------------------------------------------------------------
566 template <typename TSpace, typename TInternalInteger>
567 template <typename Vector3D>
568 inline
569 void
570 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
571 getUnitNormal( Vector3D & normal ) const
572 {
573  getNormal( normal );
574  double l = sqrt( normal[ 0 ] * normal[ 0 ]
575  + normal[ 1 ] * normal[ 1 ]
576  + normal[ 2 ] * normal[ 2 ] );
577  normal[ 0 ] /= l;
578  normal[ 1 ] /= l;
579  normal[ 2 ] /= l;
580 }
581 //-----------------------------------------------------------------------------
582 template <typename TSpace, typename TInternalInteger>
583 inline
584 void
585 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
586 getBounds( double & min, double & max ) const
587 {
588  double nx = NumberTraits<InternalInteger>::castToDouble( myState.N[ 0 ] );
589  double ny = NumberTraits<InternalInteger>::castToDouble( myState.N[ 1 ] );
590  double nz = NumberTraits<InternalInteger>::castToDouble( myState.N[ 2 ] );
591  double l = sqrt( nx*nx + ny*ny + nz*nz );
592  min = NumberTraits<InternalInteger>::castToDouble( myState.min ) / l;
593  max = NumberTraits<InternalInteger>::castToDouble( myState.max ) / l;
594 }
595 //-----------------------------------------------------------------------------
596 template <typename TSpace, typename TInternalInteger>
597 inline
598 const typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::Point &
599 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
600 minimalPoint() const
601 {
602  ASSERT( ! this->empty() );
603  return *(myState.indMin);
604 }
605 //-----------------------------------------------------------------------------
606 template <typename TSpace, typename TInternalInteger>
607 inline
608 const typename DGtal::COBANaivePlane<TSpace, TInternalInteger>::Point &
609 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
610 maximalPoint() const
611 {
612  ASSERT( ! this->empty() );
613  return *(myState.indMax);
614 }
615 
616 
617 
619 // Interface - public :
620 
625 template <typename TSpace, typename TInternalInteger>
626 inline
627 void
628 DGtal::COBANaivePlane<TSpace, TInternalInteger>::selfDisplay ( std::ostream & out ) const
629 {
630  double min, max;
631  double N[ 3 ];
632  out << "[COBANaivePlane"
633  << " axis=" << myAxis << " w=" << myWidth[ 0 ] << "/" << myWidth[ 1 ]
634  << " size=" << size() << " complexity=" << complexity() << " N=" << myState.N << ": ";
635  this->getUnitNormal( N );
636  this->getBounds( min, max );
637  out << min << " <= "
638  << N[ 0 ] << " * x + "
639  << N[ 1 ] << " * y + "
640  << N[ 2 ] << " * z "
641  << " <= " << max << " ]";
642 }
643 
648 template <typename TSpace, typename TInternalInteger>
649 inline
650 bool
651 DGtal::COBANaivePlane<TSpace, TInternalInteger>::isValid() const
652 {
653  return myG != NumberTraits< InternalInteger >::ZERO;
654 }
655 
656 
658 // Internals
660 //-----------------------------------------------------------------------------
661 template <typename TSpace, typename TInternalInteger>
662 inline
663 void
664 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
665 computeCentroidAndNormal( State & state ) const
666 {
667  if ( state.cip.empty() ) return;
668  state.centroid = state.cip.centroid();
669  ic().reduce( state.centroid );
670  switch( myAxis ){
671  case 0 :
672  state.N[ 0 ] = state.centroid[ 2 ] * myG;
673  state.N[ 1 ] = state.centroid[ 0 ];
674  state.N[ 2 ] = state.centroid[ 1 ];
675  break;
676  case 1 :
677  state.N[ 0 ] = state.centroid[ 0 ];
678  state.N[ 1 ] = state.centroid[ 2 ] * myG;
679  state.N[ 2 ] = state.centroid[ 1 ];
680  break;
681  case 2 :
682  state.N[ 0 ] = state.centroid[ 0 ];
683  state.N[ 1 ] = state.centroid[ 1 ];
684  state.N[ 2 ] = state.centroid[ 2 ] * myG;
685  break;
686  }
687 
688 }
689 //-----------------------------------------------------------------------------
690 template <typename TSpace, typename TInternalInteger>
691 inline
692 void
693 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
694 doubleCut( InternalPoint2 & grad, State & state ) const
695 {
696  // 2 cuts on the search space:
697  // Gradient.p <= cst1 - _v
698  // -Gradient.p <= cst2 + _v
699  _v = myG * ( state.ptMin[ myAxis ]
700  - state.ptMax[ myAxis ] );
701  state.cip.cut( HalfSpace( grad, myCst1 - _v ) );
702  grad.negate();
703  state.cip.cut( HalfSpace( grad, myCst2 + _v ) );
704  grad.negate();
705 }
706 
707 //-----------------------------------------------------------------------------
708 template <typename TSpace, typename TInternalInteger>
709 template <typename TInputIterator>
710 void
711 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
712 computeMinMax( State & state, TInputIterator itB, TInputIterator itE ) const
713 {
714  BOOST_CONCEPT_ASSERT(( boost::InputIterator<TInputIterator> ));
715 
716  ASSERT( itB != itE );
717  ic().getDotProduct( state.min, state.N, *itB );
718  state.max = state.min;
719  state.ptMax = state.ptMin = *itB;
720  ++itB;
721  // look for the points defining the min dot product and the max dot product
722  for ( ; itB != itE; ++itB )
723  {
724  ic().getDotProduct( _v, state.N, *itB );
725  if ( _v > state.max )
726  {
727  state.max = _v;
728  state.ptMax = *itB;
729  }
730  else if ( _v < state.min )
731  {
732  state.min = _v;
733  state.ptMin = *itB;
734  }
735  }
736 }
737 //-----------------------------------------------------------------------------
738 template <typename TSpace, typename TInternalInteger>
739 template <typename TInputIterator>
740 bool
741 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
742 updateMinMax( State & state, TInputIterator itB, TInputIterator itE ) const
743 
744 {
745  BOOST_CONCEPT_ASSERT(( boost::InputIterator<TInputIterator> ));
746 
747  bool changed = false;
748  // look for the points defining the min dot product and the max dot product
749  for ( ; itB != itE; ++itB )
750  {
751  ic().getDotProduct( _v, state.N, *itB );
752  if ( _v > state.max )
753  {
754  state.max = _v;
755  state.ptMax = *itB;
756  changed = true;
757  }
758  else if ( _v < state.min )
759  {
760  state.min = _v;
761  state.ptMin = *itB;
762  changed = true;
763  }
764  }
765  return changed;
766 }
767 //-----------------------------------------------------------------------------
768 template <typename TSpace, typename TInternalInteger>
769 bool
770 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
771 checkPlaneWidth( const State & state ) const
772 {
773  _v = ic().abs( state.N[ myAxis ] );
774  return ( ( state.max - state.min )
775  < ( _v * myWidth[ 0 ] / myWidth[ 1 ] ) );
776 }
777 //-----------------------------------------------------------------------------
778 template <typename TSpace, typename TInternalInteger>
779 void
780 DGtal::COBANaivePlane<TSpace, TInternalInteger>::
781 computeGradient( InternalPoint2 & grad, const State & state ) const
782 {
783  // computation of the gradient
784  switch( myAxis ){
785  case 0 :
786  grad[ 0 ] = state.ptMin[ 1 ] - state.ptMax[ 1 ];
787  grad[ 1 ] = state.ptMin[ 2 ] - state.ptMax[ 2 ];
788  break;
789  case 1:
790  grad[ 0 ] = state.ptMin[ 0 ] - state.ptMax[ 0 ];
791  grad[ 1 ] = state.ptMin[ 2 ] - state.ptMax[ 2 ];
792  break;
793  case 2:
794  grad[ 0 ] = state.ptMin[ 0 ] - state.ptMax[ 0 ];
795  grad[ 1 ] = state.ptMin[ 1 ] - state.ptMax[ 1 ];
796  break;
797  }
798 }
799 
800 
802 // Implementation of inline functions //
803 
804 template <typename TSpace, typename TInternalInteger>
805 inline
806 std::ostream&
807 DGtal::operator<< ( std::ostream & out,
808  const COBANaivePlane<TSpace, TInternalInteger> & object )
809 {
810  object.selfDisplay( out );
811  return out;
812 }
813 
814 // //
816 
817 
Generated on Wed Dec 19 2012 19:10:21 for DGtal by   doxygen 1.8.1.1