3d/vector.hh

Go to the documentation of this file.

/* -*- mia-c++  -*-
 *
 * This file is part of MIA - a toolbox for medical image analysis 
 * Copyright (c) Leipzig, Madrid 1999-2017 Gert Wollny
 *
 * MIA 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; either version 3 of the License, or
 * (at your option) any later version.
 *
 * 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 MIA; if not, see <http://www.gnu.org/licenses/>.
 *
 */

#ifndef __MIA_3DVECTOR_HH
#define __MIA_3DVECTOR_HH 1

#include <typeinfo>
#include <assert.h>
#include <stdexcept>
#include <math.h>
#include <complex>
#include <iostream>
#include <type_traits>

#include <mia/core/defines.hh>
#include <mia/core/type_traits.hh>
#include <mia/core/attributetype.hh>
#include <mia/3d/defines3d.hh>

NS_MIA_BEGIN

template < class T > 
class T3DVector {
public:
        T x;
        T y;
        T z;
        
        typedef T value_type; 

        T3DVector():x(T()),y(T()),z(T()){};
        
        explicit T3DVector(int dim):x(T()),y(T()),z(T()){
                assert(dim == 3);
        }
        
        T3DVector(const T3DVector<T>& other) = default; 
        
        T3DVector<T>& operator = (const T3DVector<T>& other) = default; 

        T3DVector(const T& x_,const T& y_,const T& z_):
                x(x_),y(y_),z(z_){
        }
        
        template <class in> explicit T3DVector(const T3DVector<in>& org):
                x(T(org.x)),y(T(org.y)),z(T(org.z)){
        }
        
        template <class in> 
        T3DVector<T>& operator =(const T3DVector<in>& org){
                x=org.x; y=org.y; z=org.z;
                return *this;
        }
        
        double norm2()const{
                return x * x + y * y + z * z;
        }
        
        double product() const {
                return x * y * z; 
        }
        double norm()const{
                return sqrt(norm2());
        }
        
        int size() const {
                return 3;
        }

        void fill(T v) {
                x = y = z = v;
        }

        const T operator [](size_t i) const {
                assert(i < 3); 
                switch (i) {
                case 0:return x; 
                case 1:return y; 
                case 2:return z; 
                default: 
                        throw std::logic_error("Access to vectorelement out of range"); 
                }
        }

        T& operator [](size_t i) {
                assert(i < 3); 
                switch (i) {
                case 0:return x; 
                case 1:return y; 
                case 2:return z; 
                default: 
                        throw std::logic_error("Access to vectorelement out of range"); 
                }
        }
        
        T3DVector<T>& operator +=(const T3DVector<T>& a){
                x+=a.x; y+=a.y; z+=a.z;
                return *this;
        }
        
        T3DVector<T>& operator -=(const T3DVector<T>& a){
                x-=a.x; y-=a.y; z-=a.z;
                return *this;
        }
        
        T3DVector<T>& operator *=(const double a){
                x = T(x * a); y = T(y * a); z = T(z * a);
                return *this;
        }

        T3DVector<T>& operator *=(const T3DVector<T>& a){
                x = T(x * a.x); y = T(y * a.y); z = T(z * a.z);
                return *this;
        }

        
        T3DVector<T>& operator /=(const double a){
                assert(a != 0.0);
                x = T(x/ a); y =T (y / a); z = T(z / a);
                return *this;
        }

        T3DVector operator -() const {
                return  T3DVector<T>(-x, -y, -z);
        }

        void write(std::ostream& os)const {
                os  << x << "," << y << "," << z; 
        }

        void read(std::istream& is) {
                char c; 
                
                T r,s,t; 
                is >> c;
                // if we get the opening delimiter '<' then we also expect the closing '>'
                // otherwise just read three coma separated values. 
                // could use the BOOST lexicel cast for better error handling
                if (c == '<') {
                        is >> r;
                        is >> c; 
                        if (c != ',') {
                                is.clear(std::ios::badbit);
                                return; 
                        }
                        is >> s; 
                        is >> c; 
                        if (c != ',') {
                                is.clear(std::ios::badbit);
                                return; 
                        }
                        is >> t; 
                        is >> c; 
                        if (c != '>') {
                                is.clear(std::ios::badbit);
                                return; 
                        }
                        x = r; 
                        y = s; 
                        z = t;
                }else{
                        is.putback(c); 
                        is >> r;
                        is >> c; 
                        if (c != ',') {
                                is.clear(std::ios::badbit);
                                return; 
                        }
                        is >> s; 
                        is >> c; 
                        if (c != ',') {
                                is.clear(std::ios::badbit);
                                return; 
                        }
                        is >> t; 
                        x = r; 
                        y = s; 
                        z = t;
                }
        }

        const T3DVector<T>&  xyz()const {
                return *this; 
        }
        
        const T3DVector<T> xzy()const {
                return T3DVector<T>(x,z,y); 
        }

        const T3DVector<T> yxz()const {
                return T3DVector<T>(y,x,z); 
        }

        const T3DVector<T> yzx()const {
                return T3DVector<T>(y,z,x); 
        }

        const T3DVector<T> zyx()const {
                return T3DVector<T>(z,y,x); 
        }
        
        const T3DVector<T> zxy()const {
                return T3DVector<T>(z,x,y); 
        }
        
        static T3DVector<T> _1; 

        static T3DVector<T> _0; 
        
        static const unsigned int  elements; 
};


struct EAttributeType_3d : public EAttributeType {
        
        static const int vector_3d_bit = 0x40000; 
        
        static bool is_vector3d(int type) {
                return type & vector_3d_bit; 
        }
}; 

template <typename T> 
struct attribute_type<T3DVector<T>> : public EAttributeType_3d {
        static const int value = attribute_type<T>::value | vector_3d_bit;
}; 


template <typename T> 
struct atomic_data<T3DVector<T> > {
        typedef T type; 
        static const int size; 
}; 

template <typename T> 
const int atomic_data<T3DVector<T> >::size = 3; 


template <typename T>
T3DVector<T> cross(const T3DVector<T>& a, const T3DVector<T>& b) 
{
        return T3DVector<T>( 
                a.y * b.z - b.y * a.z, 
                a.z * b.x - b.z * a.x, 
                a.x * b.y - b.x * a.y
                ); 
}


template <class T> double fabs(const T3DVector<T>& t)
{
        return t.norm();
}

template <class T> double dot(const T3DVector<T>& a, const T3DVector<T>& b)
{
        return a.x * b.x + a. y * b.y + a.z * b.z;
}


typedef T3DVector<float>    C3DFVector;

typedef T3DVector<double>   C3DDVector;

typedef T3DVector<unsigned int>   C3DBounds;



 template <class T> 
std::ostream& operator << (std::ostream& os, const T3DVector<T>& v)
{
        v.write(os);
        return os; 
}

template <class T> 
std::istream& operator >> (std::istream& is, T3DVector<T>& v)
{
        v.read(is);
        return is; 
}


template <class T> 
inline const T3DVector<T> operator +(const T3DVector<T>& a,const T3DVector<T>& b){
        T3DVector<T> tmp(a);
        tmp += b; 
        return tmp;
}

template <typename  T, typename  S>
T3DVector<T> operator +(const T3DVector<T>& a, const T3DVector<S>& b)
{
        return T3DVector<T>(a.x + b.x, a.y + b.y, a.z + b.z); 
}

                
template <class T> 
inline const T3DVector<T> operator -(const T3DVector<T>& a,const T3DVector<T>& b){
        T3DVector<T> tmp(a);
        tmp -= b; 
        return tmp;
}

template <class T> 
inline const T3DVector<T> operator *(const T3DVector<T>& a,const T3DVector<T>& b)
{
        return T3DVector<T>(b.x * a.x, b.y * a.y, b.z * a.z);
}

template <class T>
inline const T3DVector<T> operator /(const T3DVector<T>& a,double f)
{
        assert(f != T());
        T3DVector<T> tmp (a);
        tmp /= f; 
        return tmp;
}

template <class T>
inline const T3DVector<T> operator / (const T3DVector<T>& a, const T3DVector<T>& b)
{
        assert(b.x != 0.0 && b.y != 0.0 && b.z != 0.0);
        return T3DVector<T>(a.x/b.x, a.y/b.y, a.z/b.z); 
}


template <class T> 
inline const T3DVector<T> operator *(const T3DVector<T>& a, double f)
{
        T3DVector<T> tmp (a);
        tmp *= f;
        return tmp;
}


template <class T> 
inline const T3DVector<T> operator *(double f, const T3DVector<T>& a)
{
        return a * f; 
}


template <class T> 
inline const T3DVector<T> operator ^(const T3DVector<T>& a,const T3DVector<T>& b)
{
        return T3DVector<T>( a.y * b.z - b.y * a.z,
                             a.z * b.x - b.z * a.x,
                             a.x * b.y - b.x * a.y);
}

template <class T> 
inline bool operator == (const T3DVector<T>& a,const T3DVector<T>& b)
{
        return (b.x == a.x && b.y == a.y && b.z == a.z);
}

template <class T>
inline bool operator != (const T3DVector<T>& a,const T3DVector<T>& b)
{
        return ! (a == b);
}

template <class T>
bool operator <  (const T3DVector<T>& a,const T3DVector<T>& b){
        return (a.x < b.x && a.y < b.y && a.z < b.z);
}

template <class T>
bool operator <= (const T3DVector<T>& b, const T3DVector<T>& a){
        return (b.x <= a.x && b.y <= a.y && b.z <= a.z);
}

template <class T>
bool operator >  (const T3DVector<T>& b, const T3DVector<T>& a){
        return (b.x > a.x && b.y > a.y && b.z > a.z);
}

template <class T>
bool operator >= (const T3DVector<T>& b, const T3DVector<T>& a){
        return (b.x >= a.x && b.y >= a.y && b.z >= a.z);
}
template <typename T >
T3DVector<T> T3DVector<T>::_1 = T3DVector<T>(1,1,1);

template <typename T >
T3DVector<T> T3DVector<T>::_0 = T3DVector<T>(0,0,0);

template <typename T>
struct less_then<T3DVector<T> > {
        bool operator() (const T3DVector<T>& a, const T3DVector<T>& b) const{
                return a.z < b.z || 
                        (a.z ==  b.z && 
                         (a.y < b.y || (a.y == b.y && a.x < b.x))); 
        }
}; 



NS_MIA_END


#endif