#ifndef SOFA_BRANCHINGIMAGE_IMAGESAMPLER_H
#define SOFA_BRANCHINGIMAGE_IMAGESAMPLER_H

#include <BranchingImage/config.h>
#include <image/ImageSampler.h>
#include "ImageAlgorithms.h"
#include "BranchingImage.h"


namespace sofa
{

namespace component
{

namespace engine
{


/// Specialization for BranchingImage
template <class T>
struct ImageSamplerSpecialization<defaulttype::BranchingImage<T>>
{
    typedef ImageSampler<defaulttype::BranchingImage<T>> ImageSamplerT;

    typedef defaulttype::BranchingImage<SReal> DistTypes;
    typedef defaulttype::BranchingImage<unsigned int> VorTypes;

    static void init( ImageSamplerT* sampler )
    {
        sampler->addAlias( &sampler->image, "branchingImage" );
    }

    static void regularSampling( ImageSamplerT* sampler, const bool atcorners=false, const bool recursive=false )
    {
        if( !atcorners )
        {
            sampler->serr<<"ImageSamplerT::regularSampling - only at corner is implemented\n";
        }


        typedef typename ImageSamplerT::Real Real;
        typedef typename ImageSamplerT::Coord Coord;
        typedef typename ImageSamplerT::Edge Edge;
        typedef typename ImageSamplerT::Hexa Hexa;
//        typedef typename ImageSamplerT::T T;


        // get tranform and image at time t
        typename ImageSamplerT::raImage in(sampler->image);
        typename ImageSamplerT::raTransform inT(sampler->transform);
        const typename ImageSamplerT::ImageTypes::BranchingImage3D& inimg = in->imgList[sampler->time];

        // data access
        typename ImageSamplerT::waPositions pos(sampler->position);       pos.clear();
        typename ImageSamplerT::waEdges e(sampler->edges);                e.clear();
        typename ImageSamplerT::waEdges g(sampler->graphEdges);           g.clear();
        typename ImageSamplerT::waHexa h(sampler->hexahedra);             h.clear();

        const typename ImageSamplerT::ImageTypes::Dimension& dim = in->getDimension();



        unsigned index1d = 0;

        {
        std::map< unsigned, std::map<unsigned, unsigned> > hindices; // for each superimposed voxel (index1d,offset) -> hexa index

        {
        // add hexa with independant vertices
        unsigned indexVertex = 0;
        for( unsigned z=0 ; z<dim[ImageSamplerT::ImageTypes::DIMENSION_Z] ; ++z )
        for( unsigned y=0 ; y<dim[ImageSamplerT::ImageTypes::DIMENSION_Y] ; ++y )
        for( unsigned x=0 ; x<dim[ImageSamplerT::ImageTypes::DIMENSION_X] ; ++x )
        {
            for( unsigned v=0 ; v<inimg[index1d].size() ; ++v )
            {
                h.push_back( Hexa( indexVertex, indexVertex+1, indexVertex+2, indexVertex+3, indexVertex+4, indexVertex+5, indexVertex+6, indexVertex+7 ) );
                indexVertex += 8;
                hindices[index1d][v] = h.size()-1;
            }
            ++index1d;
        }
        }

        // link vertices
        index1d = 0;
        for( unsigned z=0 ; z<dim[ImageSamplerT::ImageTypes::DIMENSION_Z] ; ++z )
        for( unsigned y=0 ; y<dim[ImageSamplerT::ImageTypes::DIMENSION_Y] ; ++y )
        for( unsigned x=0 ; x<dim[ImageSamplerT::ImageTypes::DIMENSION_X] ; ++x )
        {
            for( unsigned v=0 ; v<inimg[index1d].size() ; ++v )
            {
                Hexa& hexa = h[hindices[index1d][v]];

                for( unsigned n=0 ; n<inimg[index1d][v].neighbours.size() ; ++n )
                {
                    const unsigned neighbourIndex = inimg[index1d][v].neighbours[n].index1d;

                    const unsigned neighbourOffset = inimg[index1d][v].neighbours[n].offset;

                    Hexa& neighbor = h[hindices[neighbourIndex][neighbourOffset]];

                    typename ImageSamplerT::ImageTypes::NeighbourOffset dir = in->getDirection( index1d, neighbourIndex );

                    switch( dir.connectionType() )
                    {
                        case ImageSamplerT::ImageTypes::NeighbourOffset::FACE:
                        {
                            if( neighbourIndex < index1d ) // merge only one way (ensure neighbour is only left or bottom or back) -> not enough for 26-connectivity
                            {
                                if( dir[0] ) //LEFT
                                {
                                    assert( dir[0]==-1 );
                                    mergeVertexIndex( h, hexa[0], neighbor[1] );
                                    mergeVertexIndex( h, hexa[4], neighbor[5] );
                                    mergeVertexIndex( h, hexa[7], neighbor[6] );
                                    mergeVertexIndex( h, hexa[3], neighbor[2] );
                                }
                                else if( dir[1] ) // BOTTOM
                                {
                                    assert( dir[1]==-1 );
                                    mergeVertexIndex( h, hexa[0], neighbor[3] );
                                    mergeVertexIndex( h, hexa[1], neighbor[2] );
                                    mergeVertexIndex( h, hexa[4], neighbor[7] );
                                    mergeVertexIndex( h, hexa[5], neighbor[6] );
                                }
                                else // BACK
                                {
                                    assert( dir[2]==-1 );
                                    mergeVertexIndex( h, hexa[0], neighbor[4] );
                                    mergeVertexIndex( h, hexa[1], neighbor[5] );
                                    mergeVertexIndex( h, hexa[2], neighbor[6] );
                                    mergeVertexIndex( h, hexa[3], neighbor[7] );
                                }
                            }
                            break;
                        }
                        case ImageSamplerT::ImageTypes::NeighbourOffset::EDGE: // 26-connectivity
                        {
                            // test only 9 on 12 (they will be treated by the neighbour)...
                            if( dir[0]==-1 ) // LEFT
                            {
                                if( dir[1]==-1 ) // BOTTOM
                                {
                                    assert( dir[2]==0 );
                                    mergeVertexIndex( h, hexa[0], neighbor[2] );
                                    mergeVertexIndex( h, hexa[4], neighbor[6] );
                                }
                                else if( dir[1]==1 ) // TOP
                                {
                                    assert( dir[2]==0 );
                                    mergeVertexIndex( h, hexa[3], neighbor[1] );
                                    mergeVertexIndex( h, hexa[7], neighbor[5] );
                                }
                                else if( dir[2]==-1 ) // BACK
                                {
                                    assert( dir[1]==0 );
                                    mergeVertexIndex( h, hexa[0], neighbor[5] );
                                    mergeVertexIndex( h, hexa[3], neighbor[6] );
                                }
                                else // FRONT
                                {
                                    assert( dir[2]==1 );
                                    assert( dir[1]==0 );
                                    mergeVertexIndex( h, hexa[4], neighbor[1] );
                                    mergeVertexIndex( h, hexa[7], neighbor[2] );
                                }
                            }
                            else if( dir[0]==1 ) // RIGHT
                            {
                                if( dir[1]==-1 ) // BOTTOM
                                {
                                    assert( dir[2]==0 );
                                    mergeVertexIndex( h, hexa[1], neighbor[3] );
                                    mergeVertexIndex( h, hexa[5], neighbor[7] );
                                }
                                else if( dir[1]==1 ) // TOP
                                {
                                    assert( dir[2]==0 );
//                                    mergeVertexIndex( h, hexa[2], neighbor[0] );
//                                    mergeVertexIndex( h, hexa[6], neighbor[4] );
                                }
                                else if( dir[2]==-1 ) // BACK
                                {
                                    assert( dir[1]==0 );
                                    mergeVertexIndex( h, hexa[1], neighbor[4] );
                                    mergeVertexIndex( h, hexa[2], neighbor[7] );
                                }
                                else // FRONT
                                {
                                    assert( dir[2]==1 );
                                    assert( dir[1]==0 );
//                                    mergeVertexIndex( h, hexa[5], neighbor[0] );
//                                    mergeVertexIndex( h, hexa[6], neighbor[3] );
                                }
                            }
                            else // CENTER
                            {
                                assert( dir[0]==0 );
                                assert( dir[1]!=0 );
                                assert( dir[2]!=0 );

                                if( dir[1]==-1 ) // BOTTOM
                                {
                                    if( dir[2]==-1 ) // BACK
                                    {
                                        mergeVertexIndex( h, hexa[0], neighbor[7] );
                                        mergeVertexIndex( h, hexa[1], neighbor[6] );
                                    }
                                    else // FRONT
                                    {
                                        mergeVertexIndex( h, hexa[4], neighbor[3] );
                                        mergeVertexIndex( h, hexa[5], neighbor[2] );
                                    }
                                }
                                else // TOP
                                {
                                    if( dir[2]==-1 ) // BACK
                                    {
                                        mergeVertexIndex( h, hexa[2], neighbor[5] );
                                        mergeVertexIndex( h, hexa[3], neighbor[4] );
                                    }
                                    else // FRONT
                                    {
//                                        mergeVertexIndex( h, hexa[6], neighbor[1] );
//                                        mergeVertexIndex( h, hexa[7], neighbor[0] );
                                    }
                                }
                            }
                            break;
                        }
                        case ImageSamplerT::ImageTypes::NeighbourOffset::CORNER: // 26-connectivity
                        {
                            // test only 4 on 8
                            assert( abs(dir[0])==1 && abs(dir[1])==1 && abs(dir[2])==1 );
                            if( dir[0]==-1 ) // LEFT
                            {
                                if( dir[1]==-1 ) // BOTTOM
                                {
                                    if( dir[2]==-1 ) // BACK
                                    {
                                        mergeVertexIndex( h, hexa[0], neighbor[6] );
                                    }
                                    else // FRONT
                                    {
                                        mergeVertexIndex( h, hexa[4], neighbor[2] );
                                    }
                                }
                                else // TOP
                                {
                                    if( dir[2]==-1 ) // BACK
                                    {
                                        mergeVertexIndex( h, hexa[3], neighbor[5] );
                                    }
                                    else // FRONT
                                    {
//                                        mergeVertexIndex( h, hexa[7], neighbor[1] );
                                    }
                                }
                            }
                            else // RIGHT
                            {
                                if( dir[1]==-1 ) // BOTTOM
                                {
                                    if( dir[2]==-1 ) // BACK
                                    {
                                        mergeVertexIndex( h, hexa[1], neighbor[7] );
                                    }
                                    else // FRONT
                                    {
//                                        mergeVertexIndex( h, hexa[5], neighbor[3] );
                                    }
                                }
                                else // TOP
                                {
                                    if( dir[2]==-1 ) // BACK
                                    {
//                                        mergeVertexIndex( h, hexa[2], neighbor[4] );
                                    }
                                    else // FRONT
                                    {
//                                        mergeVertexIndex( h, hexa[6], neighbor[0] );
                                    }
                                }
                            }
                            break;
                        }
                        case ImageSamplerT::ImageTypes::NeighbourOffset::ONPLACE: // 7- or 27-connectivity
                        {
                            for( unsigned w=0 ; w<8 ; ++w )
                                mergeVertexIndex( h, hexa[w], neighbor[w] );
                            break;
                        }
                        case ImageSamplerT::ImageTypes::NeighbourOffset::NOTCLOSE: // a not connected, far neighbour
                        default:
                        {
                            break;
                        }
                    }
                }
            }
            index1d++;
        }
        }

        {
            // give a continue index from 0 to max (without hole)
            unsigned continueIndex = 0;
            std::map<unsigned,unsigned> continueMap;

            for( unsigned i=0 ; i<h.size() ; ++i )
            for( unsigned j=0 ; j<8 ; ++j )
            {
                if( continueMap.find(h[i][j])==continueMap.end() ) continueMap[h[i][j]]=continueIndex++;
            }
            for( unsigned i=0 ; i<h.size() ; ++i )
            for( unsigned j=0 ; j<8 ; ++j )
            {
                h[i][j] = continueMap[h[i][j]];
            }
            pos.resize( continueIndex );
        }

        {
            std::map<unsigned,bool> alreadyAddedPos;
            std::map<Edge,bool> alreadyAddedEdge;
            index1d = 0;
            unsigned indexHexa = 0;
            for( unsigned z=0 ; z<dim[ImageSamplerT::ImageTypes::DIMENSION_Z] ; ++z )
            for( unsigned y=0 ; y<dim[ImageSamplerT::ImageTypes::DIMENSION_Y] ; ++y )
            for( unsigned x=0 ; x<dim[ImageSamplerT::ImageTypes::DIMENSION_X] ; ++x )
            {
                for( unsigned v=0 ; v<inimg[index1d].size() ; ++v )
                {
                    static const Real gap = (Real)-0.5;
                    static const Real hexaCornerGapFromCenter[8][3] = { {gap,gap,gap},{1+gap,gap,gap},{1+gap,1+gap,gap},{gap,1+gap,gap},{gap,gap,1+gap},{1+gap,gap,1+gap},{1+gap,1+gap,1+gap},{gap,1+gap,1+gap}  };

                    Hexa& hexa = h[indexHexa++];
                    for( unsigned j=0 ; j<8 ; ++j )
                    {
                        if( alreadyAddedPos.find(hexa[j])==alreadyAddedPos.end() )
                        {
                            alreadyAddedPos[hexa[j]] = true;
                            pos[hexa[j]] = Coord(x+hexaCornerGapFromCenter[j][0],y+hexaCornerGapFromCenter[j][1],z+hexaCornerGapFromCenter[j][2]);
                        }
                    }

                    Edge edge( hexa[0], hexa[1] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[1], hexa[2] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[2], hexa[3] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[3], hexa[0] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[4], hexa[5] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[5], hexa[6] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[6], hexa[7] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[7], hexa[4] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[0], hexa[4] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[1], hexa[5] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[2], hexa[6] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                    edge = Edge( hexa[3], hexa[7] );
                    if( alreadyAddedEdge.find( edge )==alreadyAddedEdge.end() )
                    {
                        alreadyAddedEdge[edge] = true;
                        alreadyAddedEdge[Edge(edge[1],edge[0])] = true;
                        e.push_back( edge );
                    }
                }

                ++index1d;
            }
        }


        if(recursive)
        {
            helper::vector<unsigned int> indices; indices.resize(pos.size()); for(unsigned int i=0; i<pos.size(); i++) indices[i]=i;
            sampler->subdivide(indices);
        }

        for(unsigned int i=0; i<pos.size(); i++) pos[i]=inT->fromImage(pos[i]);
    }

    template<class Hexas>
    static void mergeVertexIndex( Hexas& h, unsigned index0, unsigned index1 )
    {
        if( index0==index1 ) return;
        for( unsigned i=0 ; i<h.size() ; ++i )
        for( unsigned j=0 ; j<8 ; ++j )
        {
            if( h[i][j]==index1 ) h[i][j]=index0;
        }
    }

    static void uniformSampling( ImageSamplerT* sampler,const unsigned int nb=0,  const bool bias=false, const unsigned int lloydIt=100, const unsigned int method=FASTMARCHING, const unsigned int /*pmmIter*/=std::numeric_limits<unsigned int>::max(), const SReal /*pmmTol*/=10)
    {
        typedef typename ImageSamplerT::Real Real;
        typedef typename ImageSamplerT::Coord Coord;
//        typedef typename ImageSamplerT::Edge Edge;
//        typedef typename ImageSamplerT::Hexa Hexa;

        clock_t timer = clock();

        // get tranform and image at time t
        typename ImageSamplerT::raImage in(sampler->image);
        typename ImageSamplerT::raTransform inT(sampler->transform);
        const typename ImageSamplerT::ImageTypes::BranchingImage3D& inimg = in->imgList[sampler->time];
        const typename ImageSamplerT::ImageTypes* biasFactor=bias?&in.ref():NULL;

        // data access
        typename ImageSamplerT::raPositions fpos(sampler->fixedPosition);
        typename ImageSamplerT::waEdges e(sampler->edges);                e.clear();
        typename ImageSamplerT::waEdges g(sampler->graphEdges);           g.clear();
        typename ImageSamplerT::waHexa h(sampler->hexahedra);             h.clear();

        typename ImageSamplerT::ImageTypes::Dimension dim = in->getDimension();

        // init voronoi and distances
        dim[ImageSamplerT::ImageTypes::DIMENSION_S]=1;
        dim[ImageSamplerT::ImageTypes::DIMENSION_T]=1;

        typename ImageSamplerT::waVor waVor(sampler->voronoi);
        VorTypes& voronoi = waVor.wref();
        voronoi.setDimensions(dim);
        voronoi.imgList[0].cloneTopology<T> (inimg,1,0);

        typename ImageSamplerT::waDist wadist(sampler->distances);
        DistTypes& dist = wadist.wref();
        dist.setDimensions(dim);
        dist.imgList[0].cloneTopology<T> (inimg,1,-1.0);

        bimg_forCVoffT(in.ref(),c,v,off1D,t) if(t==sampler->time) if(in.ref()(off1D,v,c,t)) dist(off1D,v,c,0)=cimg_library::cimg::type<Real>::max();

        // list of seed points
        typedef typename ImageSamplerT::ImageTypes::VoxelIndex VoxelIndex;
        std::set<std::pair<Real, VoxelIndex > > trial;

        // add fixed points
        helper::vector<unsigned int> fpos_voronoiIndex;
        helper::vector<VoxelIndex> fpos_VoxelIndex;

        for(size_t i=0; i<fpos.size(); i++)
        {
            fpos_voronoiIndex.push_back(i+1);
            Coord p = inT->toImageInt(fpos[i]);
            VoxelIndex ind (dist.index3Dto1D(p[0],p[1],p[2]), 0); // take first superimposed voxel    TO DO: identify it from fine resolution
            fpos_VoxelIndex.push_back(ind);
            AddSeedPoint<Real>(trial,dist,voronoi, fpos_VoxelIndex[i],fpos_voronoiIndex[i]);
        }
        if(fpos.size())
        {
            switch(method)
            {
            case FASTMARCHING : fastMarching<Real,T>(trial,dist, voronoi, sampler->transform.getValue().getScale(),biasFactor ); break;
            case DIJKSTRA : dijkstra<Real,T>(trial,dist, voronoi, sampler->transform.getValue().getScale(), biasFactor); break;
            default : std::cerr << "Unknown Distance Field Computation Method" << std::endl; break;
            };
        }

        // farthest point sampling using geodesic distances
        helper::vector<unsigned int> pos_voronoiIndex;
        helper::vector<VoxelIndex> pos_VoxelIndex;
        while(pos_VoxelIndex.size()<nb)
        {
            Real dmax=0;  VoxelIndex indMax;
            bimg_forCVoffT(dist,c,v,off1D,t) if(dist(off1D,v,c,t)>dmax) { dmax=dist(off1D,v,c,t); indMax = VoxelIndex(off1D,v); }
            if(dmax)
            {
                pos_voronoiIndex.push_back(fpos_VoxelIndex.size()+pos_VoxelIndex.size()+1);
                pos_VoxelIndex.push_back(indMax);
                AddSeedPoint<Real>(trial,dist,voronoi, pos_VoxelIndex.back(),pos_voronoiIndex.back());
                switch(method)
                {
                case FASTMARCHING : fastMarching<Real,T>(trial,dist, voronoi, sampler->transform.getValue().getScale(),biasFactor ); break;
                case DIJKSTRA : dijkstra<Real,T>(trial,dist, voronoi, sampler->transform.getValue().getScale(), biasFactor); break;
                default : std::cerr << "Unknown Distance Field Computation Method" << std::endl; break;
                };
            }
            else break;
        }

        unsigned int it=0;
        bool converged =(it>=lloydIt)?true:false;

        while(!converged)
        {
            if(Lloyd<Real>(pos_VoxelIndex,pos_voronoiIndex,voronoi)) // one lloyd iteration
            {
                // recompute distance from scratch
                bimg_forCVoffT(dist,c,v,off1D,t) if(dist(off1D,v,c,t)!=-1) dist(off1D,v,c,t)=cimg_library::cimg::type<Real>::max();
                for(unsigned int i=0; i<fpos_voronoiIndex.size(); i++) AddSeedPoint<Real>(trial,dist,voronoi, fpos_VoxelIndex[i], fpos_voronoiIndex[i]);
                for(unsigned int i=0; i<pos_voronoiIndex.size(); i++) AddSeedPoint<Real>(trial,dist,voronoi, pos_VoxelIndex[i], pos_voronoiIndex[i]);
                switch(method)
                {
                case FASTMARCHING : fastMarching<Real,T>(trial,dist, voronoi, sampler->transform.getValue().getScale(),biasFactor ); break;
                case DIJKSTRA : dijkstra<Real,T>(trial,dist, voronoi, sampler->transform.getValue().getScale(), biasFactor); break;
                default : std::cerr << "Unknown Distance Field Computation Method" << std::endl; break;
                };
                it++; if(it>=lloydIt) converged=true;
            }
            else converged=true;
        }

        // add 3D points
        std::vector<defaulttype::Vec<3,Real> >& pos = *sampler->position.beginEdit();    pos.clear();
        for(size_t i=0; i<pos_VoxelIndex.size(); i++)
        {
            unsigned x,y,z; dist.index1Dto3D(pos_VoxelIndex[i].index1d,x,y,z);
            pos.push_back(inT->fromImage(Coord(x,y,z)));
        }
        sampler->position.endEdit();


        if(sampler->f_printLog.getValue())
        {
            std::cout<<sampler->getName()<<": sampling completed in "<< it <<" Lloyd iterations ("<< (clock() - timer) / (float)CLOCKS_PER_SEC <<"s )"<<std::endl;
        }

    }

    static void recursiveUniformSampling( ImageSamplerT* sampler,const unsigned int /*nb*/=0,  const bool /*bias*/=false, const unsigned int /*lloydIt*/=100,const unsigned int /*method*/=FASTMARCHING,  const unsigned int /*N*/=1, const unsigned int /*pmmIter*/=std::numeric_limits<unsigned int>::max(), const SReal /*pmmTol*/=10 )
    {
        sampler->serr<<"ImageSamplerT::recursiveUniformSampling is not yet immplemented for BranchingImage\n";
    }
};






} // namespace engine

} // namespace component

} // namespace sofa

#endif // SOFA_BRANCHINGIMAGE_IMAGESAMPLER_H
