00001 #ifndef __Fiber_MultiIndex_HPP
00002 #define __Fiber_MultiIndex_HPP
00003
00004 #include <eagle/vecmath.hpp>
00005 #include <eagle/Assignment.hpp>
00006 #include <eagle/Operator.hpp>
00007 #include <eagle/FixedArray.hpp>
00008 #include "Iterator.hpp"
00009
00010 namespace Fiber
00011 {
00012
00037 struct CreateFromLinear {};
00038
00039 struct Power2Alignment {};
00040
00052 template <int Dims>
00053 class MultiIndex : public MultiIndex<Dims-1>
00054 {
00055 typedef MultiIndex<Dims-1> base_t;
00056
00057 index_t idx;
00058
00059 public:
00060 typedef index_t value_type;
00061
00063 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Add&) throw()
00064 : base_t(M,D, Add() )
00065 , idx(M.idx + D.maxidx())
00066 {}
00067
00069 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Sub&) throw()
00070 : base_t(M,D, Sub() )
00071 , idx(M.idx - D.maxidx())
00072 {}
00073
00075 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Mult&) throw()
00076 : base_t(M,D, Mult() )
00077 , idx(M.idx * D.maxidx())
00078 {}
00079
00081 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Div&)
00082 : base_t(M,D, Div() )
00083 , idx(M.idx / D.maxidx())
00084 {}
00085
00086 MultiIndex(unsigned int bits, const ::Eagle::BinaryAnd&) throw()
00087 : base_t(bits, ::Eagle::BinaryAnd() )
00088 , idx( bits & (1<<(Dims-1)) ? 1 : 0 )
00089 {}
00090
00091 MultiIndex(const MultiIndex&M, const MultiIndex&D, const ::Eagle::BinaryAnd&) throw()
00092 : base_t(M,D, ::Eagle::BinaryAnd() )
00093 , idx( M.idx > D.maxidx() ? M.idx - D.maxidx() : D.maxidx() - M.idx )
00094 {}
00095
00096 MultiIndex(const MultiIndex&M, const Power2Alignment&) throw()
00097 : base_t(M, Power2Alignment() )
00098 , idx( RoundUpToNextHighestPowerOf2( M.idx) )
00099 {}
00100
00101 public:
00106 friend MultiIndex distance(const MultiIndex&M0, const MultiIndex&M1)
00107 {
00108 return MultiIndex(M0,M1, ::Eagle::BinaryAnd() );
00109 }
00110
00111 friend MultiIndex bitalign(const MultiIndex&M) throw()
00112 {
00113 return MultiIndex(M, Power2Alignment() );
00114 }
00115
00135 static MultiIndex BitIndex(unsigned int bits) throw()
00136 {
00137 return MultiIndex(bits, ::Eagle::BinaryAnd() );
00138 }
00139
00150 static MultiIndex Axis(unsigned int orientation) throw()
00151 {
00152 return BitIndex( 1 << orientation );
00153 }
00154
00170 index_t BitIndex() const throw()
00171 {
00172 return linear( MultiIndex(2) );
00173 }
00174
00179 static int log2(index_t N) throw()
00180 {
00181 if (N & (1<<(Dims-1) ) )
00182 return Dims-1;
00183
00184 return base_t::log2(N);
00185 }
00186
00206 index_t Orientation() const throw()
00207 {
00208 return log2( BitIndex() );
00209 }
00210
00211 enum
00212 {
00214 dims = Dims,
00215
00217 SIZE = Dims
00218 };
00219
00229 int getSignedOrientation(const MultiIndex&B) const
00230 {
00231 if (idx<B.idx) return +Dims;
00232 if (idx>B.idx) return -Dims;
00233 return subidx().getSignedOrientation(B.subidx() );
00234 }
00235
00237 MultiIndex() throw()
00238 : idx(0)
00239 {}
00240
00242 explicit MultiIndex(const index_t&I) throw()
00243 : base_t(I), idx(I)
00244 {}
00245
00247 MultiIndex(const MultiIndex<Dims-1>&Midx, const index_t&I) throw()
00248 : base_t(Midx), idx(I)
00249 {}
00250
00256 MultiIndex(const MultiIndex<Dims-1>&m, const index_t&Slice, int SliceDirection)
00257 {
00258 for(int i=0, k=0; i<Dims; i++)
00259 {
00260 if (i==SliceDirection)
00261 (*this)[i] = Slice;
00262 else
00263 {
00264 (*this)[i] = m[k];
00265 k++;
00266 }
00267 }
00268 }
00269
00290 MultiIndex(const index_t&lin_idx, const MultiIndex&Dimens, const CreateFromLinear&) throw()
00291 : base_t( lin_idx, Dimens, CreateFromLinear() )
00292 , idx( (lin_idx / Dimens.subidx().size()) % Dimens.maxidx() )
00293 {}
00294
00296 const MultiIndex<Dims-1>&subidx() const throw()
00297 { return *this; }
00298
00300 MultiIndex<Dims-1>&subidx() throw()
00301 { return *this; }
00302
00308 const index_t&operator[](int i) const
00309 { return (i>=Dims-1) ? idx : subidx()[i]; }
00310
00314 index_t&operator[](int i) throw()
00315 { return (i>=Dims-1) ? idx : subidx()[i]; }
00316
00317
00319 Eagle::Assignment<MultiIndex,0> operator=(const index_t&i)
00320 {
00321 return Eagle::Assignment<MultiIndex, 0>(*this, i);
00322 }
00323
00324 const index_t&maxidx() const throw()
00325 { return idx; }
00326
00327 index_t&maxidx() throw()
00328 { return idx; }
00329
00331 index_t size() const throw()
00332 {
00333 return idx*subidx().size();
00334 }
00335
00337 bool operator!=(const MultiIndex&D) const throw()
00338 {
00339 return idx != D.maxidx() || subidx() != D.subidx() ;
00340 }
00341
00343 bool operator==(const MultiIndex&M) const throw()
00344 {
00345 return idx == M.maxidx() && subidx() == M.subidx() ;
00346 }
00347
00349 MultiIndex&operator+=(const MultiIndex&D) throw()
00350 {
00351 idx += D.maxidx();
00352 subidx() += D.subidx();
00353 return *this;
00354 }
00355
00357 MultiIndex&operator-=(const MultiIndex&D) throw()
00358 {
00359 idx -= D.maxidx();
00360 subidx() -= D.subidx();
00361 return *this;
00362 }
00363
00364
00365
00366
00367
00368
00369
00370
00371
00372
00373
00374
00375
00376
00377
00378
00379
00380
00381
00387 bool inc(const MultiIndex&Dimens) throw()
00388 {
00389 idx ++;
00390 if ( idx >= Dimens.maxidx() )
00391 {
00392 bool inc = subidx().inc(Dimens);
00393 if (inc)
00394 idx = index_t(0);
00395 return inc;
00396 }
00397 return true;
00398 }
00399
00406 bool inc(const MultiIndex&Dimens, const MultiIndex&Increment) throw()
00407 {
00408 idx += Increment.maxidx();
00409 if ( idx >= Dimens.maxidx() )
00410 {
00411 bool inc = subidx().inc( Dimens, Increment );
00412 if (inc)
00413 idx = index_t(0);
00414 return inc;
00415 }
00416 return true;
00417 }
00418
00422 index_t linear(const MultiIndex&Dimens) const throw()
00423 {
00424
00425
00426
00427
00428 assert(idx >= 0 && idx < Dimens.maxidx() );
00429
00430
00431
00432
00433
00434
00435
00436 index_t L = Dimens.subidx().size();
00437 return idx * L + subidx().linear( Dimens.subidx() );
00438 }
00439
00440
00441
00442 #if defined(_IOSTREAM_) || defined(_GLIBCXX_OSTREAM) || _CPP_IOSTREAM || defined(__SGI_STL_IOSTREAM)
00443
00444 std::ostream&print(std::ostream&os, char sep) const
00445 {
00446 base_t::print(os, sep);
00447 return os << sep << idx;
00448 }
00449
00451 friend std::ostream&operator<<(std::ostream&os, const MultiIndex&MI)
00452 {
00453 return MI.print(os, 'x');
00454 }
00455 #endif
00456
00458 bool isWithin(const MultiIndex&Range) const throw()
00459 {
00460
00461 return maxidx() < Range.maxidx() &&
00462 subidx().isWithin( Range );
00463 }
00464
00466 bool operator<(const MultiIndex&Range) const throw()
00467 {
00468 return isWithin( Range );
00469 }
00470
00472 bool operator>(const MultiIndex&Range) const throw()
00473 {
00474 return !isWithin( Range );
00475 }
00476
00478 MultiIndex operator+(const MultiIndex&D) const throw()
00479 { return MultiIndex(*this, D, Add() ); }
00480
00482 MultiIndex operator-(const MultiIndex&D) const throw()
00483 { return MultiIndex(*this, D, Sub() ); }
00484
00486 MultiIndex operator*(const MultiIndex&D) const throw()
00487 { return MultiIndex(*this, D, Mult() ); }
00488
00490 MultiIndex operator/(const MultiIndex&D) const
00491 { return MultiIndex(*this, D, Div() ); }
00492
00493
00503 MultiIndex operator+(int i) const
00504 {
00505 const MultiIndex&M = *this;
00506 if (i==0) return M;
00507
00508 if (i<0)
00509 return M - MultiIndex::Axis(-i+1);
00510 else
00511 return M + MultiIndex::Axis( i-1);
00512 }
00513
00519 friend MultiIndex operator-(const MultiIndex&M, int i)
00520 {
00521 if (i==0) return M;
00522
00523 if (i<0)
00524 return M + MultiIndex::Axis(-i+1);
00525 else
00526 return M - MultiIndex::Axis( i-1);
00527 }
00528
00529 protected:
00530
00531 template <class Functor, int SuperDims>
00532 static bool ForEachRecursion(Functor&F, const MultiIndex<SuperDims>&SuperIndex,
00533 const MultiIndex&Start, const MultiIndex&End, const MultiIndex&Increment)
00534 {
00535 for(index_t idx = Start.idx; idx < End.idx; idx += Increment.idx)
00536 {
00537 MultiIndex<SuperDims+1> SupraIndex(SuperIndex, idx);
00538 if (!MultiIndex<Dims-1>::ForEachRecursion(F,SupraIndex, Start, End, Increment) )
00539 return false;
00540 }
00541 return true;
00542 }
00543
00544 public:
00545
00546 template <class Functor>
00547 static bool ForEach(Functor&F, const MultiIndex&Start, const MultiIndex&End, const MultiIndex&Increment = MultiIndex(1) )
00548 {
00549 for(index_t idx = Start.idx; idx < End.idx; idx += Increment.idx)
00550 {
00551 MultiIndex<1> SupraIndex(idx);
00552 if (!MultiIndex<Dims-1>::ForEachRecursion(F, SupraIndex, Start, End, Increment) )
00553 return false;
00554 }
00555 return true;
00556 }
00557 };
00558
00562 template <>
00563 class MultiIndex<1>
00564 {
00565 index_t idx;
00566
00567 protected:
00568 typedef index_t Index;
00569 typedef Index value_type;
00570
00571 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Add&) throw()
00572 : idx(M.idx + D.maxidx() )
00573 {}
00574
00575 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Sub&) throw()
00576 : idx(M.idx - D.maxidx() )
00577 {}
00578
00579 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Mult&) throw()
00580 : idx(M.idx * D.maxidx() )
00581 {}
00582
00583 MultiIndex(const MultiIndex&M, const MultiIndex&D, const Div&) throw()
00584 : idx(M.idx / D.maxidx() )
00585 {}
00586
00587 MultiIndex(unsigned int bits, const ::Eagle::BinaryAnd&) throw()
00588 : idx( (bits & 1) ? 1 : 0 )
00589 {}
00590
00591 MultiIndex(const MultiIndex&M, const MultiIndex&D, const ::Eagle::BinaryAnd&) throw()
00592 : idx( M.idx > D.maxidx() ? M.idx - D.maxidx() : D.maxidx() - M.idx )
00593 {}
00594
00595 MultiIndex(const MultiIndex&M, const Power2Alignment&) throw()
00596 : idx( RoundUpToNextHighestPowerOf2( M.idx) )
00597 {}
00598
00599 public:
00601 MultiIndex() throw()
00602 : idx(0)
00603 {}
00604
00606 MultiIndex(const Index&I) throw()
00607 : idx(I)
00608 {}
00609
00611 MultiIndex(const MultiIndex&Midx, const Index&I) throw()
00612 : idx(I)
00613 {}
00614
00618 MultiIndex(const Index&lin_idx, const MultiIndex&Dimens, const CreateFromLinear&)
00619 : idx(lin_idx % Dimens.maxidx() )
00620 {}
00621
00622 static int log2(index_t N) throw()
00623 {
00624 if (N&1) return 0;
00625 else return -1;
00626 }
00627
00628 int getSignedOrientation(const MultiIndex&B) const
00629 {
00630 if (idx<B.idx) return +1;
00631 if (idx>B.idx) return -1;
00632 return 0;
00633 }
00634
00635 const Index&operator[](int) const throw()
00636 { return idx; }
00637
00638 Index&operator[](int) throw()
00639 { return idx; }
00640
00641
00642 const Index&maxidx() const throw()
00643 { return idx; }
00644
00645 Index&maxidx() throw()
00646 { return idx; }
00647
00648 Index size() const throw()
00649 { return idx; }
00650
00651
00652 MultiIndex&operator+=(const MultiIndex&D) throw()
00653 {
00654 idx += D.idx;
00655 return *this;
00656 }
00657
00658 MultiIndex&operator-=(const MultiIndex&D) throw()
00659 {
00660 idx -= D.idx;
00661 return *this;
00662 }
00663
00664 MultiIndex operator+(const MultiIndex&D) throw()
00665 {
00666 return MultiIndex(idx + D.idx);
00667 }
00668
00669
00670
00671 friend MultiIndex clamp(const MultiIndex&M, const MultiIndex&D) throw()
00672 {
00673 MultiIndex retval;
00674 retval.idx = clamp(M.idx, D.idx);
00675 return retval;
00676 }
00677
00678 friend MultiIndex clamp_m1(const MultiIndex&M, const MultiIndex&D) throw()
00679 {
00680 MultiIndex retval;
00681 retval.idx = clamp(M.idx, D.idx-1);
00682 return retval;
00683 }
00684
00685
00686 bool operator==(const MultiIndex&D) const throw()
00687 { return idx == D.idx; }
00688
00690 bool operator!=(const MultiIndex&D) const throw()
00691 {
00692 return idx != D.maxidx();
00693 }
00694
00695 bool inc(const MultiIndex&Dimens, const MultiIndex&Increment) throw()
00696 {
00697 idx+= Increment.maxidx();
00698 if (idx >= Dimens.maxidx() )
00699 {
00700 idx = Index(0);
00701 return false;
00702 }
00703 return true;
00704 }
00705
00706 bool inc(const MultiIndex&Dimens) throw()
00707 {
00708 idx++;
00709 if (idx >= Dimens.maxidx() )
00710 {
00711 idx = Index(0);
00712 return false;
00713 }
00714 return true;
00715 }
00716
00717 Index linear(const MultiIndex&) const throw()
00718 { return idx; }
00719
00720 bool isWithin(const MultiIndex&Range) const throw()
00721 {
00722
00723 return maxidx() < Range.maxidx();
00724 }
00725
00726
00727 #if defined(_IOSTREAM_) || defined(_GLIBCXX_OSTREAM) || _CPP_IOSTREAM || defined(__SGI_STL_IOSTREAM)
00728
00729 std::ostream&print(std::ostream&os, char sep) const
00730
00731 { return os << idx; }
00732
00734 friend std::ostream&operator<<(std::ostream&os, const MultiIndex&MI)
00735 {
00736 return MI.print(os, 'x');
00737 }
00738 #endif
00739
00740
00741 template <class Functor, int SuperDims>
00742 static bool ForEachRecursion(Functor&F, const MultiIndex<SuperDims>&SuperIndex,
00743 const MultiIndex&Start, const MultiIndex&End, const MultiIndex&Increment)
00744 {
00745 for(index_t idx = Start.idx; idx < End.idx; idx += Increment.idx)
00746 {
00747 MultiIndex<SuperDims+1> SupraIndex(SuperIndex, idx);
00748 if (!F(SupraIndex) )
00749 return false;
00750 }
00751 return true;
00752 }
00753 };
00754
00755
00756 inline MultiIndex<1> MIndex(index_t i0)
00757 {
00758 MultiIndex<1> m = i0;
00759 return m;
00760 }
00761
00762 inline MultiIndex<2> MIndex(index_t i0, index_t i1)
00763 {
00764 MultiIndex<2> m; m = i0, i1;
00765 return m;
00766 }
00767
00768 inline MultiIndex<3> MIndex(index_t i0, index_t i1, index_t i2)
00769 {
00770 MultiIndex<3> m; m = i0, i1, i2;
00771 return m;
00772 }
00773
00774 inline MultiIndex<4> MIndex(index_t i0, index_t i1, index_t i2, index_t i3)
00775 {
00776 MultiIndex<4> m; m = i0, i1, i2, i3;
00777 return m;
00778 }
00779
00780
00781
00782
00783
00784
00785
00786 #if defined(_IOSTREAM_) || defined(_GLIBCXX_OSTREAM) || _CPP_IOSTREAM || defined(__SGI_STL_IOSTREAM)
00787 template <int Dims>
00788 inline std::ostream&operator<<(std::ostream&os, const MultiIndex<Dims>&M)
00789 {
00790 return M.print( os<<'[', ',') << ']' ;
00791 }
00792 #endif
00793
00794
00799 template <int Dims>
00800 inline FixedArray<double, Dims> div(const MultiIndex<Dims>&Divisor, const MultiIndex<Dims>&Dividend)
00801 {
00802 FixedArray<double, Dims> result;
00803
00804 for(int i=0; i<Dims; i++)
00805 result[i] = double(Divisor[i])/Dividend[i];
00806
00807 return result;
00808 }
00809
00810
00811 #if 0
00812
00816 template <int Dims>
00817 inline FixedArray<double, Dims> Modulo(double*a, double range=1.0)
00818 {
00819 FixedArray<double, Dims> m;
00820 for(int i=0; i<Dims; i++)
00821 m[i] = fmod(a[i] , range);
00822 return m;
00823 }
00824
00825
00830 template <int Dims>
00831 inline FixedArray<double, Dims> Modulo(double*a, double range=1.0)
00832 {
00833 FixedArray<double, Dims> m;
00834 for(int i=0; i<Dims; i++)
00835 m[i] = fmod(a[i] , range);
00836 return m;
00837 }
00838
00843 template <int Dims>
00844 inline FixedArray<double, Dims> Modulo(const FixedArray<double, Dims>&a, const FixedArray<double, Dims>&range)
00845 {
00846 MultiIndex<Dims, double> m;
00847 for(int i=0; i<Dims; i++)
00848 m[i] = fmod(a[i] , range[i]);
00849 return m;
00850 }
00851 #endif
00852
00853
00854 }
00855
00856 #endif // __Fiber_MultiIndex_HPP