feat/vtk-3d-view #7

Merged
gaozheng merged 301 commits from feat/vtk-3d-view into main 2026-06-27 18:43:52 +08:00
1 changed files with 166 additions and 22 deletions
Showing only changes of commit 7a294e13ab - Show all commits

View File

@ -1280,6 +1280,18 @@ struct GradStats {
};
GradStats sampleGradientMagnitude(vtkImageData* img);
// 标量分位标定P3 修可见性核心):扫该体实际体素值分布,取 2%/98% 分位作色阶/
// 不透明度的物理端点,裁掉离群。处理后体值多集中在 ±窄带、少量离群到 ±9000若按
// 全量化域(meta.vminPhys/vmaxPhys=±9249)映射 → 窄带信号落近透明区 → 整体近黑。
// 返回物理单位的 {lo, hi}(已按 quant.toPhys 反算)。前置声明,实现在 polish 段。
struct ScalarPercentiles {
double lo = 0.0, hi = 0.0; // 物理单位2% / 98% 分位)
std::size_t samples = 0;
};
ScalarPercentiles sampleScalarPercentiles(vtkImageData* img,
const geopro::core::Quant& q,
double pLo, double pHi);
geopro::core::ColorScale makeSeismicColorScale(double vmin, double vmax) {
geopro::core::ColorScale cs;
const double span = (vmax > vmin) ? (vmax - vmin) : 1.0;
@ -3059,6 +3071,10 @@ struct ViewState {
vtkRenderWindow* rw = nullptr;
double exagg = 8.0;
int lastLevel = -1;
std::string dir; // store 目录:首帧直读 level0 局部段(类 gallery),绕开 LOD 选粗层
// 预建的首帧高清段(level0 沿线中段)cmdView 已为分位标定建好,直接复用喂 mapper,
// 避免在 viewSetupDefaultFrame 内重复读盘。空则该函数再按 dir 直读。
vtkSmartPointer<vtkImageData> seedSegImg;
// 持有当前高清单图引用避免被释放mapper 仅持裸指针)。
vtkSmartPointer<vtkImageData> currentImg;
// 回调防重入:回调内部会 Render(),若 Render 又触发观察者回调会无限递归。
@ -3200,41 +3216,83 @@ std::size_t viewSetupDefaultFrame(ViewState* st, vtkRenderer* ren) {
const int totBricksX = (m.nx + brick - 1) / brick;
const int localBx = std::min(kViewDefaultLocalBricks, totBricksX);
const int bx0 = std::max(0, totBricksX / 2 - localBx / 2); // 沿线中段
// 该局部段世界 X 范围level0
(void)bx0;
(void)localBx;
(void)totBricksX;
// 体exagg 后)世界尺寸与中心 + 包围球半径。
const double wx = std::max(1.0, m.nx * m.spacing[0]);
// 该局部段世界 X 范围level0brick 列 [bx0, bx0+localBx))。把首帧相机框到这一段
// (而非整卷)是 P3 修复 #2 的关键:相机近观局部段 → C1 selectLod 选 level0 局部子区,
// C2 重组该段单图,framing 该段 → 14×796 截面 + 沿线一段充满视野(类厚 B-scan)。
// 框整卷则 selectLod 选最粗层(整条 45305 细带)、看着空白。
const double segX0 = m.origin[0] + bx0 * brick * m.spacing[0];
const double segX1 =
m.origin[0] + std::min(m.nx, (bx0 + localBx) * brick) * m.spacing[0];
// 段exagg 后)世界尺寸与中心 + 包围球半径。X 取该段宽(非整卷),Y/Z 全幅(薄轴)。
const double wx = std::max(1.0, segX1 - segX0);
const double wy = std::max(1.0, m.ny * m.spacing[1] * st->exagg);
const double wz = std::max(1.0, m.nz * m.spacing[2] * st->exagg);
const double cx = m.origin[0] + 0.5 * wx;
const double cx = 0.5 * (segX0 + segX1);
const double cy = m.origin[1] + 0.5 * wy;
const double cz = m.origin[2] + 0.5 * wz;
const double radius = 0.5 * std::sqrt(wx * wx + wy * wy + wz * wz);
// 相机从 +X 看体中心,距离 = 半径 / tan(半视角) × 余量确保整体落入视锥C1
// selectLod 不会判 empty由视距/分辨率自然选层;近观靠交互再拉近切细层。
// 相机从 +Y 看段中心(看进【X-Z 宽面=B-scan 墙】),距离 = 半径/tan(半视角)×余量。
// 段几何 = X≈12.6m 宽 × Y≈1.5m 薄(跨通道) × Z 深(exagg 后高)。exagg 只夸张深度(Z),
// Y 仍真实极薄 → 若从 +X(沿线)看只见薄前缘、近空。改从 +Y 俯看宽 X-Z 面: GPR 水平
// 分层沿 X 铺开、随 Z 叠层,这一面才读得出内部结构。整条概览靠用户滚轮拉远(细带几何必然)。
st->cam = ren->GetActiveCamera();
const double fovY = st->cam->GetViewAngle();
const double halfAngle = 0.5 * fovY * 3.14159265358979 / 180.0;
const double tanH = std::max(1e-3, std::tan(halfAngle));
const double dist = radius / tanH * 1.4; // 1.4:留余量含 aspect/边缘
st->cam->SetFocalPoint(cx, cy, cz);
st->cam->SetPosition(cx + dist, cy, cz);
st->cam->SetViewUp(0, 0, 1);
st->cam->SetPosition(cx, cy + dist, cz); // +Y 视点 → 正对 X-Z 宽面
st->cam->SetViewUp(0, 0, 1); // Z 朝上(深度向下)
ren->ResetCameraClippingRange();
// 源选层选区 + 重组单图喂 mapper。初始化场景需保证拿到首图 → 阻塞轮询到就绪。
const std::size_t blocks = viewRefreshBlocking(st);
// 首帧高清段直读(P3 修复 #2 核心):异步 LOD 源在「框一段」的视距下仍会选最粗层
// (整条 45305 细带 → 看着近黑),不可取。改为【直接从 store 读 level0 沿线中段子体】
// (与 gallery 的 buildLocalLevel0Image 同一直读路径,非 LOD 算法),喂高清 mapper —
// 保证首帧就是「全分辨率一段」的清晰块体。后续交互仍由异步源接管(用户拉远/拖动按
// 视距正常选层)。退化(读不到段)再回退异步阻塞刷新。
std::size_t blocks = 0;
{
vtkSmartPointer<vtkImageData> locImg = st->seedSegImg; // cmdView 预建,优先复用
if (locImg == nullptr && !st->dir.empty()) { // 退化:按 dir 直读
geopro::data::ChunkedVolumeStore store(st->dir);
locImg = buildLocalLevel0Image(store, m, bx0, localBx);
}
if (locImg != nullptr) {
st->currentImg = locImg;
st->lastLevel = 0;
st->mapper->SetInputData(locImg);
st->mapper->Update();
viewSyncBaseCropping(st); // 底图按该段挖空,无缝叠加
blocks = 1;
}
}
if (blocks == 0) blocks = viewRefreshBlocking(st); // 退化:回退异步源
// 框住局部段:用无参 ResetCamera按 actor 的【已 SetScale(1,exagg,exagg)】缩放
// 后包围盒框),相机角度沿用能看出结构的 Elevation/Azimuth再 Zoom 拉近填满画面。
ren->ResetCamera();
// 框住【局部段】:无参 ResetCamera 会按场景全部 actor(含常驻整卷底图,45305 长)的
// 包围盒框 → 整条细带、截面填不满 → 看着空白。改为只框高清段(currentImg=level0 沿
// 线中段子体)的包围盒,使 14×796 截面+沿线一段充满视野(类厚 B-scan);整条概览靠用户
// 滚轮拉远(细带是 1:34 几何必然,非 bug)。Z 轴按 actor 的 SetScale(1,1,exagg) 同步夸张。
if (st->currentImg != nullptr) {
double b[6];
st->currentImg->GetBounds(b); // 高清段模型坐标盒(含绝对 X origin)
b[4] *= st->exagg; // Z 下界随深度夸张
b[5] *= st->exagg; // Z 上界随深度夸张
ren->ResetCamera(b); // 只框该段 → 段充满画面
} else {
ren->ResetCamera(); // 退化(无高清段):回退全场景框
}
// 取景角度:默认相机已置 +Y 正对 X-Z 宽面。var4 的 El45/Az30 是为 gallery 的 +X /
// Y 也夸张几何调的,套到这里(+Y 视、Y 真实极薄)会让宽面强烈斜退、大片黑。改用为本
// +Y B-scan 几何调的小角度(El/Az 各 ~15-18°)留 3D 立体感而宽面仍基本正对,Zoom 拉足
// 让 X-Z 面填满。仅默认/交互/preview 取景,不动 var4 gallery 参数。
constexpr double kDefaultFrameElevation = 18.0; // 轻俯,见顶面薄边显层叠
constexpr double kDefaultFrameAzimuth = 15.0; // 轻偏,宽面仍基本正对
constexpr double kDefaultFrameZoom = 1.1; // ResetCamera 已贴合该段,只略收边距
st->cam = ren->GetActiveCamera();
st->cam->Elevation(kViewDefaultVariant.elevation); // var4 取景El18
st->cam->Azimuth(kViewDefaultVariant.azimuth); // var4 取景Az22
st->cam->Zoom(kViewDefaultVariant.zoom); // var4 取景Zoom2.0 填满画面
st->cam->Elevation(kDefaultFrameElevation);
st->cam->Azimuth(kDefaultFrameAzimuth);
st->cam->Zoom(kDefaultFrameZoom);
ren->ResetCameraClippingRange();
return blocks;
}
@ -3245,8 +3303,6 @@ int runGalleryVariant(const std::string& dir, const GalleryVariant& v,
const int winW = 1280, winH = 800;
geopro::data::ChunkedVolumeStore store(dir);
const geopro::data::StoreMeta& m = store.meta();
const double vmin = m.vminPhys, vmax = m.vmaxPhys;
const geopro::core::ColorScale cs = pickColor(v.color, vmin, vmax);
auto rw = makeOffscreenWindow(winW, winH);
vtkNew<vtkRenderer> ren;
@ -3261,6 +3317,21 @@ int runGalleryVariant(const std::string& dir, const GalleryVariant& v,
vtkSmartPointer<vtkImageData> locImg =
buildLocalLevel0Image(store, m, bx0, localBx);
// P3 传函分位标定:色阶/不透明度端点按该局部段【实际值 2%/98% 分位】裁离群。
double vmin = m.vminPhys, vmax = m.vmaxPhys;
{
const ScalarPercentiles pc =
sampleScalarPercentiles(locImg.Get(), m.quant, 0.02, 0.98);
if (pc.samples > 0) {
vmin = pc.lo;
vmax = pc.hi;
std::cout << "[gallery " << v.name << "] 传函分位标定(样本 " << pc.samples
<< "): 2%=" << vmin << " 98%=" << vmax << " (全域 ["
<< m.vminPhys << ", " << m.vmaxPhys << "])\n";
}
}
const geopro::core::ColorScale cs = pickColor(v.color, vmin, vmax);
GradStats gs;
if (v.useGradientOpacity) {
gs = sampleGradientMagnitude(locImg.Get());
@ -3455,7 +3526,38 @@ int cmdView(int argc, char** argv) {
source.setViewportHeight(winH);
source.setAspect(static_cast<double>(winW) / winH);
const auto& m = source.meta();
const double vmin = m.vminPhys, vmax = m.vmaxPhys;
// P3 首帧高清段(沿线中段 level0):直读 store 建该段单图——既作传函分位标定的基准
// (该段实际值,对比比整卷底图更punchy),又供 viewSetupDefaultFrame 直接喂高清 mapper
// (绕开异步 LOD 在「框一段」视距下仍选最粗层的问题)。同 gallery 的 buildLocalLevel0Image
// 直读路径,非 LOD 算法。读不到则回退底图/全域。
vtkSmartPointer<vtkImageData> segImg;
{
geopro::data::ChunkedVolumeStore store(dir);
const int totBx = store.bricksX(0);
const int localBx = std::min(kViewDefaultLocalBricks, totBx);
const int bx0 = std::max(0, totBx / 2 - localBx / 2);
segImg = buildLocalLevel0Image(store, m, bx0, localBx);
}
// P3 传函分位标定:色阶/不透明度端点按【实际值 2%/98% 分位】(物理域),裁离群。
// 处理后体值多集中 ±窄带、少量离群 ±9000按 meta 全域(±9249)映射 → 窄带近透明 → 全黑。
// 基准用常驻底图(整卷代表):band 较窄 → 窄带信号映到更饱和色、可见度更高;且代表整线
// 典型信号、不受单段局部高能离群影响。底图缺失再回退首帧段。退化(无样本)回退全域。
double vmin = m.vminPhys, vmax = m.vmaxPhys;
vtkImageData* pcBasis =
source.baseImage() != nullptr ? source.baseImage() : segImg.Get();
if (pcBasis != nullptr) {
const ScalarPercentiles pc =
sampleScalarPercentiles(pcBasis, m.quant, 0.02, 0.98);
if (pc.samples > 0) {
vmin = pc.lo;
vmax = pc.hi;
std::cout << "[view] 传函分位标定("
<< (source.baseImage() != nullptr ? "底图" : "局部段")
<< ",样本 " << pc.samples << "): 2%=" << vmin << " 98%=" << vmax
<< " (全域 [" << m.vminPhys << ", " << m.vmaxPhys
<< "] → 裁离群)\n";
}
}
// 配色/不透明度包络取自 var4seismic + V 形实体包络(floor/mid + opacity 作峰值)。
const geopro::core::ColorScale cs = pickColor(dv.color, vmin, vmax);
// C4默认变体(var4)开了梯度不透明度 → 从常驻底图实测梯度分布标定阈值。底图恒非空。
@ -3541,6 +3643,8 @@ int cmdView(int argc, char** argv) {
st.fpsText = fpsText.Get();
st.rw = rw.Get();
st.exagg = exagg;
st.dir = dir; // 供首帧直读 level0 局部段(绕开 LOD 选粗层)
st.seedSegImg = segImg; // cmdView 已为分位标定建好的首帧高清段,直接复用
// 相机初始定向(修复 1默认框「局部段」而非整卷。整线横截面 1:34框整卷
// 即便 exagg=8 也是一条隐形细带(看着空白);改为对准沿线中段一个 ~768 道窗口
@ -3905,6 +4009,46 @@ GradStats sampleGradientMagnitude(vtkImageData* img) {
return gs;
}
// 标量分位标定实现P3步长抽样该 VTK_SHORT 体的非空体素值,排序取 pLo/pHi 分位,
// 用 quant.toPhys 反算成物理端点。失败/退化(无样本或 lo>=hi)返回 samples=0由调用方
// 回退到 meta 全量化域。仿 sampleGradientMagnitude 的 stride/blank 处理,自适应该体值域。
ScalarPercentiles sampleScalarPercentiles(vtkImageData* img,
const geopro::core::Quant& q,
double pLo, double pHi) {
ScalarPercentiles out;
if (!img) return out;
int dims[3];
img->GetDimensions(dims);
const vtkIdType npts =
static_cast<vtkIdType>(dims[0]) * dims[1] * dims[2];
auto* arr = vtkShortArray::SafeDownCast(img->GetPointData()->GetScalars());
if (!arr || npts <= 0) return out;
const std::int16_t blank = geopro::core::ScalarVolumeI16::kBlank;
std::vector<std::int16_t> vals;
vals.reserve(static_cast<std::size_t>(npts) / 4 + 1);
// 步长抽样:~50 万样本足以代表分布(与梯度采样同档)。
const vtkIdType stride = std::max<vtkIdType>(1, npts / 500000);
for (vtkIdType id = 0; id < npts; id += stride) {
const std::int16_t v = arr->GetValue(id);
if (v == blank) continue;
vals.push_back(v);
}
if (vals.empty()) return out;
std::sort(vals.begin(), vals.end());
auto pick = [&](double p) {
const std::size_t idx = static_cast<std::size_t>(
std::min<double>(vals.size() - 1, p * (vals.size() - 1)));
return static_cast<double>(vals[idx]);
};
const double qLo = pick(pLo);
const double qHi = pick(pHi);
if (!(qLo < qHi)) return out; // 退化(常数体)→ 调用方回退全域
out.lo = q.toPhys(static_cast<std::int16_t>(std::lround(qLo)));
out.hi = q.toPhys(static_cast<std::int16_t>(std::lround(qHi)));
out.samples = vals.size();
return out;
}
// 标量不透明度V 形包络(与 makeSolidVolumeProperty 同思路floor/mid/max三图共用
// 保证唯一变量是梯度不透明度 / 光照。
void setPolishScalarOpacity(vtkVolumeProperty* prop, const geopro::core::Quant& q,