Image.h

// ===========================================================================
// Imagine++ Libraries
// Copyright (C) Imagine
// For detailed information: http://imagine.enpc.fr/software
// ===========================================================================

#include "PixelTraits.h"
#include "Interpol.h"
#include "Border.h"

namespace Imagine {

    template <typename T, int dim=2> class Image : public MultiArray<T,dim> {
        typedef MultiArray<T,dim> Base;

    public:
        typedef typename PixelTraits<T>::scalar_type scalar_type;
        typedef typename Base::iterator iterator;
        typedef typename Base::const_iterator const_iterator;

        Image() : Base() {}
        explicit Image(const Coords<dim>& sz) : Base(sz) {}
        Image(int w, int h) : Base(w,h) {}
        Image(int w, int h, int d) : Base(w,h,d) {}
        Image(T* ptr, const Coords<dim>& sz,bool handleDelete=false) : Base(ptr,sz,handleDelete) {}
        Image(T* ptr, int w,int h,bool handleDelete=false) : Base(ptr,w,h,handleDelete) {}
        Image(T* ptr, int w,int h,int d,bool handleDelete=false) : Base(ptr,w,h,d,handleDelete) {}
        Image(const Base& I) : Base(I) {}
        template <typename T2> Image(const MultiArray<T2,dim> &I) : Base(I) {}
        virtual ~Image() {}
        Image& operator=(const Base& I) { Base::operator = (I); return *this; }
        template <typename T2> Image& operator= (const MultiArray<T2,dim> &I) {
            Base::operator = (I);
            return *this;
        }
        Image clone() const { return Base::clone(); }
        Image& fill(T x) { Base::fill(x); return *this; }
        Image getSubImage(const Coords<dim> &offset, const Coords<dim> &sz) const {
            return this->getSubArray(offset,sz);
        }

        Image getSubImage(int x,int y,int w,int h) const {
            return this->getSubArray(Coords<2>(x,y),Coords<2>(w,h));
        }
        Image getSubImage(int x,int y,int z,int w,int h,int d) const {
            return this->getSubArray(Coords<3>(x,y,z),Coords<3>(w,h,d));
        }
        template <typename V> typename PixelTraits<T>::template CastPixel<V>::value_type interpolate(const FVector<V,dim>& c) const { return Imagine::interpolate(*this,c); }
        template <typename V> typename PixelTraits<T>::template CastPixel<V>::value_type interpolate(V x, V y) const { return Imagine::interpolate(*this,FVector<V,2>(x,y)); }
        template <typename V> typename PixelTraits<T>::template CastPixel<V>::value_type interpolate(V x, V y, V z) const { return Imagine::interpolate(*this,FVector<V,3>(x,y,z)); }
        T neumann(const Coords<dim>& c) const { return NeumannBorder<T,dim>()(*this,c); }
        T neumann(int x, int y) const { return NeumannBorder<T,2>()(*this,Coords<2>(x,y)); }
        T neumann(int x, int y, int z) const { return NeumannBorder<T,3>()(*this,Coords<3>(x,y,z)); }
        T dirichlet(const Coords<dim>& c, T out= T(0)) const { return DirichletBorder<T,dim>(out)(*this,c); }
        T dirichlet(int x, int y, T out= T(0)) const { return DirichletBorder<T,2>(out)(*this,Coords<2>(x,y)); }
        T dirichlet(int x, int y, int z, T out= T(0)) const { return DirichletBorder<T,3>(out)(*this,Coords<3>(x,y,z)); }
        T mirror(const Coords<dim>& c) const { return MirrorBorder<T,dim>()(*this,c); }
        T mirror(int x, int y) const { return MirrorBorder<T,2>()(*this,Coords<2>(x,y)); }
        T mirror(int x, int y, int z) const { return MirrorBorder<T,3>()(*this,Coords<3>(x,y,z)); }
        T invMirror(const Coords<dim>& c) const { return InvMirrorBorder<T,dim>()(*this,c); }
        T invMirror(int x, int y) const { return InvMirrorBorder<T,2>()(*this,Coords<2>(x,y)); }
        T invMirror(int x, int y, int z) const { return InvMirrorBorder<T,3>()(*this,Coords<3>(x,y,z)); }
        
        

        friend inline T sum(const Image& I) {
            return std::accumulate(I.begin(), I.end(), T(scalar_type(0)));
        }
        friend inline scalar_type norm2(const Image& I) {
            scalar_type n = 0;
            for (const_iterator it = I.begin() ; it != I.end() ; ++it) {
                n += (*it) * (*it);
            }
            return n;
        }
        friend inline scalar_type norm(const Image& I) {
            assert( !std::numeric_limits<T>::is_integer );
            return ::sqrt(norm2(I));
        }
        Image operator-() const {
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it3)
                *it3 = - *it1;
            return C;
        }
        Image operator+(const Image& B) const {
            assert(this->sizes()==B.sizes());
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin(), it2 = B.begin() ; it1 != this->end() ; ++it1, ++it2, ++it3)
                *it3 = *it1 + *it2;
            return C;
        }
        Image& operator+=(const Image& B) {
            assert(this->sizes()==B.sizes());
            iterator it1 = this->begin();
            for (const_iterator it2 = B.begin() ; it2 != B.end() ; ++it1, ++it2)
                *it1+=*it2;
            return *this;
        }
        Image operator-(const Image& B) const {
            assert(this->sizes()==B.sizes());
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin(), it2 = B.begin() ; it1 != this->end() ; ++it1, ++it2, ++it3)
                *it3 = *it1 - *it2;
            return C;
        }
        Image& operator-=(const Image& B) {
            assert(this->sizes()==B.sizes());
            iterator it1 = this->begin();
            for (const_iterator it2 = B.begin() ; it2 != B.end() ; ++it1, ++it2)
                *it1-=*it2;
            return *this;
        }
        Image operator+(T v) const {
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it3)
                *it3 = *it1 + v;
            return C;
        }
        friend inline Image operator+(T v, const Image& I) { return I+v; }
        Image& operator+=(T v){
            for (iterator it = this->begin() ; it != this->end() ; ++it)
                *it+=v;
            return *this;
        }
        Image operator-(T v) const {
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it3)
                *it3 = *it1 - v;
            return C;
        }
        friend inline Image operator-(T v, const Image& I) { return (-I)+v; }
        Image& operator-=(T v){
            for (iterator it = this->begin() ; it != this->end() ; ++it)
                *it-=v;
            return *this;
        }
        Image operator*(scalar_type s) const {
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it3)
                *it3 = *it1 * s;
            return C;
        }
        friend inline Image operator*(scalar_type s, const Image& I) { return I*s; }
        Image& operator*=(scalar_type s){
            for (iterator it = this->begin() ; it != this->end() ; ++it)
                *it*=s;
            return *this;
        }
        Image operator/(scalar_type s) const {
            Image C(this->sizes());
            iterator it3 = C.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it3)
                *it3 = *it1 / s;
            return C;
        }
        Image& operator/=(scalar_type s){
            for (iterator it = this->begin() ; it != this->end() ; ++it)
                *it/=s;
            return *this;
        }
        Image operator*(const Image<scalar_type,dim>& B) const {
            assert(this->sizes()==B.sizes());
            Image C(this->sizes());
            iterator it3 = C.begin();
            typename Image<scalar_type,dim>::const_iterator it2 = B.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it2, ++it3)
                *it3 = *it1 * *it2;
            return C;
        }
        Image& operator*=(const Image<scalar_type,dim>& B) {
            assert(this->sizes()==B.sizes());
            iterator it1 = this->begin();
            for (typename Image<scalar_type,dim>::const_iterator it2 = B.begin() ; it2 != B.end() ; ++it1, ++it2)
                *it1*=*it2;
            return *this;
        }
        Image operator/(const Image<scalar_type,dim>& B) const {
            assert(this->sizes()==B.sizes());
            Image C(this->sizes());
            iterator it3 = C.begin();
            typename Image<scalar_type,dim>::const_iterator it2 = B.begin();
            for (const_iterator it1 = this->begin() ; it1 != this->end() ; ++it1, ++it2, ++it3)
                *it3 = *it1 / *it2;
            return C;
        }
        Image& operator/=(const Image<scalar_type,dim>& B) {
            assert(this->sizes()==B.sizes());
            iterator it1 = this->begin();
            for (typename Image<scalar_type,dim>::const_iterator it2 = B.begin() ; it2 != B.end() ; ++it1, ++it2)
                *it1/=*it2;
            return *this;
        }
    };

}