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

*clap*

Browse Source
This commit is contained in:
papush! 2022-03-23 15:27:16 +01:00
parent 26b1175924
commit 36194f715e
10 changed files with 210 additions and 37 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

7
CMakeLists.txt
View File

@ -12,6 +12,7 @@ set(VTK_COMPONENTS
VTK::IOGeometry
VTK::IOXML
VTK::FiltersModeling
VTK::FiltersGeometry
VTK::vtksys)
set(ENABLE_VIEWER OFF CACHE BOOL "Enable the 3D viewer, depends on Qt.")
if(ENABLE_VIEWER)
@ -62,7 +63,11 @@ target_sources(pfe PRIVATE
src/remove_external_cells_filter.cc
src/remove_external_cells_filter.h
src/relaxation_filter.cc
src/relaxation_filter.h)
src/relaxation_filter.h
src/max_distance_filter.cc
src/max_distance_filter.h
src/closest_polymesh_point.cc
src/closest_polymesh_point.h)
target_link_libraries(pfe PRIVATE ${VTK_COMPONENTS})

12
data/tetrahedron.obj Normal file
View File

@ -0,0 +1,12 @@
# Blender v2.79 (sub 0) OBJ File: ''
# www.blender.org
o Solid
v 0.000000 0.769688 0.000000
v 0.000000 0.000000 -0.769688
v 0.769688 0.000000 0.000000
v 0.000000 0.000000 0.000000
s off
f 1 2 3
f 1 3 4
f 1 4 2
f 2 4 3

60
src/closest_polymesh_point.cc Normal file
View File

@ -0,0 +1,60 @@
#include "closest_polymesh_point.h"
#include "point_tris_dist.h"
#include <limits>
vtkIdType findClosestPoint(const double *point, vtkPointSet *pointSet) {
vtkIdType minPoint = 0;
double minDist2 = std::numeric_limits<double>::infinity();
for (vtkIdType i = 0; i < pointSet->GetNumberOfPoints(); i++) {
double point2[3];
pointSet->GetPoint(i, point2);
double dist2 = vtkMath::Distance2BetweenPoints(point, point2);
if (dist2 < minDist2) {
minDist2 = dist2;
minPoint = i;
}
}
return minPoint;
}
void closestPolyMeshPoint(vtkPolyData *polyMesh,
const double *point,
vtkKdTree *kdTree,
vtkStaticCellLinks *links,
double *toClosestPoint,
double *distanceToClosestPoint) {
vtkIdType closestPolyMeshPoint;
{
double dummy;
// closestPolyMeshPoint = kdTree->FindClosestPoint((double *) point, dummy);
closestPolyMeshPoint = findClosestPoint(point, polyMesh);
}
vtkIdType *cellIds = links->GetCells(closestPolyMeshPoint);
vtkIdType nCells = links->GetNcells(closestPolyMeshPoint);
if (nCells == 0) {
*distanceToClosestPoint = -1;
return;
}
double minDist = std::numeric_limits<double>::infinity();
/* For each tri the closest point belongs to. */
for (vtkIdType i = 0; i < nCells; i++) {
vtkIdType cellId = cellIds[i];
double direction[3];
double dist = pointTriangleDistance(
(double *) point, polyMesh->GetCell(cellId), direction);
/* Find the closest one. */
if (dist < minDist) {
minDist = dist;
toClosestPoint[0] = direction[0];
toClosestPoint[1] = direction[1];
toClosestPoint[2] = direction[2];
}
}
vtkMath::MultiplyScalar(toClosestPoint, minDist);
*distanceToClosestPoint = minDist;
}

17
src/closest_polymesh_point.h Normal file
View File

@ -0,0 +1,17 @@
#ifndef CLOSEST_POLYMESH_POINT_H
#define CLOSEST_POLYMESH_POINT_H
#include <vtkPolyData.h>
#include <vtkKdTree.h>
#include <vtkStaticCellLinks.h>
void closestPolyMeshPoint(vtkPolyData *polyMesh,
const double *point,
vtkKdTree *kdTree,
vtkStaticCellLinks *links,
double *toClosestPoint,
double *distanceToClosestPoint);
#endif

10
src/dihedral_angles_filter.cc
View File

@ -1,5 +1,7 @@
#include "dihedral_angles_filter.h"
#include <algorithm>
#include <limits>
#include <vtkUnstructuredGrid.h>
#include <vtkPointData.h>
#include <vtkCellData.h>
@ -36,6 +38,9 @@ vtkTypeBool DihedralAnglesFilter::RequestData(
minDihedralAngleArray->SetNumberOfTuples(input->GetNumberOfCells());
double *minDihedralAngleBase = minDihedralAngleArray->GetPointer(0);
AverageMinDegrees = 0;
MinMinDegrees = std::numeric_limits<double>::infinity();
size_t i = 0;
vtkCellIterator *it = input->NewCellIterator();
for (it->InitTraversal(); !it->IsDoneWithTraversal(); it->GoToNextCell()) {
@ -62,9 +67,14 @@ vtkTypeBool DihedralAnglesFilter::RequestData(
minDihedralAngleBase[i / 6] =
*std::min_element(dihedralAnglesBase + i,
dihedralAnglesBase + i + 6);
AverageMinDegrees += minDihedralAngleBase[i / 6];
MinMinDegrees = std::min(MinMinDegrees, minDihedralAngleBase[i / 6]);
i += 6;
}
AverageMinDegrees /= input->GetNumberOfCells();
AverageMinDegrees = AverageMinDegrees * 180. / 3.141592653589793;
MinMinDegrees = MinMinDegrees * 180. / 3.141592653589793;
cellData->AddArray((vtkAbstractArray *) dihedralAnglesArray);
cellData->AddArray((vtkAbstractArray *) minDihedralAngleArray);
return true;

4
src/dihedral_angles_filter.h
View File

@ -12,9 +12,13 @@ public:
vtkTypeBool RequestData(vtkInformation *request,
vtkInformationVector **inputVector,
vtkInformationVector *outputVector) override;
vtkGetMacro(AverageMinDegrees, double);
vtkGetMacro(MinMinDegrees, double);
protected:
~DihedralAnglesFilter() override = default;
double AverageMinDegrees;
double MinMinDegrees;
};

14
src/main.cc
View File

@ -3,7 +3,9 @@
#include "surface_points_filter.h"
#include "project_surface_points_on_poly.h"
#include "relaxation_filter.h"
#include "max_distance_filter.h"
#include <vtkGeometryFilter.h>
#include <vtkCellArrayIterator.h>
#include <vtkCellData.h>
#include <vtkDataSetReader.h>
@ -98,12 +100,24 @@ int main(int argc, char **argv) {
dihedralAnglesFilter->SetInputConnection(
relaxationFilter->GetOutputPort());
vtkNew<vtkGeometryFilter> geometryFilter;
geometryFilter->SetInputConnection(dihedralAnglesFilter->GetOutputPort());
vtkNew<MaxDistanceFilter> maxDistanceFilter;
maxDistanceFilter->SetInputConnection(0, geometryFilter->GetOutputPort());
maxDistanceFilter->SetInputConnection(1, polyMeshReader->GetOutputPort());
vtkNew<vtkXMLUnstructuredGridWriter> writer;
writer->SetInputConnection(dihedralAnglesFilter->GetOutputPort());
writer->SetDataModeToAscii();
writer->SetFileName("out.vtu");
writer->Write();
maxDistanceFilter->Update();
std::cerr << "Max distance: " << maxDistanceFilter->GetMaxDist() << "\n"
<< "Average min angle: " << dihedralAnglesFilter->GetAverageMinDegrees() << "\n"
<< "Min min angle: " << dihedralAnglesFilter->GetMinMinDegrees() << "\n";
#ifdef USE_VIEWER
/* Volume rendering properties */
vtkNew<vtkOpenGLProjectedTetrahedraMapper> volumeMapper;

63
src/max_distance_filter.cc Normal file
View File

@ -0,0 +1,63 @@
#include "max_distance_filter.h"
#include "closest_polymesh_point.h"
#include <limits>
#include <vtkPolyData.h>
#include <vtkInformation.h>
#include <vtkKdTree.h>
vtkStandardNewMacro(MaxDistanceFilter);
MaxDistanceFilter::MaxDistanceFilter() {
SetNumberOfInputPorts(2);
SetNumberOfOutputPorts(0);
}
int MaxDistanceFilter::FillInputPortInformation(int port,
vtkInformation *info) {
if (port == 0) {
vtkPolyDataAlgorithm::FillInputPortInformation(port, info);
} else if (port == 1) {
info->Set(vtkDataObject::DATA_TYPE_NAME(), "vtkPolyData");
}
return 1;
}
vtkTypeBool MaxDistanceFilter::RequestData(vtkInformation *request,
vtkInformationVector **inputVector,
vtkInformationVector *outputVector) {
(void) request;
(void) outputVector;
vtkPolyData* input1 = vtkPolyData::GetData(inputVector[0]);
vtkPolyData* input2 = vtkPolyData::GetData(inputVector[1]);
vtkNew<vtkKdTree> tree1;
tree1->BuildLocatorFromPoints(input1->GetPoints());
vtkNew<vtkKdTree> tree2;
tree2->BuildLocatorFromPoints(input2->GetPoints());
vtkNew<vtkStaticCellLinks> links1;
links1->BuildLinks(input1);
vtkNew<vtkStaticCellLinks> links2;
links2->BuildLinks(input2);
MaxDist = 0;
for (vtkIdType i = 0; i < input1->GetNumberOfPoints(); i++) {
double point[3];
input1->GetPoint(i, point);
double vec[3];
double dist;
closestPolyMeshPoint(input2, point, tree2, links2, vec, &dist);
MaxDist = std::max(dist, MaxDist);
}
for (vtkIdType i = 0; i < input2->GetNumberOfPoints(); i++) {
double point[3];
input2->GetPoint(i, point);
double vec[3];
double dist;
closestPolyMeshPoint(input1, point, tree1, links1, vec, &dist);
MaxDist = std::max(dist, MaxDist);
}
return true;
}

23
src/max_distance_filter.h Normal file
View File

@ -0,0 +1,23 @@
#ifndef MAX_DISTANCE_FILTER_H
#define MAX_DISTANCE_FILTER_H
#include <vtkPolyDataAlgorithm.h>
class MaxDistanceFilter : public vtkPolyDataAlgorithm {
public:
static MaxDistanceFilter *New();
vtkTypeMacro(MaxDistanceFilter, vtkPolyDataAlgorithm);
MaxDistanceFilter();
int FillInputPortInformation(int, vtkInformation *info) override;
vtkTypeBool RequestData(vtkInformation *request,
vtkInformationVector **inputVector,
vtkInformationVector *outputVector) override;
vtkGetMacro(MaxDist, double);
protected:
double MaxDist;
};
#endif

37
src/project_surface_points_on_poly.cc
View File

@ -1,6 +1,5 @@
#include "point_tris_dist.h"
#include "project_surface_points_on_poly.h"
#include "closest_polymesh_point.h"
#include <vtkUnstructuredGrid.h>
#include <vtkPointData.h>
@ -36,40 +35,6 @@ int ProjectSurfacePointsOnPoly::FillInputPortInformation(int port, vtkInformatio
}
static void closestPolyMeshPoint(vtkPolyData *polyMesh,
const double *point,
vtkKdTree *kdTree,
vtkStaticCellLinks *links,
double *toClosestPoint,
double *distanceToClosestPoint) {
vtkIdType closestPolyMeshPoint;
{
double dummy;
closestPolyMeshPoint = kdTree->FindClosestPoint((double *) point, dummy);
}
vtkIdType *cellIds = links->GetCells(closestPolyMeshPoint);
vtkIdType nCells = links->GetNcells(closestPolyMeshPoint);
double minDist = std::numeric_limits<double>::infinity();
/* For each tri the closest point belongs to. */
for (vtkIdType i = 0; i < nCells; i++) {
vtkIdType cellId = cellIds[i];
double direction[3];
double dist = pointTriangleDistance(
(double *) point, polyMesh->GetCell(cellId), direction);
/* Find the closest one. */
if (dist < minDist) {
minDist = dist;
toClosestPoint[0] = direction[0];
toClosestPoint[1] = direction[1];
toClosestPoint[2] = direction[2];
}
}
vtkMath::MultiplyScalar(toClosestPoint, minDist);
*distanceToClosestPoint = minDist;
}
void findPointsWithinRadius(double radius, vtkPointSet *pointSet,
const double *point, vtkIdList *points) {
for (vtkIdType i = 0; i < pointSet->GetNumberOfPoints(); i++) {