/home/oan/prosjekt/gotools/segmentation/gpl_distro/lsseg_1.0_gpl/include/Image.h

Go to the documentation of this file.

//===========================================================================
// The Level-Set Segmentation Library (LSSEG)
//
//
// Copyright (C) 2000-2005 SINTEF ICT, Applied Mathematics, Norway.
//
// This program is free software; you can redistribute it and/or          
// modify it under the terms of the GNU General Public License            
// as published by the Free Software Foundation version 2 of the License. 
//
// This program is distributed in the hope that it will be useful,        
// but WITHOUT ANY WARRANTY; without even the implied warranty of         
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the          
// GNU General Public License for more details.                           
//
// You should have received a copy of the GNU General Public License      
// along with this program; if not, write to the Free Software            
// Foundation, Inc.,                                                      
// 59 Temple Place - Suite 330,                                           
// Boston, MA  02111-1307, USA.                                           
//
// Contact information: e-mail: tor.dokken@sintef.no                      
// SINTEF ICT, Department of Applied Mathematics,                         
// P.O. Box 124 Blindern,                                                 
// 0314 Oslo, Norway.                                                     
// 
//
// Other licenses are also available for this software, notably licenses
// for:
// - Building commercial software.                                        
// - Building software whose source code you wish to keep private.        
//
//===========================================================================
//===========================================================================
//                                                                           
// File: Image.h                                                             
//                                                                           
// Created: Tue Oct 25 13:57:00 2005                                         
//                                                                           
// Author: Odd A. Andersen <Odd.Andersen@sintef.no>
//                                                                           
// Revision: $Id: Image.h,v 1.23 2006/11/22 15:23:00 oan Exp $
//                                                                           
// Description:
//                                                                           
//===========================================================================

#ifndef _IMAGE_H
#define _IMAGE_H

#include <algorithm>
#include <assert.h>
#include <vector>
#include <istream>
#include <ostream>
#include <numeric>
#include <boost/lambda/lambda.hpp>

using namespace boost;
using namespace boost::lambda;

namespace lsseg {

//===========================================================================

template<typename T>
class Image
//===========================================================================
{
public:
    typedef T value_type;

    Image() 
    {std::fill(dim_, dim_ + 4, 0);}

    Image(int x, int y, int z = 1, int chan = 1, T fill_val = 0) 
        : data_(x * y * z * chan, fill_val)
    { 
        dim_[0] = x;
        dim_[1] = y;
        dim_[2] = z;
        dim_[3] = chan;
    }
    
    virtual ~Image() {}

    template<typename U>
    Image(const Image<U>& rhs) {
        dim_[0] = rhs.dimx();
        dim_[1] = rhs.dimy();
        dim_[2] = rhs.dimz();
        dim_[3] = rhs.numChannels();
        data_.resize(size());
        assert(size() == rhs.size());
        copy(rhs.begin(), rhs.end(), &data_[0]);
    }

    Image(const Image<T>& rhs) : data_(rhs.data_) 
    { std::copy(rhs.dim_, rhs.dim_ + 4, dim_);}

    template<typename U>
    Image(const Image<U>& rhs, bool copy) {
        dim_[0] = rhs.dimx();
        dim_[1] = rhs.dimy();
        dim_[2] = rhs.dimz();
        dim_[3] = rhs.numChannels();
        data_.resize(size());
        if (copy) {
            T* target = begin();
            for (const U* it = rhs.begin(); it != rhs.end(); ++it, ++target) {
                *target = T(*it);
            }
        }
    }
    
    virtual Image<T>& operator=(const Image<T>& rhs) {
        Image<T> tmp(rhs);
        this->swap(tmp);
        return *this;
    }

    Image<T>& operator=(const T& val) {
        std::fill(data_.begin(), data_.end(), val);
        return *this;
    }

    unsigned long int size() const { return dim_[0] * dim_[1] * dim_[2] * dim_[3];}

    unsigned long int channelSize() const { return size() / numChannels();}

    unsigned long int graySize2D() const { return dimx() * dimy();}

    int dimx() const { return dim_[0];}

    int dimy() const { return dim_[1];}
    
    int dimz() const { return dim_[2];}

    int numChannels() const { return dim_[3];}

    template<typename U>
    void resize(const Image<U>& rhs) {
        resize(rhs.dimx(), rhs.dimy(), rhs.dimz(), rhs.numChannels());
    }
    
    virtual void resize(int x, int y, int z = 1, int channels = 1) 
    {
        dim_[0] = x;
        dim_[1] = y;
        dim_[2] = z;
        dim_[3] = channels;
        data_.resize(size());
    }
    
    virtual void swap(Image<T>& rhs) {
        std::swap_ranges(dim_, dim_ + 4, rhs.dim_);
        data_.swap(rhs.data_);
    }

    Image<T>& operator +=(const Image<T>& rhs) {
        assert(size_compatible(rhs));
        const T* from = rhs.begin();
        T* to;
        for (to = begin(); to != end(); ++to, ++from) {
            *to += *from;
        }
        return *this;
    }

    Image<T>& operator -=(const Image<T>& rhs) {
        assert(size_compatible(rhs));
        const T* from = rhs.begin();
        T* to;
        for (to = begin(); to != end(); ++to, ++from) {
            *to -= *from;
        }
        return *this;
    }

    Image<T>& operator +=(T a) {
        T* it;
        for (it = begin(); it != end(); ++it) {
            *it += a;
        }
        return *this;
    }

    Image<T>& operator -=(T a) {
        return operator+=(-a);
    }

    Image<T>& operator *=(T a) {
        T* it;
        for (it = begin(); it != end(); ++it) {
            *it *= a;
        }
        return *this;
    }

    Image<T>& operator /=(T a) {
        const T denom_inv = T(1) / a;
        return (*this) *= denom_inv;
    }

    T max() const { return *std::max_element(data_.begin(), data_.end()); }

    T min() const { return *std::min_element(data_.begin(), data_.end()); }

    template<typename U>
    bool spatial_compatible(const Image<U>& rhs) const {
        if (dimx() != rhs.dimx()) return false;
        if (dimy() != rhs.dimy()) return false;
        if (dimz() != rhs.dimz()) return false;
        return true;
    }


    template<typename U>
    bool size_compatible(const Image<U>& rhs) const {
        if (!spatial_compatible(rhs)) return false;
        if (numChannels() != rhs.numChannels()) return false;
        return true;
    }

    void fill(T val) { std::fill(data_.begin(), data_.end(), val); }

    void setAbsolute() {
        for (T* p = begin(); p != end(); ++p) {
            if (*p < 0) *p = -(*p);
        }
    }

    T getAverage() const
    { return std::accumulate(data_.begin(), data_.end(), T(0)) / T(size()); }

    T* begin() { return &data_[0];}

    T* end() { return &data_[0] + size();}

    const T* begin() const { return &data_[0];}

    const T* end() const { return &data_[0] + size();}

    T* channelBegin(int channel) { return &data_[0] + channel * channelSize();}

    T* channelEnd(int channel) { return channelBegin(channel) + channelSize();}

    const T* channelBegin(int channel) const { return &data_[0] + channel * channelSize();}

    const T* channelEnd(int channel) const  { return channelBegin(channel) + channelSize();}

    void superpose(const Image<T>& img) {
        assert(dimx() == img.dimx() && dimy() == img.dimy()); // at least this must be correct
        assert(img.size() <= size());
        T* to = begin();
        while (to != end()) {
            std::transform(img.begin(), img.end(), to, to, _1 + _2);
            to += img.size();
        }
    }

    void intersect(const Image<T>& img) {
        assert(size_compatible(img));
        std::transform(img.begin(), img.end(), begin(), begin(), _1 && _2);
    }

    void dumpInfo(std::ostream& os) const {
        os << "-------------\n";
        os << "X: " << data_[0] << '\n';
        os << "Y: " << data_[1] << '\n';
        os << "Z: " << data_[2] << '\n';
        os << "Channels: " << data_[3] << std::endl;
    }

    size_t indexOf(int x, int y) const {
        assert(x >= 0 && y >= 0);
        assert(x <= dimx() && y <= dimy());
        return x + dim_[0] * y;
    }

    size_t indexOf(int x, int y, int z) const {
        assert(x >= 0 && y >= 0 && z >= 0);
        assert(x <= dimx() && y <= dimy() && z <= dimz());
        return x + dim_[0] * (y + dim_[1] * z);
    }

    size_t indexOf(int x, int y, int z, int c) const {
        assert(x >= 0 && y >= 0 && z >= 0 && c >= 0);
        assert(x <= dimx() && y <= dimy() && z <= dimz() && c <= numChannels());
        return x + dim_[0] * (y + dim_[1] * (z + dim_[2] * c));
    }

    T& operator()(int x, int y) { return data_[indexOf(x, y)];}

    T& operator()(int x, int y, int z) { return data_[indexOf(x, y, z)];}

    T& operator()(int x, int y, int z, int c) { return data_[indexOf(x, y, z, c)];}

    const T& operator()(int x, int y) const { return data_[indexOf(x, y)];}

    const T& operator()(int x, int y, int z) const { return data_[indexOf(x, y, z)];}

    const T& operator()(int x, int y, int z, int c) const { return data_[indexOf(x, y, z, c)];}

    T& operator[](unsigned int ix) {
        return data_[ix];
    }

    const T& operator[](unsigned int ix) const {
        return data_[ix];
    }

    virtual void permute(const int* const perm)
    {
        Image<T> tmp_img(dim_[perm[0]], dim_[perm[1]], dim_[perm[2]], dim_[perm[3]]);
        int ix[4];

        const T* dp = begin();
        for (ix[3] = 0; ix[3] != numChannels(); ++ix[3]) {
            for (ix[2] = 0; ix[2] != dimz(); ++ix[2]) {
                for (ix[1] = 0; ix[1] != dimy(); ++ix[1]) {
                    for (ix[0] = 0; ix[0] != dimx(); ++ix[0]) {
                        const int Xperm = ix[perm[0]];
                        const int Yperm = ix[perm[1]];
                        const int Zperm = ix[perm[2]];
                        const int Cperm = ix[perm[3]];
                        tmp_img.data_[tmp_img.indexOf(Xperm, Yperm, Zperm, Cperm)] = *dp++;
                    }
                }
            }
        }
        swap(tmp_img);
        
//      Go::GoBorrowedMVGrid<4, T> tmp_grid = getGrid();
//      tmp_grid.permuteElements(perm);
//      for (int i = 0; i < 4; ++i) {
//          dim_[i] = tmp_grid.rowlength(i);
//      }
    }

    void write(std::ostream& os, bool binary=true)  const
    {
        if (binary) {
            // writing to file in binary format

            // platform-dependent stuff
            assert(sizeof(double) == 8);
            assert(sizeof(float) == 4);
            assert(sizeof(int) == 4);

            const char* dp = reinterpret_cast<const char*>(dim_);
            int dim_storage_size = 4 * sizeof(int);
            os.write(dp, dim_storage_size);

            int total_data_size = size() * sizeof(T);
            dp = reinterpret_cast<const char*>(&data_[0]);
            os.write(dp, total_data_size);

        } else {
            // writing to file in ASCII format
            for (int d = 0; d != 4; ++d) {
                os << dim_[d] << '\n';
            }
            for (size_t i = 0; i != data_.size(); ++i) {
                os << data_[i] << ' ';
            }
            os << std::endl;
        }
    }


    virtual void read(std::istream& is, bool binary=true) 
    {
        Image<T> tmp_img;

        if (binary) {
            // read in binary format

            // platform-dependent stuff
            assert(sizeof(double) == 8);
            assert(sizeof(float) == 4);
            assert(sizeof(int) == 4);

            int dim[4];
            int dim_storage_size = 4 * sizeof(int);
            char* dp = reinterpret_cast<char *>(dim);
            is.read(dp, dim_storage_size);
            tmp_img.resize(dim[0], dim[1], dim[2], dim[3]);

            int total_data_size = tmp_img.size() * sizeof(T);
            dp = reinterpret_cast<char*>(tmp_img.begin());
            is.read(dp, total_data_size);
            
        } else {
            // read in ascii format
            int dim[4];
            for (int i = 0; i != 4; ++i) {
                is >> dim[i];
            }
            tmp_img.resize(dim[0], dim[1], dim[2], dim[3]);
            T* dpointer = tmp_img.begin();
            for (size_t i = 0; i != tmp_img.size(); ++i) {
                is >> *dpointer;
                ++dpointer;
            }
        }
        swap(tmp_img);
    }

    void applyMask(const Image<char>& mask) {
        assert(spatial_compatible(mask));
        for (int ch = 0; ch  < numChannels(); ++ch) {
            const char* mp = mask.begin();
            T* ch_end = channelEnd(ch);
            for (T* pt = channelBegin(ch); pt != ch_end; ++pt, ++mp) {
                *pt = (*mp) ? *pt : 0;
            }
        }
    }

    void makeChannelImage(Image<T>& target, int ch) 
    {
        target.resize(dimx(), dimy(), dimz(), 1);
        std::copy(begin() + ch * channelSize(), 
                  begin() + (ch + 1) * channelSize(),
                  target.begin());
    }

protected:
    int dim_[4];

    std::vector<T> data_;
};

}; // end namespace lsseg

#endif // _IMAGE_H

---