ccolin/pfe
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

fix projection

Browse Source
This commit is contained in:
papush! 2022-03-14 14:30:06 +01:00
parent 6c2e3e8a89
commit 62e4cd1629
1 changed files with 42 additions and 21 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

63
src/project_surface_points_on_poly.cc
View File

@ -71,6 +71,19 @@ static void closestPolyMeshPoint(vtkPolyData *polyMesh,
}
void findPointsWithinRadius(double radius, vtkPointSet *pointSet,
const double *point, vtkIdList *points) {
for (vtkIdType i = 0; i < pointSet->GetNumberOfPoints(); i++) {
double vec[3];
vtkMath::Subtract(point, pointSet->GetPoint(i), vec);
double dist = vtkMath::Norm(vec);
if (dist < radius) {
points->InsertNextId(i);
}
}
}
static void moveSurfacePoint(vtkUnstructuredGrid *tetMesh,
vtkPolyData *polyMesh,
vtkIdType pointId,
@ -87,14 +100,15 @@ static void moveSurfacePoint(vtkUnstructuredGrid *tetMesh,
closestPolyMeshPoint(polyMesh, point, kdTree, links, direction, &distance);
vtkNew<vtkIdList> affectedPoints;
kdTree->FindPointsWithinRadius(distance * radiusScale, point, affectedPoints);
findPointsWithinRadius(distance * radiusScale, tetMesh, point, affectedPoints);
// kdTree->FindPointsWithinRadius(distance * radiusScale, point, affectedPoints);
if(distance > 0) for (vtkIdType j = 0; j < affectedPoints->GetNumberOfIds(); j++) {
vtkIdType affectedPointId = affectedPoints->GetId(j);
double affectedPoint[3];
tetMesh->GetPoint(affectedPointId, affectedPoint);
double dist2 = vtkMath::Distance2BetweenPoints(affectedPoint, point);
if(dist2 > (distance * radiusScale) * (distance * radiusScale)) continue;
// if(dist2 > (distance * radiusScale) * (distance * radiusScale)) continue;
double motion[3] = {direction[0], direction[1], direction[2]};
double factor = 1. - std::sqrt(dist2) / (distance * radiusScale);
@ -112,6 +126,29 @@ static void moveSurfacePoint(vtkUnstructuredGrid *tetMesh,
}
static void applyMotionVectors(vtkUnstructuredGrid *tetMesh,
vtkIdTypeArray *surfacePoints,
std::vector<double> &motionVectors,
const std::vector<unsigned> &numberOfAffectors) {
for (vtkIdType i = 0; i < tetMesh->GetNumberOfPoints(); i++) {
bool skip = false;
for (vtkIdType j = 0; j < surfacePoints->GetNumberOfTuples(); j++) {
if (i == surfacePoints->GetValue(j)) skip = true;
}
if (skip) continue;
if(numberOfAffectors[i] == 0) continue;
double point[3];
tetMesh->GetPoint(i, point);
double *motion = motionVectors.data() + (i * 3);
vtkMath::MultiplyScalar(motion, 1. / numberOfAffectors[i]);
vtkMath::Subtract(point, motion, point);
tetMesh->GetPoints()->SetPoint(i, point);
tetMesh->Modified();
}
}
vtkTypeBool ProjectSurfacePointsOnPoly::RequestData(
vtkInformation *request,
vtkInformationVector **inputVector,
@ -144,24 +181,8 @@ vtkTypeBool ProjectSurfacePointsOnPoly::RequestData(
numberOfAffectors,
kdTree, links);
}
for (vtkIdType i = 0; i < tetMesh->GetNumberOfPoints(); i++) {
bool skip = false;
for (vtkIdType j = 0; j < surfacePoints->GetNumberOfTuples(); j++) {
if (i == surfacePoints->GetValue(j)) skip = true;
}
if (skip) continue;
if(numberOfAffectors[i] != 0) {
double point[3];
tetMesh->GetPoint(i, point);
double *motion = motionVectors.data() + (i * 3);
vtkMath::MultiplyScalar(motion, 1. / numberOfAffectors[i]);
vtkMath::Subtract(point, motion, point);
tetMesh->GetPoints()->SetPoint(i, point);
tetMesh->Modified();
}
}
applyMotionVectors(tetMesh, surfacePoints,
motionVectors, numberOfAffectors);
return true;
}