void QuadSetTopologyModifier::removeEdgesProcess( const sofa::helper::vector<unsigned int> &indices,
        const bool removeIsolatedItems)
{
    // Note: this does not check if an edge is removed from an existing quad (it should never happen)

    if(m_container->hasEdgesInQuad()) // this method should only be called when edges exist
    {
        if(!m_container->hasQuadsAroundEdge())
            m_container->createQuadsAroundEdgeArray();

        unsigned int lastEdge = m_container->getNumberOfEdges() - 1;
        for(unsigned int i = 0; i < indices.size(); ++i, --lastEdge)
        {
            // updating the quads connected to the edge replacing the removed one:
            // for all quads connected to the last point
            for(sofa::helper::vector<unsigned int>::iterator itt = m_container->m_quadsAroundEdge[lastEdge].begin();
                itt != m_container->m_quadsAroundEdge[lastEdge].end(); ++itt)
            {
                unsigned int edgeIndex = m_container->getEdgeIndexInQuad(m_container->m_edgesInQuad[*itt], lastEdge);
                m_container->m_edgesInQuad[*itt][edgeIndex] = indices[i];
            }

            // updating the edge shell itself (change the old index for the new one)
            m_container->m_quadsAroundEdge[ indices[i] ] = m_container->m_quadsAroundEdge[ lastEdge ];
        }

        m_container->m_quadsAroundEdge.resize( m_container->getNumberOfEdges() - indices.size() );
    }

    // call the parent's method.
    EdgeSetTopologyModifier::removeEdgesProcess(indices, removeIsolatedItems);
}

void ManifoldTriangleSetTopologyModifier::updateRemovingModifications (const sofa::helper::vector< unsigned int >& edgeToBeRemoved,
        const sofa::helper::vector< unsigned int >& vertexToBeRemoved)
{
    //???
    for (unsigned int i = 0; i <vertexToBeRemoved.size(); i++)
    {
        it_modif = m_modifications.find( vertexToBeRemoved[i] );

        if(it_modif != m_modifications.end())
            m_modifications.erase( vertexToBeRemoved[i] );
    }


    for (unsigned int i = 0; i <edgeToBeRemoved.size(); i++)
    {
        for (unsigned int j = 0; j<m_modificationsEdge.size(); j++)
        {
            if (m_modificationsEdge[j] == edgeToBeRemoved[i])
            {
                m_modificationsEdge.erase(m_modificationsEdge.begin()+j);
                break;
            }
        }
    }

}

void QuadSetTopologyModifier::removeQuadsWarning( sofa::helper::vector<unsigned int> &quads)
{
    m_container->setQuadTopologyToDirty();
    /// sort vertices to remove in a descendent order
    std::sort( quads.begin(), quads.end(), std::greater<unsigned int>() );

    // Warning that these quads will be deleted
    QuadsRemoved *e=new QuadsRemoved(quads);
    addTopologyChange(e);
}

void HexahedronSetTopologyModifier::removeHexahedraWarning( sofa::helper::vector<unsigned int> &hexahedra )
{
    m_container->setHexahedronTopologyToDirty();
    /// sort vertices to remove in a descendent order
    std::sort( hexahedra.begin(), hexahedra.end(), std::greater<unsigned int>() );

    // Warning that these edges will be deleted
    HexahedraRemoved *e = new HexahedraRemoved(hexahedra);
    addTopologyChange(e);
}

void QuadSetTopologyModifier::addQuadsProcess(const sofa::helper::vector< Quad > &quads)
{
    helper::WriteAccessor< Data< sofa::helper::vector<Quad> > > m_quad = m_container->d_quad;
    m_quad.reserve(m_quad.size() + quads.size());

    for(unsigned int i=0; i<quads.size(); ++i)
    {
        addQuadProcess(quads[i]);
    }
}

/// This function is named this way because someone give the getDataFields to the one
/// that returns a dictionary of (name, value) which is not coherente with the c++
/// name of the function.
/// If you are a brave hacker, courageous enough to break backward compatibility you can
/// probably fix all this mess.
static PyObject * Base_getListOfDataFields(PyObject *self, PyObject * /*args*/) {
    Base * component = get_base(self);

    const sofa::helper::vector<BaseData*> dataFields = component->getDataFields();

    PyObject * pyList = PyList_New(dataFields.size());
    for (unsigned int i = 0; i < dataFields.size(); ++i) {
        PyList_SetItem(pyList, i, SP_BUILD_PYPTR(Data,BaseData,dataFields[i],false)) ;
    }

    return pyList;
}

static PyObject * Base_getListOfLinks(PyObject *self, PyObject * /*args*/) {
    Base * component = get_base(self);

    const sofa::helper::vector<BaseLink*> links = component->getLinks() ;

    PyObject * pyList = PyList_New(links.size());
    for (unsigned int i = 0; i < links.size(); ++i) {
        PyList_SetItem(pyList, i, SP_BUILD_PYPTR(Link, BaseLink, links[i], false)) ;
    }

    return pyList;
}

      void TriangleSetTopologyModifier::removeTrianglesWarning(sofa::helper::vector<unsigned int> &triangles)
      {
          m_container->setTriangleTopologyToDirty();


	/// sort vertices to remove in a descendent order
	std::sort( triangles.begin(), triangles.end(), std::greater<unsigned int>() );

	// Warning that these triangles will be deleted
	TrianglesRemoved *e=new TrianglesRemoved(triangles);
	addTopologyChange(e);
      }

static PyObject * Base_getDataFields(PyObject *self, PyObject * /*args*/) {
    Base * component = get_base(self);

    const sofa::helper::vector<BaseData*> dataFields = component->getDataFields();

    PyObject * pyDict = PyDict_New();
    for (size_t i=0; i<dataFields.size(); i++) {
        PyDict_SetItem(pyDict, PyString_FromString(dataFields[i]->getName().c_str()),
                       GetDataValuePython(dataFields[i]));
    }

    return pyDict;
}

      void TriangleSetTopologyModifier::removeTriangles(sofa::helper::vector< unsigned int >& triangles,
							const bool removeIsolatedEdges,
							const bool removeIsolatedPoints)
      {
	helper::ReadAccessor< Data< sofa::helper::vector<Triangle> > > m_triangle = m_container->d_triangle;
	for (unsigned int i = 0; i < triangles.size(); i++)
            {
              if( triangles[i] >= m_triangle.size())
                {
                  std::cout << "Error: TriangleSetTopologyModifier::removeTriangles: Triangle: "<< triangles[i] <<" is out of bound" << std::endl;
                  return;
                }
            }


          if (removeTrianglesPreconditions(triangles)) // Test if the topology will still fullfil the conditions if these triangles are removed.
            {
              /// add the topological changes in the queue
              removeTrianglesWarning(triangles);
              // inform other objects that the triangles are going to be removed
              propagateTopologicalChanges();
              // now destroy the old triangles.
              removeTrianglesProcess(  triangles ,removeIsolatedEdges, removeIsolatedPoints);

              m_container->checkTopology();
            }
          else
            {
              std::cout << " TriangleSetTopologyModifier::removeItems(), preconditions for removal are not fullfil. " << std::endl;
            }

      }

void HexahedronSetTopologyModifier::addHexahedraProcess(const sofa::helper::vector< Hexahedron > &hexahedra)
{
    for(unsigned int i = 0; i < hexahedra.size(); ++i)
    {
        addHexahedronProcess(hexahedra[i]);
    }
}

      void TriangleSetTopologyModifier::renumberPointsProcess( const sofa::helper::vector<unsigned int> &index,
							       const sofa::helper::vector<unsigned int> &inv_index,
							       const bool renumberDOF)
      {

	if(m_container->hasTriangles())
	{
	  if(m_container->hasTrianglesAroundVertex())
	  {
	    sofa::helper::vector< sofa::helper::vector< unsigned int > > trianglesAroundVertex_cp = m_container->m_trianglesAroundVertex;
	    for(unsigned int i=0; i<index.size(); ++i)
	    {
	      m_container->m_trianglesAroundVertex[i] = trianglesAroundVertex_cp[ index[i] ];
	    }
	  }
	  helper::WriteAccessor< Data< sofa::helper::vector<Triangle> > > m_triangle = m_container->d_triangle;

	  for(unsigned int i=0; i<m_triangle.size(); ++i)
	  {
	    m_triangle[i][0] = inv_index[ m_triangle[i][0] ];
	    m_triangle[i][1] = inv_index[ m_triangle[i][1] ];
	    m_triangle[i][2] = inv_index[ m_triangle[i][2] ];
	  }
	}

	// call the parent's method
	EdgeSetTopologyModifier::renumberPointsProcess( index, inv_index, renumberDOF );
      }

void HexahedronSetTopologyModifier::removeQuadsProcess( const sofa::helper::vector<unsigned int> &indices,
        const bool removeIsolatedEdges,
        const bool removeIsolatedPoints)
{
    if(!m_container->hasQuads()) // this method should only be called when quads exist
        return;

    if(m_container->hasQuadsInHexahedron())
    {
        if(!m_container->hasHexahedraAroundQuad())
            m_container->createHexahedraAroundQuadArray();

        unsigned int lastQuad = m_container->getNumberOfQuads() - 1;
        for(unsigned int i=0; i<indices.size(); ++i, --lastQuad)
        {
            for(sofa::helper::vector<unsigned int>::iterator itt=m_container->m_hexahedraAroundQuad[lastQuad].begin();
                itt!=m_container->m_hexahedraAroundQuad[lastQuad].end(); ++itt)
            {
                unsigned int quadIndex=m_container->getQuadIndexInHexahedron(m_container->m_quadsInHexahedron[*itt],lastQuad);
                m_container->m_quadsInHexahedron[*itt][quadIndex]=indices[i];
            }

            // updating the quad shell itself (change the old index for the new one)
            m_container->m_hexahedraAroundQuad[ indices[i] ] = m_container->m_hexahedraAroundQuad[ lastQuad ];
        }
        m_container->m_hexahedraAroundQuad.resize( m_container->m_hexahedraAroundQuad.size() - indices.size() );
    }

    // call the parent's method.
    QuadSetTopologyModifier::removeQuadsProcess( indices, removeIsolatedEdges, removeIsolatedPoints );
}

void HexahedronSetTopologyModifier::removePointsProcess(const sofa::helper::vector<unsigned int> &indices,
        const bool removeDOF)
{
    if(m_container->hasHexahedra())
    {
        if(!m_container->hasHexahedraAroundVertex())
        {
            m_container->createHexahedraAroundVertexArray();
        }

        helper::WriteAccessor< Data< sofa::helper::vector<Hexahedron> > > m_hexahedron = m_container->d_hexahedron;

        unsigned int lastPoint = m_container->getNbPoints() - 1;
        for(unsigned int i = 0; i < indices.size(); ++i, --lastPoint)
        {
            // updating the edges connected to the point replacing the removed one:
            // for all edges connected to the last point
            for(sofa::helper::vector<unsigned int>::iterator itt=m_container->m_hexahedraAroundVertex[lastPoint].begin();
                itt!=m_container->m_hexahedraAroundVertex[lastPoint].end(); ++itt)
            {
                unsigned int vertexIndex = m_container->getVertexIndexInHexahedron(m_hexahedron[*itt], lastPoint);
                m_hexahedron[*itt][ vertexIndex] = indices[i];
            }

            // updating the edge shell itself (change the old index for the new one)
            m_container->m_hexahedraAroundVertex[ indices[i] ] = m_container->m_hexahedraAroundVertex[ lastPoint ];
        }

        m_container->m_hexahedraAroundVertex.resize( m_container->m_hexahedraAroundVertex.size() - indices.size() );
    }

    // Important : the points are actually deleted from the mechanical object's state vectors iff (removeDOF == true)
    // call the parent's method.
    QuadSetTopologyModifier::removePointsProcess(  indices, removeDOF );
}

void ManifoldTriangleSetTopologyModifier::createRemovingEdgesFutureModifications (const sofa::helper::vector <unsigned int> items)
{
    EdgesInTriangle EdgesInTriangleArray;
    bool test = true;

    for (unsigned int  i = 0; i < items.size(); i++)
    {
        EdgesInTriangleArray = m_container->getEdgesInTriangle( items[i] );

        for (unsigned int j =0; j < 3 ; j++)
        {

            for (unsigned int k =0; k< m_modificationsEdge.size(); k++)
            {
                if (EdgesInTriangleArray[j] == m_modificationsEdge[k])
                {
                    test = false;
                    break;
                }
            }

            if (test)
            {
                m_modificationsEdge.push_back(EdgesInTriangleArray[j]);
            }
            test = true;
        }
    }
}

void QuadSetTopologyModifier::renumberPointsProcess( const sofa::helper::vector<unsigned int> &index,
        const sofa::helper::vector<unsigned int> &inv_index,
        const bool renumberDOF)
{
    if(m_container->hasQuads())
    {
        helper::WriteAccessor< Data< sofa::helper::vector<Quad> > > m_quad = m_container->d_quad;

        if(m_container->hasQuadsAroundVertex())
        {
            sofa::helper::vector< sofa::helper::vector< unsigned int > > quadsAroundVertex_cp = m_container->m_quadsAroundVertex;
            for(unsigned int i=0; i<index.size(); ++i)
            {
                m_container->m_quadsAroundVertex[i] = quadsAroundVertex_cp[ index[i] ];
            }
        }

        for(unsigned int i=0; i<m_quad.size(); ++i)
        {
            m_quad[i][0]  = inv_index[ m_quad[i][0]  ];
            m_quad[i][1]  = inv_index[ m_quad[i][1]  ];
            m_quad[i][2]  = inv_index[ m_quad[i][2]  ];
            m_quad[i][3]  = inv_index[ m_quad[i][3]  ];
        }
    }

    // call the parent's method
    EdgeSetTopologyModifier::renumberPointsProcess( index, inv_index, renumberDOF );
}

void HexahedronSetTopologyModifier::renumberPointsProcess( const sofa::helper::vector<unsigned int> &index, const sofa::helper::vector<unsigned int> &inv_index, const bool renumberDOF)
{
    if(m_container->hasHexahedra())
    {
        if(m_container->hasHexahedraAroundVertex())
        {
            sofa::helper::vector< sofa::helper::vector< unsigned int > > hexahedraAroundVertex_cp = m_container->m_hexahedraAroundVertex;
            for(unsigned int i=0; i<index.size(); ++i)
            {
                m_container->m_hexahedraAroundVertex[i] = hexahedraAroundVertex_cp[ index[i] ];
            }
        }

        helper::WriteAccessor< Data< sofa::helper::vector<Hexahedron> > > m_hexahedron = m_container->d_hexahedron;

        for(unsigned int i=0; i<m_hexahedron.size(); ++i)
        {
            m_hexahedron[i][0]  = inv_index[ m_hexahedron[i][0]  ];
            m_hexahedron[i][1]  = inv_index[ m_hexahedron[i][1]  ];
            m_hexahedron[i][2]  = inv_index[ m_hexahedron[i][2]  ];
            m_hexahedron[i][3]  = inv_index[ m_hexahedron[i][3]  ];
            m_hexahedron[i][4]  = inv_index[ m_hexahedron[i][4]  ];
            m_hexahedron[i][5]  = inv_index[ m_hexahedron[i][5]  ];
            m_hexahedron[i][6]  = inv_index[ m_hexahedron[i][6]  ];
            m_hexahedron[i][7]  = inv_index[ m_hexahedron[i][7]  ];
        }
    }

    // call the parent's method.
    QuadSetTopologyModifier::renumberPointsProcess( index, inv_index, renumberDOF );
}

      void TriangleSetTopologyModifier::movePointsProcess (const sofa::helper::vector <unsigned int>& id,
							   const sofa::helper::vector< sofa::helper::vector< unsigned int > >& ancestors,
							   const sofa::helper::vector< sofa::helper::vector< double > >& coefs,
							   const bool moveDOF)
      {
	(void)moveDOF;
	unsigned int nbrVertex = id.size();
	bool doublet;
	sofa::helper::vector< unsigned int > trianglesAroundVertex2Move;
	sofa::helper::vector< Triangle > trianglesArray;


	// Step 1/4 - Creating trianglesAroundVertex to moved due to moved points:
	for (unsigned int i = 0; i<nbrVertex; ++i)
	{
	  const sofa::helper::vector <unsigned int>& trianglesAroundVertex = m_container->getTrianglesAroundVertex( id[i] );

	  for (unsigned int j = 0; j<trianglesAroundVertex.size(); ++j)
	  {
	    doublet = false;
	
	    for (unsigned int k =0; k<trianglesAroundVertex2Move.size(); ++k) //Avoid double
	    {
	      if (trianglesAroundVertex2Move[k] == trianglesAroundVertex[j])
	      {
		doublet = true;
		break;
	      }
	    }
	    
	    if(!doublet)
	      trianglesAroundVertex2Move.push_back (trianglesAroundVertex[j]);
	    
	  }
	}
	
	std::sort( trianglesAroundVertex2Move.begin(), trianglesAroundVertex2Move.end(), std::greater<unsigned int>() );


	// Step 2/4 - Create event to delete all elements before moving and propagate it:
	TrianglesMoved_Removing *ev1 = new TrianglesMoved_Removing (trianglesAroundVertex2Move);
	this->addTopologyChange(ev1);
	propagateTopologicalChanges();


	// Step 3/4 - Physically move all dof:
	PointSetTopologyModifier::movePointsProcess (id, ancestors, coefs);
	
	
	// Step 4/4 - Create event to recompute all elements concerned by moving and propagate it:

	// Creating the corresponding array of Triangles for ancestors
	for (unsigned int i = 0; i<trianglesAroundVertex2Move.size(); i++)
	  trianglesArray.push_back (m_container->getTriangleArray()[ trianglesAroundVertex2Move[i] ]);

	TrianglesMoved_Adding *ev2 = new TrianglesMoved_Adding (trianglesAroundVertex2Move, trianglesArray);
	this->addTopologyChange(ev2); // This event should be propagated with global workflow
      }

void DynamicSparseGridTopologyModifier::addHexahedraProcess ( const sofa::helper::vector< Hexahedron > &hexahedra, const sofa::helper::vector< unsigned int> &indices )
{
    assert( hexahedra.size() == indices.size());

    unsigned int hexaSize = m_DynContainer->getNumberOfHexahedra(); // Get the size before adding elements
    HexahedronSetTopologyModifier::addHexahedraProcess ( hexahedra );
    helper::vector<core::topology::BaseMeshTopology::HexaID>& iirg = *m_DynContainer->idxInRegularGrid.beginEdit();

    std::map< unsigned int, core::topology::BaseMeshTopology::HexaID> &idrg2topo=*m_DynContainer->idInRegularGrid2IndexInTopo.beginEdit();
    for ( unsigned int i = 0; i < hexahedra.size(); i++ )  // For each element
    {
        iirg[hexaSize + i] = indices[i];
        idrg2topo.insert( std::make_pair ( indices[i], hexaSize + i ) );

        //TODO// init the values too ...
    }
    m_DynContainer->idInRegularGrid2IndexInTopo.endEdit();
    m_DynContainer->idxInRegularGrid.endEdit();
}

void SimpleTesselatedTetraTopologicalMapping::removeOutputTetrahedra( const sofa::helper::vector<unsigned int>& index )
{

    helper::vector< fixed_array<int, 8> >& tetrahedraMappedFromTetraData = *(tetrahedraMappedFromTetra.beginEdit());
    helper::vector<int>& tetraSourceData = *(tetraSource.beginEdit());

    int last = tetraSourceData.size() -1;
    for (unsigned int i = 0; i < index.size(); ++i)
    {
        swapOutputTetrahedra( index[i], last );
        int source = tetraSourceData[last];
        if (source != -1)
        {
			int nbt = 0;

            for (int j=0; j<8; ++j)
			{
                if (tetrahedraMappedFromTetraData[source][j] == last)
                {
                    tetrahedraMappedFromTetraData[source][j] = -1;
                }
                else if (tetrahedraMappedFromTetraData[source][j] != -1)
                    ++nbt;
			}
            if (nbt == 0) // we need to remove the source tetra
            {
//				sout << "SimpleTesselatedTetraTopologicalMapping: source tetra "<<source<<" needs to be removed."<<sendl;
                tetrahedraToRemove.insert(source);
            }
            else
            {
//			    sout << "SimpleTesselatedTetraTopologicalMapping: source tetra "<<source<<" now has "<<nbt<<" / 8 childs."<<sendl;
            }
            --last;
        }
    }

    tetraSourceData.resize( tetraSourceData.size() - index.size() );

    tetrahedraMappedFromTetra.endEdit();
    tetraSource.endEdit();
}