DirectProductArray.hpp

#ifndef __FIBER_DIRECTPRODUCT_HPP
#define __FIBER_DIRECTPRODUCT_HPP
 
#ifndef DONT_USE_PRAGMA_ONCE
#pragma once
#endif
 
#include "MemArray.hpp"
#include "CoordinateAxisArray.hpp"
 
namespace Fiber
{
using std::type_info;
 
 
template <class T, class Category, int MDIMS>
class   DirectProductMemArrayBase : public MemArrayBase<MDIMS>
{
protected:
 
        ~DirectProductMemArrayBase()
        {}
 
public:
 
#ifdef  HAVE_OLD_MEMARRAY
        /*
         * Silences the compiler of overloaded virtual hiding, which is not the case here
         */
        using MemArrayBase<MDIMS>::newMemArray;
//    using MemArrayBase<MDIMS>::createMemArray;
#endif
 
        typedef typename T::value_type component_type;
 
        typedef CoordinateAxisArray<component_type, Category> CoordinateAxisArray_t;
 
        enum
        {
                COMPONENTS = T::SIZE,
                Dims = MDIMS
        };
 
        RefPtr<CoordinateAxisArray_t>   Components[ COMPONENTS ];
 
        Eagle::Vector<int, Dims>        ComponentMap;
 
        typedef FixedArray<double, MDIMS> FloatIndex_t;
 
        DirectProductMemArrayBase(const MemBase::Creator_t&C = MemCore::NullPtr() )
        : MemBase(FiberType<T>::getFiberType()->getFiberTypeIndex(), C)
        , MemArrayBase<MDIMS>(FiberType<T>::getFiberType()->getFiberTypeIndex(), C )
        {
                for(int i=0; i<Dims; i++)
                {
                        if (i<COMPONENTS)
                                ComponentMap[ i ] = i;
                        else
                                ComponentMap[ i ] = -1;
                }
        }
 
        DirectProductMemArrayBase(const RefPtr<CoordinateAxisArray_t> C[ COMPONENTS ],
                                  const Eagle::Vector<int, Dims> &CMap,
                                  const MemBase::Creator_t&theCreator = MemCore::NullPtr() )
        : MemArrayBase<MDIMS>(FiberType<T>::getFiberType()->getFiberTypeIndex(), theCreator)
        {
                for(int i=0; i<COMPONENTS; i++)
                {
                        Components[i] = C[i];
                }
 
                for(int i=0; i<Dims; i++)
                {
                        ComponentMap[ i ] = CMap[i];
                }
        }
 
        DirectProductMemArrayBase(const DirectProductMemArrayBase&Src, const MemBase::Creator_t&theCreator = MemCore::NullPtr() )
        : MemArrayBase<MDIMS>(FiberType<T>::getFiberType()->getFiberTypeIndex(), theCreator)
        {
                for(int i=0; i<COMPONENTS; i++)
                {
                        if (Src.Components[i] )
                                Components[i] = new CoordinateAxisArray_t( *Src.Components[i] );
                }
                ComponentMap = Src.ComponentMap;
        }
 
        void    DeferredConstructor() override
        {
                MemCore::SaveRegistry<DirectProductMemArrayBase>::addInterfaces(this);
        }
 
        MemCore::memsize_t      memsize() const override
        {
                return 0;
        }
 
#ifdef  HAVE_OLD_MEMARRAY
        //TODO: was marked as override, but doesn't
        RefPtr<MemBase> newMemArray() const
        {
                return new MemArray<MDIMS, T>( this->Size() );
        }
#endif
 
        bool    copyFrom(const RefPtr<MemCore::ChunkBase>&theChunk) override
        {
                return  false;
        }
 
        RefPtr<MemArray<MDIMS,T> > createMemArray(const RefPtr<MemCore::ChunkBase>&Data, const MemBase::Creator_t&C) const
        {
                if (!Data) return MemCore::NullPtr();
                if (Data->elements() != Size().size() )
                        return MemCore::NullPtr();
 
                return new MemArray<MDIMS,T>(Data, Size(), C );
        }
 
        RefPtr<MemBase> createEquallySizedMemArray(const RefPtr<MemCore::ChunkBase>&Storage, const MemBase::Creator_t&C) const override
        {
                return createMemArray(Storage, C);
        }
 
        RefPtr<MemBase> getComponentArray(int member, const MemBase::Creator_t&C) override
        {
                return MemCore::NullPtr();
        }
 
        void    retrieve(T&retval, const MultiIndex<Dims>&I) const
        {
                for(int i=0; i<Dims; i++)
                {
                int     component = ComponentMap[ i ];
                        if (component<0 || component>=COMPONENTS)
                                continue;
 
                const RefPtr<CoordinateAxisArray_t>&that = Components[ component ];
 
                        if (that)
                                retval[ component ] = (*that)( I[i] );
                }
        }
 
        FloatIndex_t    getIndex(const T&PointCoordinates) const
        {
        FloatIndex_t    result;
                for(int i=0; i<Dims; i++)
                {
                int     component = ComponentMap[ i ];
                        if (component<0 || component>=COMPONENTS)
                                continue;
 
                        if (const RefPtr<CoordinateAxisArray_t>&that = Components[ component ] )
                        {
                                result[ component ] = that->reverse( PointCoordinates[i] );
                        }
                }
                return result;
        }
 
#if 0 // that->locate( .. ) is disabled
        bool    locate(T&pos) const
        {
                for(int i=0; i<Dims; i++)
                {
                int     component = ComponentMap[ i ];
                        if (component<0 || component>=COMPONENTS)
                        {
                                continue;
                        }
 
                const RefPtr<CoordinateAxisArray_t>&that = Components[ component ];
 
                        if (!that)
                                return false;
 
                        if (!that->locate( pos[ component ] ) )
                                return false;
                }
 
                return true;
        }
#endif
        MultiIndex<Dims> Size() const override
        {
        MultiIndex<Dims> Sz;
                for(int i=0; i<Dims; i++)
                {
                int     component = ComponentMap[ i ];
                        if (component<0 || component>=COMPONENTS)
                        {
                                Sz[ i ] = 1;
                                continue;
                        }
                const RefPtr<CoordinateAxisArray_t>&that = Components[ component ];
                        Sz[i] = that->Length();
                }
                return Sz;
        }
 
        MultiIndex<Dims> firstIndex() const
        {
        MultiIndex<Dims> I;
                for(int i=0; i<Dims; i++)
                        I[i] = 0;
 
                return I;
        }
 
        MultiIndex<Dims> lastIndex() const
        {
        MultiIndex<Dims> I;
        const MultiIndex<Dims>&Sz = this->Size();
                for(int i=0; i<Dims; i++)
                {
                index_t sz = Sz[i];
                        if (sz>1)
                                I[i] = sz - 1;
                        else
                                I[i] = 0;
                }
                return I;
//              return Size() - MIndex(1,1,1);
        }
 
//      CreativeIterator<T>* creativeIterator() override
//      {
//              return 0;
//      }
 
/*
        const type_info&getType() const override
        {
                return typeid(T);
        }
 
        const FiberTypeBase&getFiberType() const override
        {
                return FiberType<T>::getFiberType();
        }
*/
 
 
        RefPtr<MemBase> getSlice(index_t i, const MemBase::Creator_t&C) const override
        {
                return MemCore::NullPtr(); // not yet implemented...
        }
 
        void*getPtr() override
        {
                return 0;
        }
 
        HyperslabParameters& getHyperslabParameters() override
        {
        static  HyperslabParameters nix(0,0,0);
                return nix;
        }
/*
        RefPtr<UntypedCreativeIterator> getIterator() const
        {
                return MemCore::NullPtr() override;
        }
*/
};
 
struct DirectProductMemArrayCreateBaseArray {};
 
template <class T, class Category = LinearCoordinateAxisArray, int MDIMS = T::SIZE>
class   DirectProductMemArray : public DirectProductMemArrayBase<typename T::FixedArray_t, Category, MDIMS>, public TypedArray<T>
{
protected:
        ~DirectProductMemArray()
        {}
 
public:
        typedef MemArray<MDIMS, T> MemArray_t;
 
        typedef T value_type;
 
        typedef typename T::value_type component_type;
 
        typedef DirectProductMemArrayBase<typename T::FixedArray_t, Category, MDIMS> DirectProductMemArrayBase_t;
        using DirectProductMemArrayBase_t::Size;
 
        enum
        {
                COMPONENTS = DirectProductMemArrayBase_t::COMPONENTS,
                Dims       = DirectProductMemArrayBase_t::Dims
        };
 
//      T       origin;
 
        template <class TBase>
        DirectProductMemArray(const DirectProductMemArrayCreateBaseArray&,const DirectProductMemArray<TBase, Category>&D, const MemBase::Creator_t&C = MemCore::NullPtr() )
        : DirectProductMemArrayBase_t(D.Components, D.ComponentMap, C)
        , TypedArray<T>( C )
//      , origin(D.origin)
        {}
 
        DirectProductMemArray(const MemBase::Creator_t&C = MemCore::NullPtr() )
        : MemBase(FiberType<T>::getFiberType()->getFiberTypeIndex(), C) // don't use creators here, these arrays don't take significant memory
        , TypedArray<T>( C )
        {}
 
        DirectProductMemArray(const DirectProductMemArray&DPM, const MemBase::Creator_t&C = MemCore::NullPtr() )
        : MemBase(FiberType<T>::getFiberType()->getFiberTypeIndex(), C)
        , DirectProductMemArrayBase_t(DPM, C)
        , TypedArray<T>( C )
        {}
 
        T operator[](const MultiIndex<Dims>&I) const
        {
        T       retval;// = origin;
                this->retrieve( retval, I );
                return retval;
        }
 
        T       getElement(const index_t i) const
        {
        const MultiIndex<Dims> I(i, this->Size(), CreateFromLinear{} );
                return (*this)[I];
        }
 
        std::string     getElementAsString(size_t index) const override
        {
        const MultiIndex<Dims> I(index, this->Size(), CreateFromLinear{} );
                if (I.isWithin( this->Size() ))
                        return std::to_string( (*this)[I] );
 
                return {};
        }
 
        RefPtr<MemBase> copy(const MemBase::Creator_t&C = nullptr) const override
        {
                return new DirectProductMemArray(*this, C);
        }
 
 
        RefPtr<MemBase> makeMemArray(const MemBase::Creator_t&C) const override
        {
        MultiIndex<Dims> Sz = this->Size();
 
        RefPtr<MemArray_t> M = new MemArray_t( Sz, C );
        typename MemArray_t::MultiArray_t&dst = *M;
 
        MultiIndex<Dims> idx;
 
        T       val;// = origin;
                do
                {
                        this-> retrieve( val, idx);
                        dst[idx] = val;
 
/*
  cannot use retrieve(dst[idx]) here, because dst[idx] might be an ElementProxy type,
  and passing such as non-const reference doesnt work.
  would be interesting how to fix that. if it is possible at all.
 */
//                      retrieve( dst[idx], idx);
                }
                while( idx.inc( Sz ) );
 
                return M;
        }
 
        RefPtr<MemBase> newMemArraynD(const MultiIndex<Dims>&, const MemBase::Creator_t&C) const override
        {
                return MemCore::NullPtr();
        }
 
        RefPtr<MemBase> createSubMemArray(const MultiIndex<Dims>&, const MultiIndex<Dims>&, const MemBase::Creator_t&C) const override
        {
                return MemCore::NullPtr();
        }
 
#if 0
        MemCore::WeakPtr<FiberTypeBase> getFiberType() const override
        {
                return FiberType<T>::getFiberType();
        }
#endif
 
        const type_info&getType() const override
        {
                return typeid(T);
        }
 
        CreativeIterator<T>* creativeIterator() override
        {
                return nullptr;
        }
 
        const CreativeIterator<T>* creativeIterator() const override
        {
                return  nullptr;
        }
 
        RefPtr<MemCore::TypedChunk<T> > myChunk() const override
        {
                return  nullptr;
        }
 
        const void*OriginPtr() const override { return nullptr; }
        const void*DeltaPtr()  const override { return nullptr; }
 
 
        T first() const override
        {
                return (*this)[ this->firstIndex() ];
        }
 
        T last() const override
        {
                return (*this)[ this->lastIndex() ];
        }
 
};
 
 
template <class T, class DeltaT = T, int MDIMS = T::SIZE>
class   LinearDirectProductMemArray : public DirectProductMemArray<T, LinearCoordinateAxisArray, MDIMS>
{
protected:
        ~LinearDirectProductMemArray()
        {}
 
public:
        typedef DirectProductMemArray<T, LinearCoordinateAxisArray, MDIMS> Base_t;
 
        typedef typename Base_t::CoordinateAxisArray_t CoordinateAxisArray_t;
 
        using Base_t::Components;
 
        enum { COMPONENTS = T::SIZE };
 
        typedef typename T::value_type component_type;
 
        const T baseValue;
        const DeltaT deltaValue;
 
//      virtual const T&Origin() const = 0;
//      virtual const T&Delta()  const = 0;
 
 
        const void*OriginPtr() const override { return &Origin(); }
        const void*DeltaPtr()  const override { return &Delta();  }
 
        const T&Origin()      const { return baseValue; }
        const DeltaT&Delta()  const { return deltaValue; }
 
        LinearDirectProductMemArray(const index_t dims[MDIMS], const T&base, const DeltaT&delta, const MemBase::Creator_t&C = MemCore::NullPtr() )
        : MemBase(FiberType<T>::getFiberType()->getFiberTypeIndex(), C)
        , DirectProductMemArray<T, LinearCoordinateAxisArray, MDIMS>(C)
        , baseValue(base)
        , deltaValue(delta)
        {
//              puts("LinearDirectProductMemArray(const index_t dims[MDIMS], const T&base, const DeltaT&delta, const MemBase::Creator_t&C = MemCore::NullPtr() )");fflush(stdout);
 
                for(int c=0; c<COMPONENTS; c++)
                {
                        if (c<MDIMS)
                                this->Components[c] = new CoordinateAxisArray_t( base[c], delta[c], dims[c] );
                }
        }
 
        LinearDirectProductMemArray(const MultiIndex<MDIMS>&dims, const T&base, const DeltaT&delta, const MemBase::Creator_t&C = MemCore::NullPtr() )
        : MemBase(FiberType<T>::getFiberType()->getFiberTypeIndex(), C)
        , DirectProductMemArray<T, LinearCoordinateAxisArray, MDIMS>(C)
        , baseValue(base)
        , deltaValue(delta)
        {
//              puts("LinearDirectProductMemArray(const MultiIndex<MDIMS>&dims, const T&base, const DeltaT&delta, const MemBase::Creator_t&C = MemCore::NullPtr() )");fflush(stdout);
 
                for(int c=0; c<COMPONENTS; c++)
                {
                        if (c<MDIMS)
                        {
                                this->Components[c] = new CoordinateAxisArray_t( base[c], delta[c], dims[c] );
                        }
                }
        }
 
        LinearDirectProductMemArray(const LinearDirectProductMemArray&Src, const MemBase::Creator_t&C = MemCore::NullPtr() )
        : MemBase(FiberType<T>::getFiberType()->getFiberTypeIndex(), C)
        , DirectProductMemArray<T, LinearCoordinateAxisArray, MDIMS>(Src,C)
        {}
 
#if     0
        bool resize(int component, index_t newSize)
        {
                if (RefPtr<LinearArray<component_type> > LA = Components[component])
                {
                        LA->resize( newSize );
                        return true;
                }
                return false;
        }
#endif
 
        bool  getFlooredMultiIndex( const T& pos, MultiIndex<MDIMS>& index ) const
        {
//              if( info )
//              {
//              Verbose( 0 ) << "LinearDirectProductMemArray::getFlooredMultiIndex() ";
//              Verbose( 0 ) << "'crds: " << pos << ", origin: " << baseValue << ", delta: " << deltaValue;
//              }
 
                for( unsigned i = 0; i < MDIMS; ++i )
                {
                        if ( (pos[i] - baseValue[i]) / deltaValue[i] >= 0)
                        {
                                index[i] = size_t( floor( ( pos[i] - baseValue[i] ) / deltaValue[i] ) );
                        }
                        else
                        {
                                Verbose( 51 ) << "LinearDirectProductMemArray::getFlooredMultiIndex() Fragment Index " << i << " negative";
                                Verbose( 51 ) << "'crds: " << pos << ", origin: " << baseValue << ", delta: " << deltaValue;
                                return false;
                        }
 
 
                int     component = this->ComponentMap[ i ];
                        if( index[i] >= this->Components[ component ]->Length() )
                        {
                                Verbose( 51 ) << "LinearDirectProductMemArray::getFlooredMultiIndex() Fragment Index " << i << " too large: " << index[i] << ">=" << this->ComponentMap[i];
                                Verbose( 51 ) << "'crds: " << pos << ", origin: " << baseValue << ", delta: " << deltaValue;
                                return false;
                        }
                }
 
                return true;
        }
 
        bool  getNearestMultiIndex( const T& pos, MultiIndex<MDIMS>& index ) const
        {
                for( unsigned i = 0; i < MDIMS; ++i )
                {
                        if ( (pos[i] - baseValue[i]) / deltaValue[i] >= 0)
                        {
                                index[i] = size_t( round( ( pos[i] - baseValue[i] ) / deltaValue[i] ) );
                        }
                        else
                        {
                                Verbose(51) << "LinearDirectProductMemArray::getNearestMultiIndex() Fragment Index " << i << " negative";
                                Verbose(51) << "'crds: " << pos << ", origin: " << baseValue << ", delta: " << deltaValue;
                                return false;
                        }
 
                int     component = this->ComponentMap[ i ];
                        if( index[i] >= this->Components[ component ]->Length() )
                        {
                                Verbose(51) << "LinearDirectProductMemArray::getNearestMultiIndex() Fragment Index " << i << " too large: " << index[i] << ">=" << this->ComponentMap[i];
                                Verbose(51) << "'crds: " << pos << ", origin: " << baseValue << ", delta: " << deltaValue;
                                return false;
                        }
                }
 
                return true;
        }
};
 
} /* namespace Fiber */
 
#endif /* __FIBER_DIRECTPRODUCT_HPP */