From ebf1e0929ddff89ab014bfb8ea6d202722f218f8 Mon Sep 17 00:00:00 2001 From: gaozheng Date: Sun, 7 Jun 2026 21:51:21 +0800 Subject: [PATCH] =?UTF-8?q?feat(render):=20dd=5Fvoxel=20=E4=BD=93=E7=BB=98?= =?UTF-8?q?=E5=88=B6(IDW->vtkImageData->GPU=20RayCast)=20+=20=E4=BA=A4?= =?UTF-8?q?=E4=BA=92=E5=88=87=E7=89=87?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Co-Authored-By: Claude Opus 4.8 (1M context) --- src/app/CMakeLists.txt | 2 + src/app/main.cpp | 168 +++++++++++++++++++++++- src/data/repo/LocalSampleRepository.cpp | 10 ++ src/data/repo/LocalSampleRepository.hpp | 4 + src/render/CMakeLists.txt | 4 +- src/render/actors/VoxelActor.cpp | 99 ++++++++++++++ src/render/actors/VoxelActor.hpp | 26 ++++ tests/CMakeLists.txt | 4 +- tests/render/test_voxel_build.cpp | 56 ++++++++ 9 files changed, 368 insertions(+), 5 deletions(-) create mode 100644 src/render/actors/VoxelActor.cpp create mode 100644 src/render/actors/VoxelActor.hpp create mode 100644 tests/render/test_voxel_build.cpp diff --git a/src/app/CMakeLists.txt b/src/app/CMakeLists.txt index cd9cf42..571de88 100644 --- a/src/app/CMakeLists.txt +++ b/src/app/CMakeLists.txt @@ -4,7 +4,9 @@ find_package(VTK REQUIRED COMPONENTS GUISupportQt RenderingOpenGL2 + RenderingVolumeOpenGL2 InteractionStyle + InteractionWidgets FiltersGeometry FiltersModeling ) diff --git a/src/app/main.cpp b/src/app/main.cpp index 7a6fe55..760ddfe 100644 --- a/src/app/main.cpp +++ b/src/app/main.cpp @@ -35,10 +35,23 @@ #include "CameraPreset.hpp" #include "Scene.hpp" #include "actors/GridContourActor.hpp" +#include "actors/VoxelActor.hpp" + +#include "algo/IInterpolator.hpp" +#include "algo/IdwInterpolator.hpp" + +#include +#include +#include #include #include +#include +#include +#include #include +#include +#include namespace { @@ -71,6 +84,93 @@ std::string readPem(const std::string& path) return ss.str(); } +// dd_voxel 构建结果：体绘制 volume + 其底层 vtkImageData（供切片 widget）+ 网格信息。 +struct VoxelBuildResult { + vtkSmartPointer volume; + vtkSmartPointer image; + int nx = 0, ny = 0, nz = 0; + double vmin = 0.0, vmax = 0.0; +}; + +// dd_voxel：合并两条交叉剖面散点 → 局部化 → IDW 体素 → GPU 体绘制。 +// 与 tools/validate_voxel.py 同逻辑（power=2, maxDist≈4, dxy=1m, dz=0.5m）。 +VoxelBuildResult buildVoxelFromScatters( + const std::vector& scatters, + const geopro::core::ColorScale& cs) +{ + using geopro::core::GridSpec; + using geopro::core::IdwInterpolator; + using geopro::core::PointSet; + + constexpr double kDxy = 1.0; // 水平步长(米) + constexpr double kDz = 0.5; // 垂向步长(米) + constexpr double kPower = 2.0; // IDW 幂 + constexpr double kMaxDist = 4.0; // 超距留空(NaN) + + // 合并两剖面点为一组 (X=projX, Y=projY, Z=z)。 + std::vector X, Y, Z, V; + for (const auto& s : scatters) { + const size_t n = s.v.size(); + for (size_t p = 0; p < n; ++p) { + // projX/projY 为 GIS 平面坐标；z 为高程；v 为值。容错缺失字段。 + if (p < s.projX.size() && p < s.projY.size() && p < s.z.size()) { + X.push_back(s.projX[p]); + Y.push_back(s.projY[p]); + Z.push_back(s.z[p]); + V.push_back(s.v[p]); + } + } + } + if (V.empty()) return VoxelBuildResult{}; + + // 局部化：三轴各减最小值（规避大坐标；ox=oy=oz=0）。 + const double minX = *std::min_element(X.begin(), X.end()); + const double minY = *std::min_element(Y.begin(), Y.end()); + const double minZ = *std::min_element(Z.begin(), Z.end()); + PointSet pts; + double extX = 0.0, extY = 0.0, extZ = 0.0; + for (size_t p = 0; p < V.size(); ++p) { + const double lx = X[p] - minX, ly = Y[p] - minY, lz = Z[p] - minZ; + pts.x.push_back(lx); + pts.y.push_back(ly); + pts.z.push_back(lz); + pts.v.push_back(V[p]); + extX = std::max(extX, lx); + extY = std::max(extY, ly); + extZ = std::max(extZ, lz); + } + + GridSpec spec; + spec.ox = spec.oy = spec.oz = 0.0; + spec.dx = spec.dy = kDxy; + spec.dz = kDz; + spec.nx = static_cast(extX / kDxy) + 1; + spec.ny = static_cast(extY / kDxy) + 1; + spec.nz = static_cast(extZ / kDz) + 1; + spec.power = kPower; + spec.maxDist = kMaxDist; + + const auto vol = IdwInterpolator().interpolate(pts, spec); + + // 算 vmin/vmax（忽略 NaN）。 + double vmin = std::numeric_limits::max(); + double vmax = std::numeric_limits::lowest(); + for (double v : vol.data()) { + if (std::isnan(v)) continue; + vmin = std::min(vmin, v); + vmax = std::max(vmax, v); + } + if (vmin > vmax) { vmin = 0.0; vmax = 1.0; } // 全 NaN 兜底 + + VoxelBuildResult r; + r.image = vtkSmartPointer::New(); + r.volume = geopro::render::buildVoxel( + vol, cs, spec.ox, spec.oy, spec.oz, spec.dx, spec.dy, spec.dz, vmin, vmax, r.image); + r.nx = spec.nx; r.ny = spec.ny; r.nz = spec.nz; + r.vmin = vmin; r.vmax = vmax; + return r; +} + // 在给定 QMainWindow 上构建 M1 工作台：ADS 三栏 + 对象树 → 渲染联动 + 属性面板。 // repo 生命周期须覆盖到事件循环结束（由调用方保证）。 // 相机模式：默认二维俯视。 @@ -118,18 +218,30 @@ void buildWorkbench(QMainWindow& window, geopro::data::LocalSampleRepository& re cameraGroup->addAction(act2D); cameraGroup->addAction(act3D); act2D->setChecked(true); // 默认二维 + viewToolBar->addSeparator(); + auto* actVoxel = viewToolBar->addAction(QStringLiteral("dd_voxel")); centerLayout->addWidget(viewToolBar); centerLayout->addWidget(vtkWidget, 1); + // 交互切片 widget：dd_voxel 时对体素 vtkImageData 加可拖切面（dd_slice）。 + // 需 interactor（render window 已设），切回非体素视图时 Off。用 shared_ptr 跨 lambda 共享。 + auto sliceWidget = std::make_shared>(); + + auto hideSlice = [sliceWidget]() { + if (sliceWidget->Get()) (*sliceWidget)->Off(); + }; + QObject::connect(act2D, &QAction::triggered, vtkWidget, - [cameraMode, rendererPtr, renderWindowPtr]() { + [cameraMode, rendererPtr, renderWindowPtr, hideSlice]() { *cameraMode = CameraMode::Top2D; + hideSlice(); // 切回剖面视图：关闭体素切片 geopro::render::applyTop2D(rendererPtr); renderWindowPtr->Render(); }); QObject::connect(act3D, &QAction::triggered, vtkWidget, - [cameraMode, rendererPtr, renderWindowPtr]() { + [cameraMode, rendererPtr, renderWindowPtr, hideSlice]() { *cameraMode = CameraMode::Free3D; + hideSlice(); geopro::render::applyFree3D(rendererPtr); renderWindowPtr->Render(); }); @@ -182,6 +294,58 @@ void buildWorkbench(QMainWindow& window, geopro::data::LocalSampleRepository& re QObject::connect(tree, &QTreeWidget::itemClicked, tree, [renderDataset](QTreeWidgetItem* it, int) { renderDataset(it); }); + // 当前体素 vtkImageData（dd_slice 切面 widget 需复用，保活到下次重建）。 + auto voxelImage = std::make_shared>(); + + // dd_voxel：读两条交叉剖面 → IDW 体素 → GPU 体绘制 + 默认 3D 视角 + 可拖切片。 + QObject::connect(actVoxel, &QAction::triggered, vtkWidget, + [&repo, scene, rendererPtr, renderWindowPtr, cameraMode, act3D, propLabel, + sliceWidget, voxelImage]() { + const auto scatters = repo.loadVoxelScatters(); + const auto cs = repo.loadColorScale("grid1"); + auto result = buildVoxelFromScatters(scatters, cs); + if (!result.volume) { + propLabel->setText(QStringLiteral("dd_voxel: 无可用散点数据")); + return; + } + + // 关闭旧切片（避免引用已释放的旧 image）。 + if (sliceWidget->Get()) (*sliceWidget)->Off(); + + scene->clear(); + rendererPtr->AddViewProp(result.volume); + *voxelImage = result.image; + + // 体素默认 3D 视角（同步工具条勾选态）。 + *cameraMode = CameraMode::Free3D; + act3D->setChecked(true); + geopro::render::applyFree3D(rendererPtr); + rendererPtr->ResetCamera(); + renderWindowPtr->Render(); // 先渲一帧确保 interactor 就绪 + + // 交互切片 widget：对体素 vtkImageData 加 X 轴可拖切面（dd_slice）。 + vtkRenderWindowInteractor* interactor = renderWindowPtr->GetInteractor(); + if (interactor) { + if (!sliceWidget->Get()) + *sliceWidget = vtkSmartPointer::New(); + (*sliceWidget)->SetInteractor(interactor); + (*sliceWidget)->SetInputData(result.image); + (*sliceWidget)->SetPlaneOrientationToXAxes(); + (*sliceWidget)->SetSliceIndex(result.nx / 2); + (*sliceWidget)->On(); + } + + renderWindowPtr->Render(); + propLabel->setText( + QStringLiteral("dd_voxel 体素\n网格: %1 x %2 x %3\nvmin / vmax: %4 / %5\n" + "(拖动切面查看 dd_slice)") + .arg(result.nx) + .arg(result.ny) + .arg(result.nz) + .arg(result.vmin) + .arg(result.vmax)); + }); + // 默认渲染第一个 DS，让窗口一打开就有图。 if (auto* first = tree->topLevelItemCount() > 0 ? tree->topLevelItem(0) : nullptr) { // 递归找到首个带 dsId 的项。 diff --git a/src/data/repo/LocalSampleRepository.cpp b/src/data/repo/LocalSampleRepository.cpp index 6dde8ca..e684b2a 100644 --- a/src/data/repo/LocalSampleRepository.cpp +++ b/src/data/repo/LocalSampleRepository.cpp @@ -20,6 +20,9 @@ constexpr const char* kDsId = "grid1"; constexpr const char* kGridFile = u8"剖面网格数据1.txt"; constexpr const char* kColorScaleFile = u8"剖面网格数据的色阶数据1.txt"; constexpr const char* kScatterFile = u8"剖面原数据1.txt"; +// dd_voxel：两条交叉剖面散点（约 78° 相交），合并为体素插值输入。 +constexpr const char* kVoxelScatterFile1 = u8"剖面原数据1.txt"; +constexpr const char* kVoxelScatterFile2 = u8"剖面原数据2.txt"; constexpr const char* kAnomalyFile = u8"剖面网格数据1——对应的异常圈定数据.txt"; // 校验 dsId，未知则抛错（输入边界验证）。 @@ -89,4 +92,11 @@ std::vector LocalSampleRepository::loadAnomalies(const std::string& dsI return parseAnomalies(readFile(kAnomalyFile)); } +std::vector LocalSampleRepository::loadVoxelScatters() { + std::vector out; + out.push_back(parseScatter(readFile(kVoxelScatterFile1))); + out.push_back(parseScatter(readFile(kVoxelScatterFile2))); + return out; +} + } // namespace geopro::data diff --git a/src/data/repo/LocalSampleRepository.hpp b/src/data/repo/LocalSampleRepository.hpp index 2cc7b08..026ddc8 100644 --- a/src/data/repo/LocalSampleRepository.hpp +++ b/src/data/repo/LocalSampleRepository.hpp @@ -17,6 +17,10 @@ public: geopro::core::ColorScale loadColorScale(const std::string& dsId) override; std::vector loadAnomalies(const std::string& dsId) override; + // dd_voxel 输入：读两条交叉剖面散点（剖面原数据1.txt + 剖面原数据2.txt）， + // 返回两者，供体素插值合并。具体类专有方法（不进 IDatasetRepository 接口）。 + std::vector loadVoxelScatters(); + private: // 用 QFile（路径 = fromUtf8(dir)+fromUtf8(name)）读全文件并转 UTF-8 std::string。 // 读失败抛 std::runtime_error。 diff --git a/src/render/CMakeLists.txt b/src/render/CMakeLists.txt index dffef63..4662654 100644 --- a/src/render/CMakeLists.txt +++ b/src/render/CMakeLists.txt @@ -1,6 +1,6 @@ -find_package(VTK REQUIRED COMPONENTS CommonCore CommonDataModel FiltersGeometry FiltersModeling RenderingOpenGL2 InteractionStyle) +find_package(VTK REQUIRED COMPONENTS CommonCore CommonDataModel FiltersGeometry FiltersModeling RenderingOpenGL2 RenderingVolumeOpenGL2 InteractionStyle InteractionWidgets) add_library(geopro_render STATIC - Scene.cpp ColorLutBuilder.cpp CameraPreset.cpp actors/GridContourActor.cpp) + Scene.cpp ColorLutBuilder.cpp CameraPreset.cpp actors/GridContourActor.cpp actors/VoxelActor.cpp) target_include_directories(geopro_render PUBLIC ${CMAKE_CURRENT_SOURCE_DIR}) target_link_libraries(geopro_render PUBLIC geopro_core ${VTK_LIBRARIES}) target_compile_features(geopro_render PUBLIC cxx_std_17) diff --git a/src/render/actors/VoxelActor.cpp b/src/render/actors/VoxelActor.cpp new file mode 100644 index 0000000..63bf16b --- /dev/null +++ b/src/render/actors/VoxelActor.cpp @@ -0,0 +1,99 @@ +#include "actors/VoxelActor.hpp" + +#include +#include + +#include +#include +#include +#include +#include +#include +#include + +namespace geopro::render { + +namespace { + +// 颜色/不透明度传递函数采样级数。 +constexpr int kTransferSamples = 64; +// 体绘制最大不透明度（[vmin,vmax] 线性 0→kMaxOpacity）。 +constexpr double kMaxOpacity = 0.15; + +// NaN/留空格的哨兵：落在 [vmin,vmax] 之外，传递函数把它映射为完全透明。 +double sentinel(double vmin) { return vmin - 1.0; } + +} // namespace + +vtkSmartPointer buildVoxel(const geopro::core::ScalarVolume& vol, + const geopro::core::ColorScale& cs, + double ox, double oy, double oz, + double dx, double dy, double dz, + double vmin, double vmax, + vtkSmartPointer& outImage) +{ + const int nx = vol.nx(), ny = vol.ny(), nz = vol.nz(); + + // vmin/vmax 退化兜底，避免传递函数区间为零。 + if (vmin >= vmax) vmax = vmin + 1.0; + const double blank = sentinel(vmin); + + auto img = vtkSmartPointer::New(); + img->SetDimensions(nx, ny, nz); + img->SetOrigin(ox, oy, oz); + img->SetSpacing(dx, dy, dz); + + vtkNew sc; + sc->SetName("v"); + sc->SetNumberOfTuples(static_cast(nx) * ny * nz); + // 点序 i 最快、j 次之、k 最慢（匹配 vtkImageData 与 ScalarVolume::idx）。 + for (int k = 0; k < nz; ++k) + for (int j = 0; j < ny; ++j) + for (int i = 0; i < nx; ++i) { + const double v = vol.at(i, j, k); + const vtkIdType id = (static_cast(k) * ny + j) * nx + i; + sc->SetValue(id, std::isnan(v) ? blank : v); // NaN → 哨兵 + } + img->GetPointData()->SetScalars(sc); + outImage = img; + + // 颜色传递函数：在 [vmin,vmax] 按 ColorScale 采样若干 RGB 点。 + vtkNew color; + for (int t = 0; t < kTransferSamples; ++t) { + const double val = vmin + (vmax - vmin) * t / (kTransferSamples - 1); + const auto c = cs.colorAt(val); + color->AddRGBPoint(val, c.r / 255.0, c.g / 255.0, c.b / 255.0); + } + + // 不透明度传递函数：哨兵 → 0（透明）；[vmin,vmax] 线性递增到 kMaxOpacity。 + vtkNew opacity; + opacity->AddPoint(blank, 0.0); + opacity->AddPoint(vmin, 0.0); + opacity->AddPoint(vmax, kMaxOpacity); + + vtkNew mapper; + mapper->SetInputData(img); + + vtkNew prop; + prop->SetColor(color); + prop->SetScalarOpacity(opacity); + prop->SetInterpolationTypeToLinear(); + prop->ShadeOff(); + + auto volume = vtkSmartPointer::New(); + volume->SetMapper(mapper); + volume->SetProperty(prop); + return volume; +} + +vtkSmartPointer buildVoxel(const geopro::core::ScalarVolume& vol, + const geopro::core::ColorScale& cs, + double ox, double oy, double oz, + double dx, double dy, double dz, + double vmin, double vmax) +{ + vtkSmartPointer ignored; + return buildVoxel(vol, cs, ox, oy, oz, dx, dy, dz, vmin, vmax, ignored); +} + +} // namespace geopro::render diff --git a/src/render/actors/VoxelActor.hpp b/src/render/actors/VoxelActor.hpp new file mode 100644 index 0000000..4a840b9 --- /dev/null +++ b/src/render/actors/VoxelActor.hpp @@ -0,0 +1,26 @@ +#pragma once +#include +#include +#include +#include "model/Field.hpp" +#include "model/ColorScale.hpp" +namespace geopro::render { + +// 把 core 规则标量体（IDW 输出，含 NaN 留空）转 vtkImageData，再建 GPU 光线投射体绘制。 +// 颜色按 ColorScale 在 [vmin,vmax] 采样；NaN/留空格 → 不透明度 0（透明）。 +// 返回 vtkVolume（由调用方加入 renderer）。 +vtkSmartPointer buildVoxel(const geopro::core::ScalarVolume& vol, + const geopro::core::ColorScale& cs, + double ox, double oy, double oz, + double dx, double dy, double dz, + double vmin, double vmax); + +// 同上，但额外通过 out 参数暴露内部 vtkImageData，供 main 建交互切片 widget。 +vtkSmartPointer buildVoxel(const geopro::core::ScalarVolume& vol, + const geopro::core::ColorScale& cs, + double ox, double oy, double oz, + double dx, double dy, double dz, + double vmin, double vmax, + vtkSmartPointer& outImage); + +} // namespace geopro::render diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 703dcdf..bd16b57 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -57,8 +57,10 @@ endif() # render 层：ColorLutBuilder（core ColorScale -> vtkLookupTable）。 # 需 vtkLookupTable（VTK::CommonCore）；geopro_render 已 PUBLIC 传递其余 VTK 组件。 -find_package(VTK REQUIRED COMPONENTS CommonCore) +find_package(VTK REQUIRED COMPONENTS CommonCore CommonDataModel) target_sources(geopro_tests PRIVATE render/test_color_lut.cpp) +# dd_voxel：buildVoxel(ScalarVolume->vtkImageData->GPU 体绘制) 构建不崩 + dims 正确。 +target_sources(geopro_tests PRIVATE render/test_voxel_build.cpp) target_link_libraries(geopro_tests PRIVATE geopro_render ${VTK_LIBRARIES}) vtk_module_autoinit(TARGETS geopro_tests MODULES ${VTK_LIBRARIES}) diff --git a/tests/render/test_voxel_build.cpp b/tests/render/test_voxel_build.cpp new file mode 100644 index 0000000..53abbbb --- /dev/null +++ b/tests/render/test_voxel_build.cpp @@ -0,0 +1,56 @@ +#include + +#include +#include + +#include + +#include "actors/VoxelActor.hpp" +#include "model/ColorScale.hpp" +#include "model/Field.hpp" + +using namespace geopro::core; + +// buildVoxel 用小 ScalarVolume 构建：不崩、返回非空 volume、image dims/spacing/origin 正确。 +TEST(VoxelBuild, BuildsVolumeFromSmallScalarVolume) { + // 2x3x4 小体，填一些值 + 一个 NaN（验证哨兵透明分支）。 + ScalarVolume vol(2, 3, 4); + double n = 0.0; + for (int k = 0; k < 4; ++k) + for (int j = 0; j < 3; ++j) + for (int i = 0; i < 2; ++i) vol.at(i, j, k) = n++; + vol.at(1, 1, 1) = std::numeric_limits::quiet_NaN(); + + ColorScale cs; + cs.addStop(0.0, Rgba{0, 0, 255, 255}); + cs.addStop(24.0, Rgba{255, 0, 0, 255}); + + vtkSmartPointer image; + auto volume = geopro::render::buildVoxel( + vol, cs, /*ox*/ 1.0, /*oy*/ 2.0, /*oz*/ 3.0, + /*dx*/ 1.0, /*dy*/ 1.0, /*dz*/ 0.5, /*vmin*/ 0.0, /*vmax*/ 24.0, image); + + ASSERT_NE(volume.GetPointer(), nullptr); + ASSERT_NE(image.GetPointer(), nullptr); + + int dims[3]; + image->GetDimensions(dims); + EXPECT_EQ(dims[0], 2); + EXPECT_EQ(dims[1], 3); + EXPECT_EQ(dims[2], 4); + + double sp[3]; + image->GetSpacing(sp); + EXPECT_DOUBLE_EQ(sp[0], 1.0); + EXPECT_DOUBLE_EQ(sp[1], 1.0); + EXPECT_DOUBLE_EQ(sp[2], 0.5); + + double org[3]; + image->GetOrigin(org); + EXPECT_DOUBLE_EQ(org[0], 1.0); + EXPECT_DOUBLE_EQ(org[1], 2.0); + EXPECT_DOUBLE_EQ(org[2], 3.0); + + // NaN 格应被替换为哨兵（vmin-1 = -1），不应残留 NaN。 + EXPECT_FALSE(std::isnan(image->GetScalarComponentAsDouble(1, 1, 1, 0))); +}