From 83a712a772ba873f60889c2d47d5632b2d05120a Mon Sep 17 00:00:00 2001
From: Christian Zimmermann <christian.zimmermann@ur.de>
Date: Mon, 27 Mar 2017 19:29:51 +0200
Subject: [PATCH] finished first draft of index mapping (TO BE TESTED AND
 REFACTORED)

---
 src/base_def.h               | 11 +++++
 src/ma_functional.cc         | 67 +++++++++++++++--------------
 src/ma_functional.h          | 41 ++++++++++++++----
 src/multi_array.cc           |  4 +-
 src/multi_array.h            |  5 ++-
 src/multi_array_header.h     |  2 +
 src/multi_array_operation.cc | 81 +++++++++++++++++++++++++++---------
 src/multi_array_operation.h  | 20 +++++----
 src/range_base.h             |  3 +-
 src/single_range.cc          | 22 +++++-----
 src/single_range.h           |  8 ++--
 11 files changed, 176 insertions(+), 88 deletions(-)

diff --git a/src/base_def.h b/src/base_def.h
index 6552a4b..8e77542 100644
--- a/src/base_def.h
+++ b/src/base_def.h
@@ -108,6 +108,10 @@ namespace MultiArrayTools
     // multi_array_operation.h
     template <typename T, class Range>
     class ConstMultiArrayOperationRoot;
+
+    // multi_array_operation.h
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    class MultiArrayOperationMap;
     
     // multi_array_operation.h
     template <typename T, class Operation, class... MAOps>
@@ -132,6 +136,13 @@ namespace MultiArrayTools
     // manipulator.h
     template <typename T>
     class BinReader;
+
+    // ma_functional.h
+    template <class InRange, class OutRange>
+    class IndexMapFunction;
+
+    // ma_functional.h
+    class vec3d2Function;
 }
 
 #endif
diff --git a/src/ma_functional.cc b/src/ma_functional.cc
index c86247a..5d5e3b8 100644
--- a/src/ma_functional.cc
+++ b/src/ma_functional.cc
@@ -4,49 +4,48 @@
 
 namespace MultiArrayTools
 {
-    namespace
+
+    template <class InRange, class OutRange>
+    IndexMapFunction<InRange,OutRange>::
+    IndexMapFunction(const MultiArrayBase<typename OutRange::IndexType,InRange>& ma,
+		     const OutRange& outRange,
+		     const Name& inName, const Name& outName) : mMap(ma, inName),
+								mOutRange(new OutRange( outRange )),
+								mOIndex(mOutRange->begin())
     {
-
-	template <size_t N>
-	struct MapEvaluation
-	{
-	    template <class OutIndex, class MapTuple>
-	    static void eval(OutIndex& oi, const MapTuple& mt)
-	    {
-		oi.template getIndex<N>() = std::get<N>(mt);
-		MapEvaluation<N-1>::eval(oi, mt);
-	    }
-	};
-	
-	template <>
-	struct MapEvaluation<0>
-	{
-	    template <class OutIndex, class MapTuple>
-	    static void eval(OutIndex& oi, const MapTuple& mt)
-	    {
-		oi.template getIndex<0>() = std::get<0>(mt);
-	    }
-	};
-	
-    } // anonymous namespace
-
+	mOIndex.name(outName);
+    }
     
-    template <class OutIndex, class... Maps>
-    void IndexMapFunction<OutIndex,Maps...>::linkIndicesTo(IndefinitIndexBase* target)
+    template <class InRange, class OutRange>
+    void IndexMapFunction<InRange,OutRange>::linkIndicesTo(IndefinitIndexBase* target)
     {
-	/*!!!!*/
+	mMap.linkIndicesTo(target);
     }
 
-    template <class OutIndex, class... Maps>
-    void IndexMapFunction<OutIndex,Maps...>::eval()
+    template <class InRange, class OutRange>
+    void IndexMapFunction<InRange,OutRange>::eval() const
     {
-	MapEvaluation<sizeof...(Maps)-1>::eval(mOIndex, mMap);
+	mOIndex.copyPos( mMap.get() );
     }
 
-    template <class OutIndex, class... Maps>
-    IndefinitIndexBase& IndexMapFunction<OutIndex,Maps...>::index()
+    template <class InRange, class OutRange>
+    IndefinitIndexBase& IndexMapFunction<InRange,OutRange>::index() const
     {
 	return mOIndex;
     }    
-    
+
+
+    /*
+    vec3d2Function::vec3d2Function(std::shared_ptr<OutRange>& outRange) : mOutRange(outRange),
+									  out(mOutRange->begin()) {} 
+    */
+    vec3d2Function::OutIndex vec3d2Function::operator()(const InIndex& i) const
+    {
+	OutSubIndex& osi = out.template getIndex<0>();
+	
+	osi.atMeta( i.template getIndex<0>().getMetaPos() * i.template getIndex<0>().getMetaPos() +
+		    i.template getIndex<1>().getMetaPos() * i.template getIndex<1>().getMetaPos() +
+		    i.template getIndex<2>().getMetaPos() * i.template getIndex<2>().getMetaPos() );
+	return out;
+    }
 }
diff --git a/src/ma_functional.h b/src/ma_functional.h
index ca87bb6..4b86de6 100644
--- a/src/ma_functional.h
+++ b/src/ma_functional.h
@@ -12,25 +12,52 @@
 namespace MultiArrayTools
 {
 
-    // Maps... are ConstMultiArrayOperationRoots where the corresponding MultiArray defines the map
-    template <class OutIndex, class... Maps>
+    // Map is a ConstMultiArrayOperationRoot where the corresponding MultiArray defines the map
+    template <class InRange, class OutRange>
     class IndexMapFunction
     {
     public:
 
-	typedef std::tuple<Maps...> MapType;
+	typedef OutRange OR;
+	typedef typename OutRange::IndexType OutIndex;
+
+	IndexMapFunction(const MultiArrayBase<typename OutRange::IndexType,InRange>& ma,
+			 const OutRange& outRange,
+			 const Name& inName, const Name& outName);
 	
 	void linkIndicesTo(IndefinitIndexBase* target);
-	void eval();
+	void eval() const;
 	
-	IndefinitIndexBase& index();
+	IndefinitIndexBase& index() const;
 	
     private:
 	
-	MapType mMap;
-	OutIndex mOIndex;
+	ConstMultiArrayOperationRoot<typename OutRange::IndexType,InRange> mMap;
+	std::shared_ptr<OutRange> mOutRange;
+	mutable typename OutRange::IndexType mOIndex;	
     };
 
+    class vec3d2Function
+    {
+    public:
+
+	typedef SingleIndex<int,RangeType::SPACE> CoordIndex;
+	typedef MultiIndex<CoordIndex,CoordIndex,CoordIndex> InIndex;
+	typedef SingleIndex<size_t,RangeType::DISTANCE> OutSubIndex;
+	typedef MultiIndex<OutSubIndex> OutIndex;
+	
+	vec3d2Function() = default;
+	//vec3d2Function(std::shared_ptr<OutRange>& outRange);
+	vec3d2Function(const vec3d2Function& in) = default;
+	vec3d2Function& operator=(const vec3d2Function& in) = default;
+	
+	OutIndex operator()(const InIndex& i) const;
+
+    private:
+	//std::shared_ptr<OutRange> mOutRange;
+	mutable OutIndex out;
+    };
+    
 }
 
 #include "ma_functional.cc"
diff --git a/src/multi_array.cc b/src/multi_array.cc
index 779a24f..91c8ab0 100644
--- a/src/multi_array.cc
+++ b/src/multi_array.cc
@@ -563,12 +563,12 @@ namespace MultiArrayTools
 
     template <typename T, class Range, class Function>
     FunctionalMultiArray<T,Range,Function>::FunctionalMultiArray(const Range& range) :
-	MultiArrayBase<T>(range), mFunc() {}
+	MultiArrayBase<T,Range>(range), mFunc() {}
 
     template <typename T, class Range, class Function>
     FunctionalMultiArray<T,Range,Function>::FunctionalMultiArray(const Range& range,
 								 const Function& func) :
-	MultiArrayBase<T>(range), mFunc(func) {}
+	MultiArrayBase<T,Range>(range), mFunc(func) {}
 
     template <typename T, class Range, class Function>
     const T& FunctionalMultiArray<T,Range,Function>::operator[](const typename Range::IndexType& i) const
diff --git a/src/multi_array.h b/src/multi_array.h
index 0fb9f4d..fe2ffa2 100644
--- a/src/multi_array.h
+++ b/src/multi_array.h
@@ -22,6 +22,8 @@ namespace MultiArrayTools
     {
     public:
 
+	typedef T value_type;
+	
 	class const_iterator : public std::iterator<std::random_access_iterator_tag,T>
 	{
 	public:
@@ -234,7 +236,7 @@ namespace MultiArrayTools
     };
 
     template <typename T, class Range, class Function>
-    class FunctionalMultiArray : public MultiArrayBase<T>
+    class FunctionalMultiArray : public MultiArrayBase<T,Range>
     {
     public:
 	typedef MultiArrayBase<T,Range> MAB;
@@ -252,7 +254,6 @@ namespace MultiArrayTools
     protected:
 	mutable T mVal;
 	Function mFunc;
-	// FUNCTION !!!!
     };
     
 }
diff --git a/src/multi_array_header.h b/src/multi_array_header.h
index bcac116..9b50760 100644
--- a/src/multi_array_header.h
+++ b/src/multi_array_header.h
@@ -14,6 +14,7 @@
 #include "slice.h"
 #include "manipulator.h"
 #include "range_transformer.h"
+#include "ma_functional.h"
 
 namespace MultiArrayTools
 {
@@ -38,6 +39,7 @@ namespace MultiArrayTools
     typedef SingleRange<double,RangeType::ANY> GenericFR;
     typedef SingleRange<size_t,RangeType::LORENTZ> LorentzR;
     typedef SingleRange<int,RangeType::SPACE> Space1dNR;
+    typedef SingleRange<size_t,RangeType::DISTANCE> DistanceNR;
     typedef SingleRange<int,RangeType::MOMENTUM> Mom1dNR;
     typedef SingleRange<size_t, RangeType::ENSEMBLE> EnsR;
     typedef MultiRange<Space1dNR,Space1dNR,Space1dNR> Space3dNR;
diff --git a/src/multi_array_operation.cc b/src/multi_array_operation.cc
index 6fadf0b..d649d3f 100644
--- a/src/multi_array_operation.cc
+++ b/src/multi_array_operation.cc
@@ -175,6 +175,13 @@ namespace MultiArrayTools
 				     MAOps...>(cop, *ind, begin, end, *this, mao...);
     }
 
+    template <typename T, class Range>
+    template<class TotalInRange, class InRange, class OutRange>
+    MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,Range>
+    MultiArrayOperationRoot<T,Range>::map(const IndexMapFunction<InRange,OutRange>& imf)
+    {
+	return MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,Range>(*this, imf);
+    }
     
     template <typename T, class Range>
     template <class MAOp>
@@ -517,10 +524,11 @@ namespace MultiArrayTools
      *   MultiArrayOperationMap     *
      ********************************/
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::
-    MultiArrayOperationMap(MultiArrayOperationRoot<T,OutRange>& root, const MapFunction mf) :
-	MultiArrayOperationBase<T>(),
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::
+    MultiArrayOperationMap(MultiArrayOperationRoot<T,TotalRange>& root,
+			   const IndexMapFunction<InRange,OutRange>& mf) :
+	MutableMultiArrayOperationBase<T>(),
 	mMF(mf),
 	mRoot(root)
 	//mIndex(mArrayRef.beginIndex()),
@@ -531,11 +539,41 @@ namespace MultiArrayTools
 	//mIndex.name(nm);
     }
 	
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    MultiArrayOperationMap& MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::
-    operator=(const ConstMultiArrayOperationRoot<T,InRange>& in)
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>&
+    MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::
+    operator=(const MultiArrayOperationRoot<T,TotalInRange>& in)
     {
-	mIndex = dynamic_cast<typename InRange::IndexType&>( in.getIndex() );
+	mIndex = dynamic_cast<typename TotalInRange::IndexType const&>( in.index() );
+	MAOB::mIibPtr = &mIndex;
+	mNm = in.name();
+	mIndex.name(mNm); // to be sure...
+	mIndex.setPos( mIndex.max() );
+	typename TotalInRange::IndexType endIndex = mIndex;
+
+	// Implement Map Functions !!!!
+	mRoot.linkIndicesTo( &mMF.index() );
+	mMF.linkIndicesTo( &mIndex );
+
+	MultiArray<size_t,TotalRange> cnt(mRoot->range());
+	auto cnto = cnt(mRoot.name(), true);
+	cnto.linkIndicesTo( &mMF.index() );
+	
+	for(mIndex.setPos(0), mMF.eval(); mIndex != endIndex; ++mIndex, mMF.eval()){
+	    get() += in.get();
+	    ++cnto.get();
+	}
+	// CHECK whether T / size_t mixture works!!
+	mRoot /= cnto;
+	return *this;
+    }
+
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>&
+    MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::
+    operator=(const ConstMultiArrayOperationRoot<T,TotalInRange>& in)
+    {
+	mIndex = dynamic_cast<typename InRange::IndexType&>( in.index() );
 	MAOB::mIibPtr = &mIndex;
 	mNm = in.name();
 	mIndex.name(mNm); // to be sure...
@@ -556,39 +594,42 @@ namespace MultiArrayTools
 	}
 	// CHECK whether T / size_t mixture works!!
 	mRoot /= cnto;
+	return *this;
     }
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    size_t MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::argNum() const
+    
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    size_t MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::argNum() const
     {
 	return 1;
     }
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    IndefinitIndexBase* MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    IndefinitIndexBase* MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::
     getLinked(const std::string& name) const
     {
 	return mRoot.getLinked(name);
     }
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    void MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::linkIndicesTo(IndefinitIndexBase* target) const
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    void MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::
+    linkIndicesTo(IndefinitIndexBase* target) const
     {
 	mRoot.linkIndicesTo(target);
     }
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    void MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::setInternalLinks() const
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    void MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::setInternalLinks() const
     { }
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    const T& MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::get() const
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    const T& MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::get() const
     {
 	return mRoot.get();
     }
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
-    T& MultiArrayOperationMap<T,MapFunction,InRange,OutRange>::get()
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
+    T& MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,TotalRange>::get()
     {
 	return mRoot.get();
     }
diff --git a/src/multi_array_operation.h b/src/multi_array_operation.h
index 896207a..bc79481 100644
--- a/src/multi_array_operation.h
+++ b/src/multi_array_operation.h
@@ -18,7 +18,7 @@ namespace MultiArrayTools
     {
     public:
 
-	typedef T ValType;
+	typedef T value_type;
 	
 	MultiArrayOperationBase() /*{ CHECK; }*/ = default;
 	virtual ~MultiArrayOperationBase();
@@ -97,6 +97,9 @@ namespace MultiArrayTools
 		 size_t end,
 		 const MAOps&... mao) const;
 
+	template<class TotalInRange, class InRange, class OutRange>
+	MultiArrayOperationMap<T,InRange,TotalInRange,OutRange,Range>
+	map(const IndexMapFunction<InRange,OutRange>& imf);
 	
 	template <class MAOp>
 	auto operator+(const MAOp& sec) -> decltype(operator()(std::plus<T>(), sec));
@@ -238,15 +241,18 @@ namespace MultiArrayTools
 	Name mNm;
     };
 
-    template <typename T, class MapFunction, class InRange, class OutRange>
+    template <typename T, class InRange, class TotalInRange, class OutRange, class TotalRange>
     class MultiArrayOperationMap : public MutableMultiArrayOperationBase<T>
     {
     public:
 	typedef MultiArrayOperationBase<T> MAOB;
 	
-	MultiArrayOperationMap(MultiArrayOperationRoot<T,OutRange>& root, const MapFunction mf);
-	MultiArrayOperationMap& operator=(const ConstMultiArrayOperationRoot<T,InRange>& in);
+	MultiArrayOperationMap(MultiArrayOperationRoot<T,TotalRange>& root,
+			       const IndexMapFunction<InRange,OutRange>& mf);
 
+	MultiArrayOperationMap& operator=(const MultiArrayOperationRoot<T,TotalInRange>& in);
+	MultiArrayOperationMap& operator=(const ConstMultiArrayOperationRoot<T,TotalInRange>& in);
+	
 	virtual size_t argNum() const override;
 
 	virtual IndefinitIndexBase* getLinked(const std::string& name) const override;
@@ -258,9 +264,9 @@ namespace MultiArrayTools
 	virtual T& get() override;
 	// !!!!
     protected:
-	MapFunction mMF;
-	MultiArrayOperationRoot<T,OutRange>& mRoot;
-	mutable IndexType mIndex; // Index of incoming range
+	IndexMapFunction<InRange,OutRange> mMF;
+	MultiArrayOperationRoot<T,TotalRange>& mRoot;
+	mutable typename TotalInRange::IndexType mIndex; // Index of incoming range
 	Name mNm; // Name of incoming range
     };
     
diff --git a/src/range_base.h b/src/range_base.h
index a96c88d..6c63d13 100644
--- a/src/range_base.h
+++ b/src/range_base.h
@@ -21,7 +21,8 @@ namespace MultiArrayTools
 	LORENTZ = 4,
 	SPIN = 5,
 	ENSEMBLE = 6,
-	VALUE_ERROR = 7
+	VALUE_ERROR = 7,
+	DISTANCE = 8
     };
 
     class MultiRangeType
diff --git a/src/single_range.cc b/src/single_range.cc
index 19c29ff..aab76a8 100644
--- a/src/single_range.cc
+++ b/src/single_range.cc
@@ -96,38 +96,38 @@ namespace MultiArrayTools
 
     // 
     
-    SingleRange<size_t,RangeType::SPACE>::SingleRange(size_t ext) :
-	RangeBase<SingleIndex<size_t,RangeType::SPACE> >(),
+    SingleRange<size_t,RangeType::DISTANCE>::SingleRange(size_t ext) :
+	RangeBase<SingleIndex<size_t,RangeType::DISTANCE> >(),
 	mExt(ext) {}
     
-    size_t SingleRange<size_t,RangeType::SPACE>::get(size_t pos) const
+    size_t SingleRange<size_t,RangeType::DISTANCE>::get(size_t pos) const
     {
 	return pos;
     }
 
-    size_t SingleRange<size_t,RangeType::SPACE>::getMeta(size_t metaPos) const
+    size_t SingleRange<size_t,RangeType::DISTANCE>::getMeta(size_t metaPos) const
     {
 	return metaPos;
     }
     
-    size_t SingleRange<size_t,RangeType::SPACE>::size() const
+    size_t SingleRange<size_t,RangeType::DISTANCE>::size() const
     {
 	return mExt;
     }
     
-    MultiRangeType SingleRange<size_t,RangeType::SPACE>::type() const
+    MultiRangeType SingleRange<size_t,RangeType::DISTANCE>::type() const
     {
-	return MultiRangeType(RangeType::SPACE);
+	return MultiRangeType(RangeType::DISTANCE);
     }
    
-    SingleIndex<size_t,RangeType::SPACE> SingleRange<size_t,RangeType::SPACE>::begin() const
+    SingleIndex<size_t,RangeType::DISTANCE> SingleRange<size_t,RangeType::DISTANCE>::begin() const
     {
-	return SingleIndex<size_t,RangeType::SPACE>(this, 0);
+	return SingleIndex<size_t,RangeType::DISTANCE>(this, 0);
     }
     
-    SingleIndex<size_t,RangeType::SPACE> SingleRange<size_t,RangeType::SPACE>::end() const
+    SingleIndex<size_t,RangeType::DISTANCE> SingleRange<size_t,RangeType::DISTANCE>::end() const
     {
-	return SingleIndex<size_t,RangeType::SPACE>(this, size());
+	return SingleIndex<size_t,RangeType::DISTANCE>(this, size());
     }
 
     // 
diff --git a/src/single_range.h b/src/single_range.h
index 8971e2e..b1ca52e 100644
--- a/src/single_range.h
+++ b/src/single_range.h
@@ -113,10 +113,10 @@ namespace MultiArrayTools
     };
 
     template <>
-    class SingleRange<size_t,RangeType::SPACE> : public RangeBase<SingleIndex<size_t,RangeType::SPACE> >
+    class SingleRange<size_t,RangeType::DISTANCE> : public RangeBase<SingleIndex<size_t,RangeType::DISTANCE> >
     {
     public:
-	typedef typename RangeBase<SingleIndex<int,RangeType::SPACE> >::IndexType IndexType;
+	typedef typename RangeBase<SingleIndex<size_t,RangeType::DISTANCE> >::IndexType IndexType;
 	
 	DEFAULT_MEMBERS(SingleRange);
 
@@ -129,8 +129,8 @@ namespace MultiArrayTools
 
 	virtual MultiRangeType type() const override;
 	
-	SingleIndex<size_t,RangeType::SPACE> begin() const override;
-	SingleIndex<size_t,RangeType::SPACE> end() const override;
+	SingleIndex<size_t,RangeType::DISTANCE> begin() const override;
+	SingleIndex<size_t,RangeType::DISTANCE> end() const override;
 	
     protected:
 	size_t mExt;