diff --git a/tsd/src/tsd/app/Core.cpp b/tsd/src/tsd/app/Core.cpp index a9912986..69d52e02 100644 --- a/tsd/src/tsd/app/Core.cpp +++ b/tsd/src/tsd/app/Core.cpp @@ -93,6 +93,8 @@ void Core::parseCommandLine(std::vector &args) importerType = tsd::io::ImporterType::USD2; else if (arg == "-vtp") importerType = tsd::io::ImporterType::VTP; + else if (arg == "-vtu") + importerType = tsd::io::ImporterType::VTU; else if (arg == "-xyzdp") importerType = tsd::io::ImporterType::XYZDP; else if (arg == "-volume") diff --git a/tsd/src/tsd/io/importers.hpp b/tsd/src/tsd/io/importers.hpp index 53d149ce..567e6b85 100644 --- a/tsd/src/tsd/io/importers.hpp +++ b/tsd/src/tsd/io/importers.hpp @@ -38,6 +38,7 @@ void import_TRK(Scene &scene, const char *filename, LayerNodeRef location = {}); void import_USD(Scene &scene, const char *filename, LayerNodeRef location = {}); void import_USD2(Scene &scene, const char *filename, LayerNodeRef location = {}); void import_VTP(Scene &scene, const char *filepath, LayerNodeRef location = {}); +void import_VTU(Scene &scene, const char *filepath, LayerNodeRef location); void import_XYZDP(Scene &scene, const char *filename, LayerNodeRef location = {}); // Spatial field importers // @@ -93,6 +94,7 @@ enum class ImporterType USD, USD2, VTP, + VTU, XYZDP, VOLUME, TSD, diff --git a/tsd/src/tsd/io/importers/import_VTU.cpp b/tsd/src/tsd/io/importers/import_VTU.cpp index d90baca7..50e90e8f 100644 --- a/tsd/src/tsd/io/importers/import_VTU.cpp +++ b/tsd/src/tsd/io/importers/import_VTU.cpp @@ -5,28 +5,66 @@ #define TSD_USE_VTK 1 #endif +#include "tsd/core/Logging.hpp" +#include "tsd/core/algorithms/computeScalarRange.hpp" #include "tsd/io/importers.hpp" #include "tsd/io/importers/detail/importer_common.hpp" -#include "tsd/core/Logging.hpp" #if TSD_USE_VTK // vtk #include +#include #include #include #include #include +#include +#include #include +#include #include #include #include #endif -// std -#include -#include namespace tsd::io { #if TSD_USE_VTK + +static bool isSurfaceCell(int type) +{ + return type == VTK_TRIANGLE || type == VTK_QUAD || type == VTK_POLYGON; +} + +static bool isVolumeCell(int type) +{ + return type == VTK_TETRA || type == VTK_VOXEL || type == VTK_HEXAHEDRON + || type == VTK_WEDGE || type == VTK_PYRAMID; +} + +static vtkSmartPointer loadVTUGrid(const char *filepath) +{ + auto reader = vtkSmartPointer::New(); + auto legacyReader = vtkSmartPointer::New(); + + vtkUnstructuredGrid *grid = nullptr; + if (reader->CanReadFile(filepath)) { + reader->SetFileName(filepath); + reader->Update(); + grid = reader->GetOutput(); + } else { + legacyReader->SetFileName(filepath); + legacyReader->Update(); + grid = legacyReader->GetOutput(); + } + + if (!grid) { + logError("[import_VTU] failed to load .vtu file '%s'", filepath); + return nullptr; + } + + return grid; +} + static ArrayRef makeFloatArray1D( Scene &scene, vtkDataArray *array, vtkIdType count) { @@ -46,101 +84,360 @@ static ArrayRef makeFloatArray1D( return arr; } -SpatialFieldRef import_VTU(Scene &scene, const char *filepath) +static bool isScalarArray(const ArrayRef &a) +{ + if (!a) + return false; + + const auto type = a->elementType(); + return !anari::isObject(type) && anari::componentsOf(type) == 1; +} + +static ArrayRef firstScalarArray(const std::vector &arrays) { - vtkUnstructuredGrid *grid{nullptr}; + for (const auto &a : arrays) { + if (isScalarArray(a)) + return a; + } - // Read .vtu file - vtkSmartPointer reader = - vtkSmartPointer::New(); + return {}; +} - vtkSmartPointer legacyReader = - vtkSmartPointer::New(); +// Build a triangle surface from surface cells in the grid +static void createSurfaceFromGrid(Scene &scene, + vtkUnstructuredGrid *grid, + const char *filepath, + LayerNodeRef location) +{ + auto filename = fileOf(filepath); - if (reader->CanReadFile(filepath)) { - reader->SetFileName(filepath); - reader->Update(); - grid = reader->GetOutput(); - } else { - // legacy reader has no CanReadFile - legacyReader->SetFileName(filepath); - legacyReader->Update(); - grid = legacyReader->GetOutput(); + // Extract surface cells into a vtkPolyData for triangulation + auto polyData = vtkSmartPointer::New(); + polyData->SetPoints(grid->GetPoints()); + + auto polys = vtkSmartPointer::New(); + std::vector surfaceCellIndices; + + for (vtkIdType i = 0; i < grid->GetNumberOfCells(); ++i) { + if (!isSurfaceCell(grid->GetCellType(i))) + continue; + surfaceCellIndices.push_back(i); + vtkCell *cell = grid->GetCell(i); + vtkIdType npts = cell->GetNumberOfPoints(); + std::vector pts(npts); + for (vtkIdType j = 0; j < npts; ++j) + pts[j] = cell->GetPointId(j); + polys->InsertNextCell(npts, pts.data()); } - if (!grid) { - logError("[import_VTU] failed to load .vtu file '%s'", filepath); - return {}; + polyData->SetPolys(polys); + polyData->GetPointData()->ShallowCopy(grid->GetPointData()); + + auto *srcCellData = grid->GetCellData(); + auto *dstCellData = polyData->GetCellData(); + for (int a = 0; a < srcCellData->GetNumberOfArrays(); ++a) { + vtkDataArray *srcArr = srcCellData->GetArray(a); + if (!srcArr) + continue; + auto filtered = vtkSmartPointer::Take(srcArr->NewInstance()); + filtered->SetName(srcArr->GetName()); + filtered->SetNumberOfComponents(srcArr->GetNumberOfComponents()); + filtered->SetNumberOfTuples( + static_cast(surfaceCellIndices.size())); + for (vtkIdType i = 0; i < (vtkIdType)surfaceCellIndices.size(); ++i) + filtered->SetTuple(i, surfaceCellIndices[i], srcArr); + dstCellData->AddArray(filtered); } - auto field = - scene.createObject(tokens::spatial_field::unstructured); - field->setName(fileOf(filepath).c_str()); + // Triangulate + auto triangleFilter = vtkSmartPointer::New(); + triangleFilter->SetInputData(polyData); + triangleFilter->Update(); - vtkIdType numPoints = grid->GetNumberOfPoints(); + vtkPolyData *triangleMesh = triangleFilter->GetOutput(); + + vtkPoints *points = triangleMesh->GetPoints(); + if (!points || points->GetNumberOfPoints() == 0) + return; + + vtkCellArray *triangles = triangleMesh->GetPolys(); + if (!triangles || triangles->GetNumberOfCells() == 0) + return; + + // Build ANARI representation + const vtkIdType numPoints = points->GetNumberOfPoints(); + auto vertexArray = scene.createArray(ANARI_FLOAT32_VEC3, numPoints); + auto *vertexPtr = vertexArray->mapAs(); + double p[3]; + for (vtkIdType i = 0; i < numPoints; ++i) { + points->GetPoint(i, p); + vertexPtr[i] = tsd::math::float3(static_cast(p[0]), + static_cast(p[1]), + static_cast(p[2])); + } + vertexArray->unmap(); + + auto idList = vtkSmartPointer::New(); + const vtkIdType numCells = triangles->GetNumberOfCells(); + std::vector triangleIndices; + triangleIndices.reserve(numCells * 3); + + triangles->InitTraversal(); + while (triangles->GetNextCell(idList)) { + if (idList->GetNumberOfIds() != 3) + continue; + for (int i = 0; i < 3; ++i) + triangleIndices.push_back(static_cast(idList->GetId(i))); + } + + auto indexArray = + scene.createArray(ANARI_UINT32_VEC3, triangleIndices.size() / 3); + indexArray->setData(triangleIndices.data()); + + // Vertex data + std::vector vertexDataArrays; + vtkPointData *ptData = triangleMesh->GetPointData(); + for (int i = 0; i < ptData->GetNumberOfArrays(); ++i) { + vtkDataArray *arr = ptData->GetArray(i); + if (!arr) + continue; + std::string name = arr->GetName() ? arr->GetName() : "unnamed_point_array"; + auto a = makeArray1DFromVTK(scene, arr, "[import_VTU]"); + a->setName(name.c_str()); + vertexDataArrays.push_back(a); + } + + // Face data + std::vector faceDataArrays; + vtkCellData *clData = triangleMesh->GetCellData(); + for (int i = 0; i < clData->GetNumberOfArrays(); ++i) { + vtkDataArray *arr = clData->GetArray(i); + if (!arr) + continue; + std::string name = arr->GetName() ? arr->GetName() : "unnamed_cell_array"; + auto a = makeArray1DFromVTK(scene, arr, "[import_VTU]"); + a->setName(name.c_str()); + faceDataArrays.push_back(a); + } + + // Assemble ANARI surfaces + auto mesh = + scene.createObject(tokens::geometry::triangle); + mesh->setName(("vtu_mesh | " + std::string(filename)).c_str()); + mesh->setParameterObject("vertex.position", *vertexArray); + mesh->setParameterObject("primitive.index", *indexArray); + + for (size_t i = 0; i < vertexDataArrays.size(); ++i) + mesh->setParameterObject( + "vertex.attribute" + std::to_string(i), *vertexDataArrays[i]); + + for (size_t i = 0; i < faceDataArrays.size(); ++i) + mesh->setParameterObject( + "primitive.attribute" + std::to_string(i), *faceDataArrays[i]); + + auto mat = scene.createObject( + tokens::material::physicallyBased); + mat->setName(("vtu_material | " + std::string(filename)).c_str()); + + if (auto colorArray = firstScalarArray(vertexDataArrays); colorArray) { + auto colorRange = tsd::core::computeScalarRange(*colorArray); + mat->setParameterObject( + "baseColor", *makeDefaultColorMapSampler(scene, colorRange)); + } else if (auto colorArray = firstScalarArray(faceDataArrays); colorArray) { + auto colorRange = tsd::core::computeScalarRange(*colorArray); + mat->setParameterObject( + "baseColor", *makeDefaultColorMapSampler(scene, colorRange)); + } else if (!vertexDataArrays.empty() || !faceDataArrays.empty()) { + logWarning( + "[import_VTU] no scalar point/cell arrays found for colormapping"); + } + + scene.insertChildObjectNode(location, + scene.createSurface( + ("vtu_surface | " + std::string(filename)).c_str(), mesh, mat)); +} + +// Build an unstructured spatial field from volume cells only +static SpatialFieldRef createFieldFromVolumeCells( + Scene &scene, vtkUnstructuredGrid *grid, const char *filepath) +{ vtkIdType numCells = grid->GetNumberOfCells(); - // --- Write vertex positions --- + // Collect volume cell indices + std::vector volumeCellIndices; + for (vtkIdType i = 0; i < numCells; ++i) { + if (isVolumeCell(grid->GetCellType(i))) + volumeCellIndices.push_back(i); + } + + if (volumeCellIndices.empty()) + return {}; + + vtkIdType numPoints = grid->GetNumberOfPoints(); + vtkIdType numVolumeCells = static_cast(volumeCellIndices.size()); + + auto field = scene.createObject( + tokens::spatial_field::unstructured); + field->setName(fileOf(filepath).c_str()); + + // Vertices auto vertexArray = scene.createArray(ANARI_FLOAT32_VEC3, numPoints); - auto *vertexData = vertexArray->mapAs(); + auto *vertexData = vertexArray->mapAs(); for (vtkIdType i = 0; i < numPoints; ++i) { double *pt = grid->GetPoint(i); - vertexData[i] = float3(pt[0], pt[1], pt[2]); + vertexData[i] = tsd::math::float3(pt[0], pt[1], pt[2]); } vertexArray->unmap(); field->setParameterObject("vertex.position", *vertexArray); - // --- Write cell indices (flattened connectivity) --- + // Cells std::vector connectivity; std::vector cellIndex; - for (vtkIdType i = 0; i < numCells; ++i) { - vtkCell *cell = grid->GetCell(i); + std::vector cellTypes; + + for (vtkIdType idx : volumeCellIndices) { + vtkCell *cell = grid->GetCell(idx); int n = cell->GetNumberOfPoints(); + int vtkType = grid->GetCellType(idx); + cellIndex.push_back(static_cast(connectivity.size())); - for (int j = 0; j < n; ++j) - connectivity.push_back(static_cast(cell->GetPointId(j))); + + if (vtkType == VTK_VOXEL) { + // Voxels need vertex reordering to become hexahedra: + // swap vertices 2<->3 and 6<->7 + std::vector pts(n); + for (int j = 0; j < n; ++j) + pts[j] = cell->GetPointId(j); + if (n >= 8) { + std::swap(pts[2], pts[3]); + std::swap(pts[6], pts[7]); + } + for (int j = 0; j < n; ++j) + connectivity.push_back(static_cast(pts[j])); + cellTypes.push_back(static_cast(VTK_HEXAHEDRON)); + } else { + for (int j = 0; j < n; ++j) + connectivity.push_back(static_cast(cell->GetPointId(j))); + cellTypes.push_back(static_cast(vtkType)); + } } + auto indexArray = scene.createArray(ANARI_UINT32, connectivity.size()); indexArray->setData(connectivity.data()); field->setParameterObject("index", *indexArray); + auto cellIndexArray = scene.createArray(ANARI_UINT32, cellIndex.size()); cellIndexArray->setData(cellIndex.data()); field->setParameterObject("cell.index", *cellIndexArray); - // --- Write cell types --- - auto cellTypesArray = scene.createArray(ANARI_UINT8, numCells); - auto *cellTypes = cellTypesArray->mapAs(); - for (vtkIdType i = 0; i < numCells; ++i) - cellTypes[i] = static_cast(grid->GetCellType(i)); - cellTypesArray->unmap(); + auto cellTypesArray = scene.createArray(ANARI_UINT8, cellTypes.size()); + cellTypesArray->setData(cellTypes.data()); field->setParameterObject("cell.type", *cellTypesArray); - // --- Write point data arrays --- + // Vertex data vtkPointData *pointData = grid->GetPointData(); - uint32_t numPointArrays = pointData->GetNumberOfArrays(); - for (uint32_t i = 0; i < std::min(1u, numPointArrays); ++i) { + for (int i = 0; i < std::min(1, pointData->GetNumberOfArrays()); ++i) { vtkDataArray *array = pointData->GetArray(i); + if (!array) + continue; auto a = makeFloatArray1D(scene, array, numPoints); field->setParameterObject("vertex.data", *a); } - // --- Write cell data arrays --- + // Cell data vtkCellData *cellData = grid->GetCellData(); - uint32_t numCellArrays = cellData->GetNumberOfArrays(); - for (uint32_t i = 0; i < std::min(1u, numCellArrays); ++i) { + for (int i = 0; i < std::min(1, cellData->GetNumberOfArrays()); ++i) { vtkDataArray *array = cellData->GetArray(i); + if (!array) + continue; auto a = makeFloatArray1D(scene, array, numCells); - field->setParameterObject("cell.data", *a); + + // Filter to volume cells only + auto filtered = scene.createArray(ANARI_FLOAT32, numVolumeCells); + auto *src = a->mapAs(); + auto *dst = filtered->mapAs(); + for (vtkIdType j = 0; j < numVolumeCells; ++j) + dst[j] = src[volumeCellIndices[j]]; + a->unmap(); + filtered->unmap(); + + field->setParameterObject("cell.data", *filtered); } return field; } + +// Full-scene importer: surfaces + volumes +void import_VTU(Scene &scene, const char *filepath, LayerNodeRef location) +{ + auto grid = loadVTUGrid(filepath); + if (!grid) + return; + + vtkIdType numCells = grid->GetNumberOfCells(); + bool hasSurface = false; + bool hasVolume = false; + int unsupportedCount = 0; + + for (vtkIdType i = 0; i < numCells; ++i) { + int type = grid->GetCellType(i); + if (isSurfaceCell(type)) + hasSurface = true; + else if (isVolumeCell(type)) + hasVolume = true; + else + unsupportedCount++; + } + + if (unsupportedCount > 0) { + logWarning("[import_VTU] %d cells with unsupported types skipped", + unsupportedCount); + } + + auto root = scene.insertChildTransformNode( + location ? location : scene.defaultLayer()->root(), + tsd::math::mat4(tsd::math::identity), + fileOf(filepath).c_str()); + + if (hasSurface) + createSurfaceFromGrid(scene, grid, filepath, root); + + if (hasVolume) { + auto field = createFieldFromVolumeCells(scene, grid, filepath); + if (field) { + tsd::math::float2 valueRange = field->computeValueRange(); + auto [inst, volume] = scene.insertNewChildObjectNode( + root, tokens::volume::transferFunction1D); + volume->setName(fileOf(filepath).c_str()); + volume->setParameterObject("value", *field); + volume->setParameter("valueRange", ANARI_FLOAT32_BOX1, &valueRange); + } + } +} + +// Spatial field-only importer (backward compatible, used by -volume) +SpatialFieldRef import_VTU(Scene &scene, const char *filepath) +{ + auto grid = loadVTUGrid(filepath); + if (!grid) + return {}; + + return createFieldFromVolumeCells(scene, grid, filepath); +} + #else + +void import_VTU(Scene &scene, const char *filepath, LayerNodeRef location) +{ + logError("[import_VTU] VTK not enabled in TSD build."); +} + SpatialFieldRef import_VTU(Scene &scene, const char *filepath) { logError("[import_VTU] VTK not enabled in TSD build."); return {}; } + #endif -} // namespace tsd +} // namespace tsd::io diff --git a/tsd/src/tsd/io/importers/import_file.cpp b/tsd/src/tsd/io/importers/import_file.cpp index 73748872..f0d396ca 100644 --- a/tsd/src/tsd/io/importers/import_file.cpp +++ b/tsd/src/tsd/io/importers/import_file.cpp @@ -86,6 +86,8 @@ void import_file(Scene &scene, tsd::io::import_USD2(scene, file.c_str(), root); } else if (f.first == ImporterType::VTP) tsd::io::import_VTP(scene, file.c_str(), root); + else if (f.first == ImporterType::VTU) + tsd::io::import_VTU(scene, file.c_str(), root); else if (f.first == ImporterType::XYZDP) tsd::io::import_XYZDP(scene, file.c_str(), root); else if (f.first == ImporterType::VOLUME)