DGtal  0.6.devel
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
SeparableMetricHelper.h
1 
17 #pragma once
18 
34 #if defined(SeparableMetricHelper_RECURSES)
35 #error Recursive header files inclusion detected in SeparableMetricHelper.h
36 #else // defined(SeparableMetricHelper_RECURSES)
37 
38 #define SeparableMetricHelper_RECURSES
39 
40 #if !defined SeparableMetricHelper_h
41 
42 #define SeparableMetricHelper_h
43 
45 // Inclusions
46 #include <iostream>
47 #include <cmath>
48 #include "DGtal/base/Common.h"
49 #include "DGtal/kernel/NumberTraits.h"
50 #include "DGtal/kernel/CBoundedInteger.h"
51 #include "DGtal/base/BasicFunctors.h"
53 
54 namespace DGtal
55 {
56 
58  // template class SeparableMetricHelper
74  template <typename TPoint, typename TInternalValue, DGtal::uint32_t tp>
75  struct SeparableMetricHelper
76  {
77  // ----------------------- Standard services ------------------------------
78 
79  typedef TInternalValue InternalValue;
80  typedef typename TPoint::Coordinate Abscissa;
81  typedef TPoint Point;
82 
83 
84  BOOST_CONCEPT_ASSERT(( CBoundedInteger<Abscissa> ));
85  BOOST_CONCEPT_ASSERT(( CBoundedInteger<TInternalValue> ));
86 
91  static const DGtal::uint32_t p = tp;
92 
93 
102  double getApproxValue( const InternalValue & aInternalValue ) const
103  {
104  return std::pow( NumberTraits<InternalValue>::castToDouble(aInternalValue),
105  (double) 1.0 / p);
106  }
107 
118  InternalValue F ( const Abscissa pos, const Abscissa ci, const InternalValue hi ) const
119  {
120  return std::pow( abs(NumberTraits<Abscissa>::castToDouble(pos - ci)),
121  (double)p) + hi;
122  }
123 
134  InternalValue reversedF ( const Abscissa pos, const Abscissa ci, const InternalValue hi ) const
135  {
136  return hi - std::pow( abs((double )pos - ci) , (double)p);
137  }
138 
139 
148  InternalValue power ( const Abscissa pos ) const
149  {
150  return ( InternalValue ) std::pow ( (double)pos, (double)p );
151  }
152 
153 
165  Abscissa reversedSep ( const Abscissa i, const InternalValue hi,
166  const Abscissa j, const InternalValue hj ) const
167  {
168  ASSERT(false && "Not-Yet-Implemented");
169  }
170 
182  Abscissa Sep ( const Abscissa i, const InternalValue hi,
183  const Abscissa j, const InternalValue hj ) const
184  {
185  ASSERT(false && "Not-Yet-Implemented");
186  }
187 
188 
189  enum Closest { FIRST=0, SECOND=1, BOTH=2};
190 
191 
203  Closest closest(const Point &origin,
204  const Point &first,
205  const Point &second) const
206  {
207  InternalValue a=NumberTraits<InternalValue>::ZERO,
208  b=NumberTraits<InternalValue>::ZERO;
209 
210  for(typename Point::Dimension i=0; i < Point::dimension; i++)
211  {
212  a += power(abs (origin[i] - first[i]));
213  b += power(abs (origin[i] - second[i]));
214  }
215  if (a<b)
216  return FIRST;
217  else
218  if (a>b)
219  return SECOND;
220  else
221  return BOTH;
222  }
223 
224 
241  Abscissa binarySearchHidden(const Abscissa &udim,
242  const Abscissa &vdim,
243  const InternalValue &nu,
244  const InternalValue &nv,
245  const Abscissa &lower,
246  const Abscissa &upper) const
247  {
248  ASSERT( (nu + (InternalValue) std::pow( (double)abs( udim - lower), (double) p)) <=
249  (nv + (InternalValue) std::pow( (double)abs( vdim - lower), (double)p)));
250 
251  //Recurrence stop
252  if ( (upper - lower) <= NumberTraits<Abscissa>::ONE)
253  return lower;
254 
255  Abscissa mid = (lower + upper)/2;
256  InternalValue nuUpdated = nu + (InternalValue) std::pow( (double)abs( udim - mid ), (double)p);
257  InternalValue nvUpdated = nv + (InternalValue) std::pow( (double)abs( vdim - mid ), (double)p);
258 
259  //Recursive call
260  if ( nuUpdated < nvUpdated)
261  return binarySearchHidden(udim,vdim,nu,nv,mid,upper);
262  else
263  return binarySearchHidden(udim,vdim,nu,nv,lower,mid);
264 
265  }
266 
285  bool hiddenBy(const Point &u,
286  const Point &v,
287  const Point &w,
288  const Point &startingPoint,
289  const Point &endPoint,
290  const typename Point::UnsignedComponent dim) const
291  {
292 
293  //Interval bound for the binary search
294  Abscissa lower = startingPoint[dim];
295  Abscissa upper = endPoint[dim];
296 
297  //Partial norm computation
298  // (sum_{i!=dim} |u_i-v_i|^p
299  InternalValue nu = NumberTraits<InternalValue>::ZERO;
300  InternalValue nv = NumberTraits<InternalValue>::ZERO;
301  InternalValue nw = NumberTraits<InternalValue>::ZERO;
302  for(Dimension i = 0 ; i < Point::dimension ; i++)
303  if (i != dim)
304  {
305  nu += ( InternalValue ) std::pow ( (double)abs(u[i] - startingPoint[i] ) , (double)p);
306  nv += ( InternalValue ) std::pow ( (double)abs(v[i] - startingPoint[i] ) , (double)p);
307  nw += ( InternalValue ) std::pow ( (double)abs(w[i] - startingPoint[i] ) , (double)p);
308  }
309 
310  //Intersection of voronoi boundary
311 
312 
313  //Optimization if vw lies before starting
314  if ((nv + (InternalValue) std::pow( (double)abs( v[dim] - lower), (double) p)) >
315  (nw + (InternalValue) std::pow( (double)abs( w[dim] - lower), (double)p)))
316  { //trace.endBlock();
317  return true;
318  }
319  //Optimization if vw lies before starting
320  if ((nu + (InternalValue) std::pow( (double)abs( u[dim] - lower), (double) p)) >
321  (nv + (InternalValue) std::pow( (double)abs( v[dim] - lower), (double)p)))
322  { //trace.endBlock();
323  return false;
324  }
325 
326  //Binary search
327  Abscissa uv = binarySearchHidden(u[dim],v[dim],nu,nv,lower,upper);
328  Abscissa vw = binarySearchHidden(v[dim],w[dim],nv,nw,lower,upper);
329 #ifdef DEBUG_VERBOSE
330  trace.info() << "Midpoint (u,v) ="<< uv<< " Midpoint (v,w) ="<<vw<<std::endl;
331 #endif
332 
333  //trace.endBlock();
334  if ( uv > vw )
335  return true;
336  else
337  return false;
338 
339  }
340 
341  }; // end of class SeparableMetricHelper
342 
343  // ------------------------------------------------------------------------
344  // ----------------------- Specializations ------------------------------
345  // ------------------------------------------------------------------------
346 
351  template <typename TPoint, typename TInternalValue>
352  struct SeparableMetricHelper<TPoint, TInternalValue, 2>
353  {
354  typedef TInternalValue InternalValue;
355  typedef typename TPoint::Coordinate Abscissa;
356  typedef TPoint Point;
357 
358  static const DGtal::uint32_t p = 2;
359 
360  inline double getApproxValue ( const InternalValue & aInternalValue ) const
361  {
362  return ( double ) sqrt ( NumberTraits<InternalValue>::castToDouble(aInternalValue) );
363  }
364 
365  inline InternalValue F ( const Abscissa pos,
366  const Abscissa ci,
367  const InternalValue hi ) const
368  {
369  return ( pos - ci ) * ( pos - ci ) + hi;
370  }
371 
372  inline InternalValue reversedF ( const Abscissa pos,
373  const Abscissa ci,
374  const InternalValue hi ) const
375  {
376  return hi - ( pos - ci ) * ( pos - ci ) ;
377  }
378 
379 
380  inline Abscissa Sep ( const Abscissa i, const InternalValue hi,
381  const Abscissa j, const InternalValue hj ) const
382  {
383  if ( ( ( j*j - i*i ) + hj - hi ) / ( 2* ( j - i ) ) >= 0)
384  return (Abscissa)( ( j*j - i*i ) + hj - hi ) / ( 2* ( j - i ) );
385  else
386  return (Abscissa)( ( j*j - i*i ) + hj - hi ) / ( 2* ( j - i ) ) -1;
387 
388  }
389 
390  inline Abscissa reversedSep ( const Abscissa i, const InternalValue hi,
391  const Abscissa j, const InternalValue hj ) const
392  {
393  return ( ( i*i -j*j ) + hj - hi ) / ( 2* ( i - j ) );
394  }
395 
396  inline InternalValue power ( const Abscissa i ) const
397  {
398  return (InternalValue) (i*i);
399  }
400  enum Closest { FIRST=0, SECOND=1, BOTH=2};
401 
402 
403  Closest closest(const Point &origin,
404  const Point &first,
405  const Point &second) const
406  {
407  InternalValue a=NumberTraits<InternalValue>::ZERO,
408  b=NumberTraits<InternalValue>::ZERO;
409 
410  for(typename Point::Dimension i=0; i < Point::dimension; i++)
411  {
412  a += (origin[i] - first[i])*(origin[i] - first[i]);
413  b += (origin[i] - second[i])*(origin[i] - second[i]);
414  }
415  if (a<b)
416  return FIRST;
417  else
418  if (a>b)
419  return SECOND;
420  else
421  return BOTH;
422  }
423 
438  bool hiddenBy(const Point &u,
439  const Point &v,
440  const Point &w,
441  const Point &startingPoint,
442  const Point &/*endPoint*/,
443  const typename Point::UnsignedComponent dim) const
444  {
445  //decide if (a,c) hide b in the lines (startingPoint, dim)
446 
447  Abscissa a,b, c;
448 
449  a = v[dim] - u[dim];
450  b = w[dim] - v[dim];
451  c = a + b;
452 
453  Abscissa d2_v=0, d2_u=0 ,d2_w=0;
454 
455  for(Dimension i = 0 ; i < Point::dimension ; i++)
456  if (i != dim)
457  {
458  d2_u += (u[i] - startingPoint[i] ) *(u[i] - startingPoint[i] );
459  d2_v += (v[i] - startingPoint[i] ) *(v[i] - startingPoint[i] );
460  d2_w += (w[i] - startingPoint[i] ) *(w[i] - startingPoint[i] );
461  }
462 
463  return (c * d2_v - b*d2_u - a*d2_w - a*b*c) > 0 ;
464  }
465 
466  };
467 
472  template <typename TPoint, typename TInternalValue>
473  struct SeparableMetricHelper<TPoint, TInternalValue, 1>
474  {
475 
476  typedef TInternalValue InternalValue;
477  static const DGtal::uint32_t p = 1;
478  typedef typename TPoint::Coordinate Abscissa;
479  typedef TPoint Point;
480 
481 
482  inline double getApproxValue ( const InternalValue & aInternalValue ) const
483  {
484  return ( double ) aInternalValue;
485  }
486 
487  inline InternalValue F ( const Abscissa pos,
488  const Abscissa ci,
489  const InternalValue hi ) const
490  {
491  return ( InternalValue ) ( ((long int) pos - ci)>=0 ? ((long int) pos - ci) : - ((long int) pos - ci) ) + hi;
492  }
493 
494  inline InternalValue reversedF ( const Abscissa pos,
495  const Abscissa ci,
496  const InternalValue hi ) const
497  {
498  return ( InternalValue ) hi - abs ( pos - ci );
499  }
500 
501 
502  inline Abscissa Sep ( const Abscissa i, const InternalValue hi,
503  const Abscissa j, const InternalValue hj ) const
504  {
505  if (hj >= hi + j - i)
506  return NumberTraits<Abscissa>::max();
507  if (hi > hj + j - i)
508  return NumberTraits<Abscissa>::min();
509  return (int)((hj - hi + j + i) / 2);
510  }
511 
512  inline Abscissa reversedSep ( const Abscissa i, const InternalValue hi,
513  const Abscissa j, const InternalValue hj ) const
514  {
515  if (hj <= hi - j + i)
516  return NumberTraits<Abscissa>::max();
517  if (hi < hj - j + i)
518  return NumberTraits<Abscissa>::min();
519  return (hi + i - hj + j ) / 2;
520  }
521 
522 
523  inline InternalValue power ( const Abscissa i ) const
524  {
525  return (InternalValue) abs(i);
526  }
527 
528  enum Closest { FIRST=0, SECOND=1, BOTH=2};
529 
530 
531  Closest closest(const Point &origin,
532  const Point &first,
533  const Point &second) const
534  {
535  InternalValue a=(origin-first).norm(Point::L_1),
536  b=(origin-second).norm(Point::L_1);
537 
538  if (a<b)
539  return FIRST;
540  else
541  if (a>b)
542  return SECOND;
543  else
544  return BOTH;
545  }
546 
561  bool hiddenBy(const Point &u,
562  const Point &v,
563  const Point &w,
564  const Point &startingPoint,
565  const Point &endPoint,
566  const typename Point::UnsignedComponent dim) const
567  {
568  ASSERT(false && "Not-Yet-Implemented");
569  Abscissa uv, vw;
570 
571 
572  }
573 
574  }; // end of class SeparableMetricHelper
575 
580  template <typename TPoint, typename TInternalValue>
581  struct SeparableMetricHelper<TPoint, TInternalValue, 0>
582  {
583  typedef TInternalValue InternalValue;
584  static const DGtal::uint32_t p = 0;
585  typedef typename TPoint::Coordinate Abscissa;
586  typedef TPoint Point;
587 
588 
589  inline double getApproxValue ( const InternalValue & aInternalValue ) const
590  {
591  return ( double ) aInternalValue;
592  }
593 
594  inline InternalValue F ( const Abscissa pos, const Abscissa ci,
595  const InternalValue hi ) const
596  {
597  return ( InternalValue )
598  max( (Abscissa) (((long int)pos - ci) >= 0 ? ((long int)pos - ci) :
599  -((long int)pos - ci)), (Abscissa) hi);
600  }
601 
602  inline InternalValue reversedF ( const Abscissa pos,
603  const Abscissa ci,
604  const InternalValue hi ) const
605  {
606  ASSERT(false && "Not-Implemented");
607  }
608 
609 
610  inline Abscissa Sep ( const Abscissa i, const InternalValue hi,
611  const Abscissa j, const InternalValue hj ) const
612  {
613  if (hi <= hj)
614  return max ((Abscissa)(i + hj), (Abscissa)(i + j) / 2);
615  else
616  return min ((Abscissa)(j - hi), (Abscissa)(i + j) / 2);
617  }
618 
619  inline Abscissa reversedSep ( const Abscissa i, const InternalValue hi,
620  const Abscissa j, const InternalValue hj ) const
621  {
622  ASSERT(false && "Not-Implemented");
623  }
624 
625 
626 
627  inline InternalValue power ( const Abscissa i ) const
628  {
629  return (InternalValue) abs(i);
630  }
631 
632  enum Closest { FIRST=0, SECOND=1, BOTH=2};
633 
634 
635  Closest closest(const Point &origin,
636  const Point &first,
637  const Point &second) const
638  {
639  InternalValue a=(origin-first).norm(Point::L_infty),
640  b=(origin-second).norm(Point::L_infty);
641 
642  if (a<b)
643  return FIRST;
644  else
645  if (a>b)
646  return SECOND;
647  else
648  return BOTH;
649  }
650 
665  bool hiddenBy(const Point &a,
666  const Point &b,
667  const Point &c,
668  const Point &startingPoint,
669  const Point &endPoint,
670  const typename Point::UnsignedComponent dim) const
671  {
672  ASSERT(false && "Not-Yet-Implemented");
673  }
674  }; // end of class SeparableMetricHelper
675 
676 
677 } // namespace DGtal
678 
679 
680 // //
682 
683 #endif // !defined SeparableMetricHelper_h
684 
685 #undef SeparableMetricHelper_RECURSES
686 #endif // else defined(SeparableMetricHelper_RECURSES)
Generated on Wed Dec 19 2012 19:10:34 for DGtal by   doxygen 1.8.1.1