MC_v3d.cpp

Go to the documentation of this file.

#include <MC_string_helper.hpp>
#include <MC_v3d.hpp>
#include <cmath>

#include <MC_v3d_vector.hpp>
#include <MC_matrix.hpp>
#include <MC_quaternion.hpp>
#include <MC_segment.hpp>


namespace mesh_conv
{

MC_v3d::MC_v3d()
{v[0]=0;v[1]=0;v[2]=0;}

MC_v3d::MC_v3d(const double& x,const double& y,const double& z)
{v[0]=x;v[1]=y;v[2]=z;}

MC_v3d::~MC_v3d(){}

MC_v3d::MC_v3d(const MC_v3d_vector& vec)
{
    if(vec.size()==1)
        *this=vec.first();
    else
    {std::cout<<"Error in MC_v3d::MC_v3d(const MC_v3d_vector& "<<vec<<"), size must be 1"<<std::endl;exit(-1);}
}

const double& MC_v3d::operator()(const int& k_dim) const
{
    if(k_dim<0 || k_dim>2)
    {std::cout<<"Error in const MC_v3d("<<k_dim<<"), size must be 0,1 or 2"<<std::endl;exit(-1);}
    return v[k_dim];
}
double& MC_v3d::operator()(const int& k_dim)
{
    if(k_dim<0 || k_dim>2)
    {std::cout<<"Error in MC_v3d("<<k_dim<<"), size must be 0,1 or 2"<<std::endl;exit(-1);}
    return v[k_dim];
}
const double& MC_v3d::operator[](const int& k_dim) const
{
    if(k_dim<0 || k_dim>2)
    {std::cout<<"Error in const MC_v3d["<<k_dim<<"], size must be 0,1 or 2"<<std::endl;exit(-1);}
    return v[k_dim];
}
double& MC_v3d::operator[](const int& k_dim)
{
    if(k_dim<0 || k_dim>2)
    {std::cout<<"Error in MC_v3d["<<k_dim<<"], size must be 0,1 or 2"<<std::endl;exit(-1);}
    return v[k_dim];
}
MC_v3d::MC_v3d(const MC_v3d& _v)
{v[0]=_v[0];v[1]=_v[1];v[2]=_v[2];}

MC_v3d::MC_v3d(const double* _v)
{v[0]=_v[0];v[1]=_v[1];v[2]=_v[2];}

MC_v3d MC_v3d::operator-() const
{return MC_v3d(-v[0],-v[1],-v[2]);}


MC_v3d operator+(const MC_v3d& vec,const double& to_add)
{return MC_v3d(vec[0]+to_add,vec[1]+to_add,vec[2]+to_add);}
MC_v3d operator-(const MC_v3d& vec,const double& to_sub)
{return MC_v3d(vec[0]-to_sub,vec[1]-to_sub,vec[2]-to_sub);}
MC_v3d operator*(const MC_v3d& vec,const double& to_mult)
{return MC_v3d(vec[0]*to_mult,vec[1]*to_mult,vec[2]*to_mult);}
MC_v3d operator*(const double& to_mult,const MC_v3d& vec)
{return vec*to_mult;}

MC_v3d operator/(const MC_v3d& vec,const double& to_subdiv)
{
    double epsilon=0.000000001;
    if(fabs(to_subdiv)<epsilon)
    {std::cout<<"Error in MC_v3d("<<vec<<"/"<<to_subdiv<<", divide by zero"<<std::endl;exit(-1);}

    return MC_v3d(vec[0]/to_subdiv,vec[1]/to_subdiv,vec[2]/to_subdiv);
}
MC_v3d operator+(const MC_v3d& vec,const MC_v3d& to_add)
{return MC_v3d(vec[0]+to_add[0],vec[1]+to_add[1],vec[2]+to_add[2]);}
MC_v3d operator-(const MC_v3d& vec,const MC_v3d& to_sub)
{return MC_v3d(vec[0]-to_sub[0],vec[1]-to_sub[1],vec[2]-to_sub[2]);}

MC_v3d& MC_v3d::operator+=(const double& to_add)
{
    v[0]+=to_add;v[1]+=to_add;v[2]+=to_add;
    return *this;
}

MC_v3d& MC_v3d::operator-=(const double& to_sub)
{
    v[0]-=to_sub;v[1]-=to_sub;v[2]-=to_sub;
    return *this;
}
MC_v3d& MC_v3d::operator*=(const double& to_mult)
{
    v[0]*=to_mult;v[1]*=to_mult;v[2]*=to_mult;
    return *this;
}
MC_v3d& MC_v3d::operator/=(const double& to_subdiv)
{
    double epsilon=0.000000001;
    if(fabs(to_subdiv)<epsilon)
    {std::cout<<"Error in MC_v3d("<<*this<<"/="<<to_subdiv<<", divide by zero"<<std::endl;exit(-1);}

    v[0]/=to_subdiv;v[1]/=to_subdiv;v[2]/=to_subdiv;
    return *this;
}
MC_v3d& MC_v3d::operator+=(const MC_v3d& to_add){
    v[0]+=to_add[0];v[1]+=to_add[1];v[2]+=to_add[2];
    return *this;
}
MC_v3d& MC_v3d::operator-=(const MC_v3d& to_sub){
    v[0]-=to_sub[0];v[1]-=to_sub[1];v[2]-=to_sub[2];return *this;
}


std::ostream& operator<<(std::ostream& stream,const MC_v3d& v)
{
    stream<<"("<<v[0]<<","<<v[1]<<","<<v[2]<<")";
    return stream;
}

MC_v3d::MC_v3d(const std::string& s)
{
    std::vector <std::string> v_s=mesh_conv::MC_string_tokenizer::tokenize(s," (),;");

    int number_ok=0;
    for(unsigned int k=0;number_ok<3 && k<v_s.size();k++)
    {
        bool is_converted=false;
        double temp_value=MC_string_converter::value_of<double>(v_s[k],&is_converted);
        if(is_converted==true)
            v[number_ok++]=temp_value;
    }

}


    bool operator>>(std::istream& stream,MC_v3d& vec)
    {
        int current_position=0;
        std::string temp_string;
        while(current_position<3 && stream.good())
        {
            if(stream>>temp_string)
            {
                std::vector <std::string> v_token=MC_string_tokenizer::tokenize(temp_string,",;() ");

                for(unsigned int k=0;current_position<3 && k<v_token.size();++k)
                {
                    bool is_converted=false;
                    double temp_value=MC_string_converter::value_of<double>(v_token[k],&is_converted);
                    if(is_converted==true)
                        vec[current_position++]=temp_value;
                }
            }
        }
        return 0;
    }

    MC_v3d& MC_v3d::scale(const double& sx,const double& sy,const double& sz)
    {v[0]*=sx;v[1]*=sy;v[2]*=sz;return *this;}
    MC_v3d MC_v3d::mask(const double& sx,const double& sy,const double& sz) const
    {return MC_v3d(v[0]*sx,v[1]*sy,v[2]*sz);}
    MC_v3d& MC_v3d::scale(const MC_v3d& scaling){v[0]*=scaling[0];v[1]*=scaling[1];v[2]*=scaling[2];return *this;}
    MC_v3d MC_v3d::mask(const MC_v3d& scaling) const{return MC_v3d(v[0]*scaling[0],v[1]*scaling[1],v[2]*scaling[2]);}
    MC_v3d& MC_v3d::scale(const double& scaling){return (*this)*=scaling;}
    MC_v3d MC_v3d::scale(const double& scaling) const{return (*this)*scaling;}

    MC_v3d MC_v3d::scaled(const MC_v3d& scaling) const
    {return MC_v3d(v[0]*scaling[0],v[1]*scaling[1],v[2]*scaling[2]);}
    MC_v3d MC_v3d::scaled(const double& scaling) const
    {return MC_v3d(v[0]*scaling,v[1]*scaling,v[2]*scaling);}


    double MC_v3d::dot(const MC_v3d& vec) const
    {return v[0]*vec[0]+v[1]*vec[1]+v[2]*vec[2];}
    MC_v3d MC_v3d::cross(const MC_v3d& vec) const
    {return MC_v3d(v[1]*vec[2]-v[2]*vec[1],v[2]*vec[0]-v[0]*vec[2],v[0]*vec[1]-v[1]*vec[0]);}
    double MC_v3d::norm() const
    {return sqrt(dot(*this));}
    MC_v3d MC_v3d::normalized() const
    {
        double epsilon=0.000000001;
        double n=norm();
        if(n<epsilon)
        {//std::cout<<"Warning in MC_v3d::normalized(), norm is "<<n<<std::endl;
            return MC_v3d(0,0,1);}
        return (*this)/n;
    }

    std::pair <MC_v3d,MC_v3d_vector> MC_v3d::normalized_with_gradient() const
    {
        MC_v3d_vector nabla=MC_v3d_vector::zeros(3);

        double norm=sqrt(v[0]*v[0]+v[1]*v[1]+v[2]*v[2]);
        double epsilon=0.00000001;
        if(norm<epsilon)
            return std::pair <MC_v3d,MC_v3d_vector> (MC_v3d(0,0,1),nabla);

        MC_v3d n_normalized=MC_v3d(v[0]/norm,v[1]/norm,v[2]/norm);
        nabla[0]=1/norm*(MC_v3d(1,0,0)-n_normalized*v[0]/norm);
        nabla[1]=1/norm*(MC_v3d(0,1,0)-n_normalized*v[1]/norm);
        nabla[2]=1/norm*(MC_v3d(0,0,1)-n_normalized*v[2]/norm);

        return std::pair<MC_v3d,MC_v3d_vector> (n_normalized,nabla);
    }



    MC_v3d& MC_v3d::operator*=(const MC_matrix& M)
    {M.internal_product(this); return *this;}

    MC_v3d::MC_v3d(const MC_matrix& M)
    {
        if( ( (M.size_1()==3 || M.size_1()==4) && M.size_2()==1) ||
            ( (M.size_2()==3 || M.size_2()==4) && M.size_1()==1) )
        {v[0]=M(0);v[1]=M(1);v[2]=M(2);}
    }

    MC_v3d& MC_v3d::operator=(const MC_matrix& M)
    {*this=MC_v3d(M);return *this;}

    MC_v3d operator*(const MC_quaternion& q,const MC_v3d& vec)
    {
        MC_quaternion q_b=q.conjugated();
        MC_quaternion v_quat;
        v_quat.set_v3d(vec);

        MC_quaternion result = q*v_quat*q_b;
        return MC_v3d(result[0],result[1],result[2]);
    }
    MC_v3d& MC_v3d::operator*=(const MC_quaternion& q)
    {
        MC_quaternion q_b=q.conjugated();
        MC_quaternion v_quat(*this,0);
        v_quat = q*v_quat*q_b;
        v[0]=v_quat[0];v[1]=v_quat[1];v[2]=v_quat[2];
        return *this;
    }
    const double* MC_v3d::pointer() const{return v;}

    bool operator==(const MC_v3d& a1,const MC_v3d& a2)
    {
        MC_v3d_less L;
        if(L(a1,a2)==false && L(a2,a1)==false)
            return true;
        return false;
    }
    bool operator!=(const MC_v3d& a1,const MC_v3d& a2)
    {return !(a1==a2);}

    double MC_v3d::area(const MC_v3d& v0,const MC_v3d& v1)
    {
        return (v0.cross(v1)).norm()/2.0;
    }

    double MC_v3d::cotan(const MC_v3d& u0,const MC_v3d& u1)
    {
       double div=(u0.cross(u1)).norm();

       double epsilon=0.00000001;
       if(div<epsilon)
       {std::cout<<"Error in V_3D::cotan(), points are aligned"<<std::endl;return 0.0;}

       return (u0.dot(u1))/(div);
    }
    double* MC_v3d::pointer_unprotected() {return &v[0];}

    double MC_v3d::angle(const MC_v3d& u0,const MC_v3d& u1)
    {
        double epsilon=0.00001;
        double n0=u0.norm(),n1=u1.norm();
        if(n0>epsilon && n1>epsilon)
        {
            double dot=u0.dot(u1)/(n0*n1);
            if(dot>=1.0)
                return 0.0;
            if(dot<=-1.0)
                return M_PI;
            return acos(dot);
        }
        else
            return 0.0;
    }

    MC_v3d MC_v3d::project_on_plane(const MC_v3d& normal) const
    {
        double epsilon=0.00001;
        if(normal.norm()<epsilon)
            return *this;
        MC_v3d n=normal.normalized();
        return (*this)-(*this).dot(n)*n;
    }

    MC_v3d::MC_v3d(const MC_double_vector& vec)
    {
        if(vec.size()!=3)
        {std::cout<<"Error in MC_v3d::MC_v3d(MC_double_vector), size ="<<vec.size()<<"!=3"<<std::endl;exit(-1);}

        for(int k_dim=0;k_dim<3;++k_dim)
            v[k_dim]=vec[k_dim];
    }

    std::string MC_v3d::to_string() const
    {
        std::stringstream stream;
        stream<<*this;
        return stream.str();
    }

    MC_v3d MC_v3d::project_on_line(const MC_v3d& dir) const
    {
        MC_v3d u=dir.normalized();
        return (*this).dot(u)*u;
    }
    MC_v3d MC_v3d::project_on_line(const MC_segment& seg) const
    {
        MC_v3d u=seg.unit_vector();
        return seg[0]+((*this)-seg[0]).dot(u)*u;
    }
    void MC_v3d::set_zero()
    {v[0]=0;v[1]=0;v[2]=0;}
}