-
Notifications
You must be signed in to change notification settings - Fork 310
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
[SofaBaseTopology] Fix SparseGridTopology initialization with an input mesh #670
Changes from 2 commits
7f55b08
6458f7f
92cd417
c9a86bd
0f2e89f
0c92a17
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
Original file line number | Diff line number | Diff line change |
---|---|---|
|
@@ -273,22 +273,18 @@ void SparseGridTopology::buildAsFinest( ) | |
// seqPoints, seqHexahedra.getValue(), nbPoints | ||
if(_filename.length() > 4 && _filename.compare(_filename.length()-4, 4, ".obj")==0) | ||
{ | ||
// std::cout << "SparseGridTopology: using mesh "<<_filename<<std::endl; | ||
buildFromTriangleMesh(_filename); | ||
} | ||
else if(_filename.length() > 6 && _filename.compare(_filename.length()-6, 6, ".trian")==0) | ||
{ | ||
// std::cout << "SparseGridTopology: using mesh "<<_filename<<std::endl; | ||
buildFromTriangleMesh(_filename); | ||
} | ||
else if(_filename.length() > 4 && _filename.compare(_filename.length()-4, 4, ".stl")==0) | ||
{ | ||
// std::cout << "SparseGridTopology: using mesh "<<_filename<<std::endl; | ||
buildFromTriangleMesh(_filename); | ||
} | ||
else if(_filename.length() > 4 && _filename.compare(_filename.length()-4, 4, ".raw")==0) | ||
{ | ||
// std::cout << "SparseGridTopology: using mesh "<<_filename<<std::endl; | ||
_usingMC = true; | ||
|
||
buildFromRawVoxelFile(_filename); | ||
|
@@ -785,13 +781,9 @@ void SparseGridTopology::buildFromTriangleMesh(const std::string& filename) | |
SReal xMin, xMax, yMin, yMax, zMin, zMax; | ||
computeBoundingBox(mesh->getVertices(), xMin, xMax, yMin, yMax, zMin, zMax); | ||
|
||
// increase the box a little | ||
Vector3 diff ( xMax-xMin, yMax - yMin, zMax - zMin ); | ||
diff /= 100.0; | ||
|
||
_min.setValue(Vector3( xMin - diff[0], yMin - diff[1], zMin - diff[2] )); | ||
_max.setValue(Vector3( xMax + diff[0], yMax + diff[1], zMax + diff[2] )); | ||
sout << "BBox size: " << (_max.getValue() - _min.getValue()) << sendl; | ||
_min.setValue(Vector3( xMin, yMin, zMin )); | ||
_max.setValue(Vector3( xMax, yMax, zMax )); | ||
msg_info() << "BBox size: " << (_max.getValue() - _min.getValue()); | ||
} | ||
|
||
// if cellWidth is given, update n | ||
|
@@ -800,7 +792,7 @@ void SparseGridTopology::buildFromTriangleMesh(const std::string& filename) | |
SReal w = _cellWidth.getValue(); | ||
Vector3 diff = _max.getValue() - _min.getValue(); | ||
setN(Vec3i((int)ceil(diff[0] / w)+1, (int)ceil(diff[1] / w)+1, (int)ceil(diff[2] / w)+1)); | ||
sout << "Grid size: " << n.getValue() << sendl; | ||
msg_info() << "Grid size: " << n.getValue(); | ||
} | ||
|
||
_regularGrid->setSize(getNx(),getNy(),getNz()); | ||
|
@@ -811,7 +803,7 @@ void SparseGridTopology::buildFromTriangleMesh(const std::string& filename) | |
|
||
buildFromRegularGridTypes(_regularGrid, regularGridTypes); | ||
|
||
sout << "Mesh Loaded" << sendl; | ||
msg_info() << "Mesh Loaded"; | ||
|
||
delete mesh; | ||
} | ||
|
@@ -840,38 +832,45 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
const helper::vector< Vector3 >& vertices = mesh->getVertices(); | ||
const size_t vertexSize = vertices.size(); | ||
helper::vector< int > verticesHexa(vertexSize); | ||
helper::vector< bool > facetDone; | ||
const Vector3 delta = (regularGrid->getDx() + regularGrid->getDy() + regularGrid->getDz()) / 2; | ||
|
||
// Compute the grid element for each mesh vertex | ||
for (size_t i = 0; i < vertexSize; ++i) | ||
{ | ||
const Vector3& vertex = vertices[i]; | ||
int index = regularGrid->findHexa(vertex); | ||
|
||
if (index > 0) | ||
regularGridTypes[index] = BOUNDARY; | ||
|
||
// Case where 'findHexa' did not find the right hexa | ||
// Here we test the case where the point is close the surface (delta /2) | ||
// Useful when the point is on the boundary | ||
if (index == -1) | ||
{ | ||
Vector3 vertex2 = vertex; | ||
const Vector3 gmin = regularGrid->getMin(); | ||
const Vector3 gmax = regularGrid->getMax(); | ||
|
||
if (vertex2[0] < gmin[0]) | ||
vertex2[0] = gmin[0]; | ||
else if (vertex2[0] > gmax[0]) | ||
if (vertex2[0] <= gmin[0]) | ||
vertex2[0] = gmin[0] + delta[0]; | ||
else if (vertex2[0] >= gmax[0]) | ||
vertex2[0] = gmax[0] - delta[0]; | ||
|
||
if (vertex2[1] < gmin[1]) | ||
vertex2[1] = gmin[1]; | ||
else if (vertex2[1] > gmax[1]) | ||
if (vertex2[1] <= gmin[1]) | ||
vertex2[1] = gmin[1] + delta[1]; | ||
else if (vertex2[1] >= gmax[1]) | ||
vertex2[1] = gmax[1] - delta[1]; | ||
|
||
if (vertex2[2] < gmin[2]) | ||
vertex2[2] = gmin[2]; | ||
else if (vertex2[2] > gmax[2]) | ||
if (vertex2[2] <= gmin[2]) | ||
vertex2[2] = gmin[2] + delta[2]; | ||
else if (vertex2[2] >= gmax[2]) | ||
vertex2[2] = gmax[2] - delta[2]; | ||
|
||
index = regularGrid->findHexa(vertex2); | ||
|
||
if(index == -1) | ||
msg_error() << "vertex "<<i<<" not found in hexahedral topology"; | ||
} | ||
|
||
verticesHexa[i] = index; | ||
|
@@ -897,6 +896,7 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
const Vector3 i0 = regularGrid->getCubeCoordinate(c0); | ||
const Vector3 i1 = regularGrid->getCubeCoordinate(c1); | ||
const Vector3 i2 = regularGrid->getCubeCoordinate(c2); | ||
|
||
const Vector3 iMin( std::min(i0[0],std::min(i1[0],i2[0])), | ||
std::min(i0[1],std::min(i1[1],i2[1])), | ||
std::min(i0[2],std::min(i1[2],i2[2]))); | ||
|
@@ -912,13 +912,16 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
{ | ||
// if already inserted discard | ||
unsigned int index = regularGrid->getCubeIndex(x,y,z); | ||
|
||
if(regularGridTypes[index]==BOUNDARY) | ||
continue; | ||
|
||
Hexa c = regularGrid->getHexaCopy(index); | ||
|
||
CubeCorners corners; | ||
for(int k=0; k<8; ++k) | ||
corners[k] = regularGrid->getPoint( c[k] ); | ||
|
||
#ifdef SOFA_NEW_HEXA | ||
Vector3 cubeDiagonal = corners[6] - corners[0]; | ||
#else | ||
|
@@ -932,6 +935,7 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
|
||
// Scale the triangle to the unit cube matching | ||
float points[3][3]; | ||
|
||
for (unsigned short w=0; w<3; ++w) | ||
{ | ||
points[0][w] = (float) ((A[w]-cubeCenter[w])/cubeDiagonal[w]); | ||
|
@@ -952,6 +956,7 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
} | ||
} | ||
|
||
|
||
// TODO: regarder les cellules pleines, et les ajouter | ||
vector<bool> alreadyTested(regularGrid->getNbHexahedra(),false); | ||
std::stack< Vec3i > seed; | ||
|
@@ -963,6 +968,7 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
launchPropagationFromSeed(Vec3i(0,y,z), regularGrid, regularGridTypes, alreadyTested,seed ); | ||
launchPropagationFromSeed(Vec3i(regularGrid->getNx()-2,y,z), regularGrid, regularGridTypes, alreadyTested,seed ); | ||
} | ||
|
||
// y==0 and y=ny-2 | ||
for(int x=0; x<regularGrid->getNx()-1; ++x) | ||
{ | ||
|
@@ -971,6 +977,7 @@ void SparseGridTopology::voxelizeTriangleMesh(helper::io::Mesh* mesh, | |
launchPropagationFromSeed(Vec3i(x,0,z), regularGrid, regularGridTypes, alreadyTested,seed ); | ||
launchPropagationFromSeed(Vec3i(x,regularGrid->getNy()-2,z), regularGrid, regularGridTypes, alreadyTested,seed ); | ||
} | ||
|
||
// z==0 and z==Nz-2 | ||
for(int y=0; y<regularGrid->getNy()-1; ++y) | ||
{ | ||
|
@@ -1419,8 +1426,8 @@ void SparseGridTopology::buildVirtualFinerLevels() | |
_virtualFinerLevels[0]->_fillWeighted.setValue( _fillWeighted.getValue() ); | ||
_virtualFinerLevels[0]->init(); | ||
|
||
sout<<"SparseGridTopology "<<getName()<<" buildVirtualFinerLevels : "; | ||
sout<<"("<<newnx<<"x"<<newny<<"x"<<newnz<<") -> "<< _virtualFinerLevels[0]->getNbHexahedra() <<" elements , "; | ||
msg_info()<<"SparseGridTopology "<<getName()<<" buildVirtualFinerLevels : " << msgendl; | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. There is no needs for msgendl at he end of message it is only to separate lines within a message. |
||
msg_info()<<"("<<newnx<<"x"<<newny<<"x"<<newnz<<") -> "<< _virtualFinerLevels[0]->getNbHexahedra() <<" elements , " << msgendl; | ||
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. It is better (more efficient, more clear) to have a single message with multiple lines in it instead of two message. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Here is an example of one message without line breaks. msg_info() <<"("<<newnx<<"x"<<newny<<"x"<<newnz<<") -> "<<_virtualFinerLevels[0]->getNbHexahedra() <<" elements , "
<<"("<<_virtualFinerLevels[i]->getNx()<<"x"<<_virtualFinerLevels[i]->getNy()
<<"x"<<_virtualFinerLevels[i]->getNz()<<") -> "<< _virtualFinerLevels[i]->getNbHexahedra() <<" elements , " If you want linebreaks into a message...add msgendl where you want the break. There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. Does my latest commit do what you expected @damienmarchal ? |
||
|
||
for(int i=1; i<nb; ++i) | ||
{ | ||
|
@@ -1431,11 +1438,9 @@ void SparseGridTopology::buildVirtualFinerLevels() | |
|
||
_virtualFinerLevels[i]->init(); | ||
|
||
sout<<"("<<_virtualFinerLevels[i]->getNx()<<"x"<<_virtualFinerLevels[i]->getNy()<<"x"<<_virtualFinerLevels[i]->getNz()<<") -> "<< _virtualFinerLevels[i]->getNbHexahedra() <<" elements , "; | ||
msg_info()<<"("<<_virtualFinerLevels[i]->getNx()<<"x"<<_virtualFinerLevels[i]->getNy()<<"x"<<_virtualFinerLevels[i]->getNz()<<") -> "<< _virtualFinerLevels[i]->getNbHexahedra() <<" elements , "<<msgendl; | ||
} | ||
|
||
sout<<sendl; | ||
|
||
this->setFinerSparseGrid(_virtualFinerLevels[nb-1].get()); | ||
} | ||
|
||
|
@@ -1507,11 +1512,14 @@ int SparseGridTopology::findNearestCube(const Vector3& pos, SReal& fx, SReal &fy | |
|
||
helper::fixed_array<int,6> SparseGridTopology::findneighboorCubes( int indice ) | ||
{ | ||
sout<<"SparseGridTopology::findneighboorCubes : "<<indice<<" -> "<<_indicesOfCubeinRegularGrid[indice]<<sendl; | ||
sout<<_indicesOfRegularCubeInSparseGrid[ _indicesOfCubeinRegularGrid[indice] ] <<sendl; | ||
msg_info()<<"SparseGridTopology::findneighboorCubes : "<<indice<<" -> "<<_indicesOfCubeinRegularGrid[indice]; | ||
msg_info()<<_indicesOfRegularCubeInSparseGrid[ _indicesOfCubeinRegularGrid[indice] ]; | ||
|
||
helper::fixed_array<int,6> result; | ||
Vector3 c = _regularGrid->getCubeCoordinate( _indicesOfCubeinRegularGrid[indice] ); | ||
sout<<c<<sendl; | ||
|
||
msg_info()<<c; | ||
|
||
result[0] = c[0]<=0 ? -1 : _indicesOfRegularCubeInSparseGrid[ _regularGrid->getCubeIndex( (int)c[0]-1,(int)c[1],(int)c[2] )]; | ||
result[1] = c[0]>=getNx()-2 ? -1 : _indicesOfRegularCubeInSparseGrid[ _regularGrid->getCubeIndex( (int)c[0]+1,(int)c[1],(int)c[2] )]; | ||
result[2] = c[1]<=0 ? -1 : _indicesOfRegularCubeInSparseGrid[ _regularGrid->getCubeIndex( (int)c[0],(int)c[1]-1,(int)c[2] )]; | ||
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Not sure about this Delta value, If I understand well. If one of the direction is very big and the other are small. You could miss some points in the "big" resolution direction in your test below.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
this could explain the crashing scene, but I don't get why actually. Why could I miss some points in the "big" resolution direction?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
not sure to remember exactly what mean Dx, Dy and Dz as they are 3D.
But in a general way in 1 dimension, if for example: Dx = 1000 and Dy = Dz = 10.
delta = 1020/2 = 510.
So when looking at next point in Dx direction, if the scale is ~1000, your delta of 510 is too small.