MC_grid_3d_scalar.cpp

Go to the documentation of this file.

#include <MC_grid_3d_scalar.hpp>
#include <MC_potential_sphere.hpp>

namespace mesh_conv
{
    MC_grid_3d_scalar::MC_grid_3d_scalar()
            :dimension_internal(MC_v3d(1,1,1)),offset_internal(MC_v3d(0,0,0))
    {}

    MC_grid_3d_scalar::MC_grid_3d_scalar(const MC_int_vector& given_size,const MC_v3d& real_size,const MC_v3d& offset)
            :dimension_internal(real_size),offset_internal(offset)
    {
        potential=MC_grid_3d<float>(given_size);
    }
    const MC_int_vector& MC_grid_3d_scalar::size() const
    {return potential.size();}
    bool MC_grid_3d_scalar::is_within_bound(const MC_v3d& x) const
    {
        MC_v3d y=x-offset_internal;
        if(y[0]>=0 && y[0]<dimension_internal[0])
            if(y[1]>=0 && y[1]<dimension_internal[1])
                if(y[2]>=0 && y[2]<dimension_internal[2])
                    return true;
        return false;
    }
    bool MC_grid_3d_scalar::is_within_bound_voxel(const MC_int_vector& x) const
    {
        MC_int_vector dim=size();
        if(x[0]>=0 && x[0]<dim[0])
            if(x[1]>=0 && x[1]<dim[1])
                if(x[2]>=0 && x[2]<dim[2])
                    return true;
        return false;
    }

    void MC_grid_3d_scalar::fill(const float& value)
    {
        potential.fill(value);
    }
    const MC_v3d& MC_grid_3d_scalar::dimension() const
    {
        return dimension_internal;
    }
    MC_v3d& MC_grid_3d_scalar::dimension()
    {
        return dimension_internal;
    }

    const MC_v3d& MC_grid_3d_scalar::offset() const
    {
        return offset_internal;
    }
    MC_v3d& MC_grid_3d_scalar::offset()
    {
        return offset_internal;
    }

    MC_v3d MC_grid_3d_scalar::convert_real_to_voxel_index(const MC_v3d& x) const
    {
        MC_v3d u;
        for(unsigned int k_dim=0;k_dim<3;++k_dim)
            u[k_dim]=(x[k_dim]-offset_internal[k_dim])/dimension_internal[k_dim]*size()[k_dim];
        return u;
    }
    MC_v3d MC_grid_3d_scalar::convert_voxel_index_to_real(const MC_v3d& u) const
    {
        MC_v3d x;
        for(unsigned int k_dim=0;k_dim<3;++k_dim)
            x[k_dim]=u[k_dim]/static_cast<double>(size()[k_dim])*dimension_internal[k_dim]+offset_internal[k_dim];
        return x;
    }

    const float& MC_grid_3d_scalar::get(const MC_int_vector& u) const
    {
        return potential.get(u);
    }
    float& MC_grid_3d_scalar::get(const MC_int_vector& u)
    {
        return potential.get(u);
    }
    void MC_grid_3d_scalar::set(const MC_int_vector& u,const float& value)
    {
        potential.get(u)=value;
    }
    float MC_grid_3d_scalar::operator()(const MC_v3d& x) const
    {
        return linear_interpolation(x);
    }

    float MC_grid_3d_scalar::linear_interpolation(const MC_v3d& x) const
    {
        if(is_within_bound(x)==false)
            return 0.0;



        MC_v3d u=convert_real_to_voxel_index(x);
        MC_int_vector u0=MC_int_vector(static_cast<int>(u[0]),
                                       static_cast<int>(u[1]),
                                       static_cast<int>(u[2])
                                       );
        MC_v3d alpha=u-MC_v3d(u0[0],u0[1],u0[2]);



        //special case at the boundaries
        if(u0[0]==0 || u0[0]==size()[0]-1 || u0[1]==0 || u0[1]==size()[1]-1 || u0[2]==0 || u0[2]==size()[2]-1)
            return get(u0);

        //linear interpolation
        float x000=get(u0+MC_int_vector(0,0,0));
        float x001=get(u0+MC_int_vector(0,0,1));
        float x010=get(u0+MC_int_vector(0,1,0));
        float x011=get(u0+MC_int_vector(0,1,1));
        float x100=get(u0+MC_int_vector(1,0,0));
        float x101=get(u0+MC_int_vector(1,0,1));
        float x110=get(u0+MC_int_vector(1,1,0));
        float x111=get(u0+MC_int_vector(1,1,1));

        float y=x000*(1-alpha[0])*(1-alpha[1])*(1-alpha[2]) +
                x001*(1-alpha[0])*(1-alpha[1])*(  alpha[2]) +
                x010*(1-alpha[0])*(  alpha[1])*(1-alpha[2]) +
                x011*(1-alpha[0])*(  alpha[1])*(  alpha[2]) +
                x100*(  alpha[0])*(1-alpha[1])*(1-alpha[2]) +
                x101*(  alpha[0])*(1-alpha[1])*(  alpha[2]) +
                x110*(  alpha[0])*(  alpha[1])*(1-alpha[2]) +
                x111*(  alpha[0])*(  alpha[1])*(  alpha[2]) ;
        return y;
    }


    void MC_grid_3d_scalar::add_potential(const implicit::MC_potential_sphere& sphere_implicit)
    {
        const MC_v3d& center=sphere_implicit.sphere().center();
        const double& radius=sphere_implicit.sphere().radius();

        MC_v3d center_voxel=convert_real_to_voxel_index(center);
        MC_int_vector center_voxel_integer=MC_int_vector(static_cast<int>(center_voxel[0]+0.5),static_cast<int>(center_voxel[1]+0.5),static_cast<int>(center_voxel[2]+0.5));


        // ******************************* //
        // TO DO: calculer la fonction de potentiel sur certains voxels et l'ajouter a la grille
        //
        // pour tout les voxels du voisinage
        //   get(voxel) += sphere_implicit(position_du_voxel_dans_l'espace)
        // ******************************* //

    }

}