#include "analysis/dihedral_angles_filter.h" #include "analysis/max_distance_filter.h" #include "fitting/remove_external_cells_filter.h" #include "fitting/project_surface_points_on_poly.h" #include "fitting/relaxation_filter.h" #include "surface_points_filter.h" #include #include #include #include #include #include #include #include #include #include #include #include #include #include #include #ifdef USE_VIEWER #include #include #include #include #include #include #include #include #include #include #endif // USE_VIEWER #include vtkAlgorithm *readerFromFileName(const char *fileName) { std::string extension = vtksys::SystemTools::GetFilenameLastExtension(fileName); if (extension == ".vtk") { auto reader = vtkDataSetReader::New(); reader->SetFileName(fileName); return {reader}; } else if (extension == ".vtu") { auto reader = vtkXMLUnstructuredGridReader::New(); reader->SetFileName(fileName); return {reader}; } else if (extension == ".vtp") { auto reader = vtkXMLPolyDataReader::New(); reader->SetFileName(fileName); return {reader}; } else if (extension == ".obj") { auto reader = vtkOBJReader::New(); reader->SetFileName(fileName); return {reader}; } throw std::runtime_error("Invalid file extension."); } static void usage(char **argv) { std::cerr << "Usage: " << argv[0] << " analyze|fit tetmesh polydata [radiusScale relaxationIterCount]" << std::endl; } int main(int argc, char **argv) { if (argc < 4) { usage(argv); return EXIT_FAILURE; } bool fit = false; if (std::string(argv[1]) == "fit") { fit = true; } else if (std::string(argv[1]) != "analyze") { usage(argv); return EXIT_FAILURE; } if ((fit && argc > 6) || (!fit && argc != 4)) { usage(argv); return EXIT_FAILURE; } double radiusScale = 2.; int iterCount = 1; if(argc > 4) { radiusScale = std::stod(argv[4]); } if (argc > 5) { iterCount = std::stoi(argv[5]); } auto tetMeshReader = readerFromFileName(argv[2]); auto polyMeshReader = readerFromFileName(argv[3]); vtkNew removeExternalCellsFilter; removeExternalCellsFilter->SetInputConnection(0, tetMeshReader->GetOutputPort()); removeExternalCellsFilter->SetInputConnection(1, polyMeshReader->GetOutputPort()); vtkNew surfacePointsFilter; surfacePointsFilter->SetInputConnection( removeExternalCellsFilter->GetOutputPort()); vtkNew project; project->SetRadiusScale(radiusScale); project->SetInputConnection(0, surfacePointsFilter->GetOutputPort()); project->SetInputConnection(1, polyMeshReader->GetOutputPort()); vtkNew relaxationFilter; relaxationFilter->SetIterCount(iterCount); relaxationFilter->SetInputConnection(project->GetOutputPort()); vtkNew dihedralAnglesFilter; if (fit) { dihedralAnglesFilter->SetInputConnection( relaxationFilter->GetOutputPort()); } else { dihedralAnglesFilter->SetInputConnection( tetMeshReader->GetOutputPort()); } vtkNew geometryFilter; geometryFilter->SetInputConnection(dihedralAnglesFilter->GetOutputPort()); vtkNew maxDistanceFilter; maxDistanceFilter->SetInputConnection(0, geometryFilter->GetOutputPort()); maxDistanceFilter->SetInputConnection(1, polyMeshReader->GetOutputPort()); vtkNew writer; writer->SetInputConnection(dihedralAnglesFilter->GetOutputPort()); writer->SetDataModeToAscii(); writer->SetFileName("out.vtu"); writer->Write(); maxDistanceFilter->Update(); std::cerr << "Max distance: " << maxDistanceFilter->GetMaxDist() << "\n" << "Max dist point id: " << maxDistanceFilter->GetMaxId() << " in " << (maxDistanceFilter->GetMaxInput() == 0 ? "tetMesh\n" : "polyMesh\n") << "Average min angle: " << dihedralAnglesFilter->GetAverageMinDegrees() << "\n" << "Min min angle: " << dihedralAnglesFilter->GetMinMinDegrees() << "\n"; #ifdef USE_VIEWER /* Volume rendering properties */ vtkNew volumeMapper; volumeMapper->SetInputConnection(externalPointsFilter->GetOutputPort()); vtkNew volume; volume->SetMapper(volumeMapper); vtkNew transferFunction; transferFunction->AddPoint(-1, 0); transferFunction->AddPoint(1, 1); volume->GetProperty()->SetScalarOpacity(transferFunction); volume->GetProperty()->SetColor(transferFunction); /* Renderer */ vtkNew colors; std::array bkg{{26, 51, 102, 255}}; colors->SetColor("BkgColor", bkg.data()); vtkNew renderer; renderer->AddVolume(volume); renderer->SetBackground(colors->GetColor3d("BkgColor").GetData()); renderer->ResetCamera(); renderer->GetActiveCamera()->Zoom(1.5); vtkNew renderWindow; renderWindow->SetSize(300, 300); renderWindow->AddRenderer(renderer); renderWindow->SetWindowName("PFE"); vtkNew renderWindowInteractor; renderWindowInteractor->SetRenderWindow(renderWindow); renderWindow->Render(); renderWindowInteractor->Start(); #endif tetMeshReader->Delete(); polyMeshReader->Delete(); return EXIT_SUCCESS; }