feat(core): GPR 结构化建体 buildGprVolume(X/Z 落格 + Y 向 1D 线性插值 → int16 量化体)

- 新增 GprSurvey 规则化建体输入模型(放 core/model 保持 geopro_core 自洽,避免 core->io 反向依赖)
- buildGprVolume: X/Z 取最近道/采样落格,仅跨通道 Y 做 1D 线性插值,边界外不外推
- int16 量化用值域中点为 offset 对称铺满 ~64000 码位,两端留余量不撞 int16/kBlank
- 整型乘积索引走 size_t
This commit is contained in:
gaozheng 2026-06-23 10:45:06 +08:00
parent 9874af77ee
commit a9e8eb9d5c
6 changed files with 239 additions and 0 deletions

View File

@ -8,6 +8,7 @@ add_library(geopro_core STATIC
model/ColorScale.cpp model/ColorScale.cpp
algo/IdwInterpolator.cpp algo/IdwInterpolator.cpp
algo/VolumeBuilder.cpp algo/VolumeBuilder.cpp
algo/GprVolumeBuilder.cpp
) )
target_include_directories(geopro_core PUBLIC ${CMAKE_CURRENT_SOURCE_DIR}) target_include_directories(geopro_core PUBLIC ${CMAKE_CURRENT_SOURCE_DIR})

View File

@ -0,0 +1,112 @@
#include "algo/GprVolumeBuilder.hpp"
#include <cmath>
#include <limits>
namespace geopro::core {
namespace {
// 把世界坐标轴值映射到最近整数索引;落在 [0, n-1] 外返回 -1标记越界
int nearestIndex(double world, double origin, double step, int n) {
if (n <= 0 || step == 0.0) return -1;
long idx = std::lround((world - origin) / step);
if (idx < 0 || idx >= n) return -1;
return static_cast<int>(idx);
}
} // namespace
BuiltI16 buildGprVolume(const GprSurvey& s, const GridSpec& spec) {
BuiltI16 out;
out.origin = {spec.ox, spec.oy, spec.oz};
out.spacing = {spec.dx, spec.dy, spec.dz};
// 1. 物理值域(扫 values跳过 NaN
double vmin = std::numeric_limits<double>::infinity();
double vmax = -std::numeric_limits<double>::infinity();
for (double v : s.values) {
if (std::isnan(v)) continue;
if (v < vmin) vmin = v;
if (v > vmax) vmax = v;
}
if (!(vmin <= vmax)) { // 无有效值:退化为 [0,0]。
vmin = 0.0;
vmax = 0.0;
}
out.vminPhys = vmin;
out.vmaxPhys = vmax;
// 2. 量化scale=(vmax-vmin)/64000 把整个物理值域铺满 int16 的 ~64000 个码位
// -32000..+32000故 offset 取值域中点使量化对称——既用满 64000 码位,
// 又两端各留 ~700 余量不撞 int16 边界(±32767)与 kBlank(INT16_MIN)。
// vmax==vmin 时 scale=1。
out.quant.scale = (vmax > vmin) ? (vmax - vmin) / 64000.0 : 1.0;
out.quant.offset = 0.5 * (vmin + vmax);
// 3. 分配体(构造即填 0origin/spacing 已设。
out.vol = ScalarVolumeI16(spec.nx, spec.ny, spec.nz);
const int nChan = static_cast<int>(s.channelY.size());
// 4. 逐网格点落值。
for (int gk = 0; gk < spec.nz; ++gk) {
const double worldZ = spec.oz + gk * spec.dz;
const int sIdx = nearestIndex(worldZ, s.z0, s.dz, s.samples);
for (int gj = 0; gj < spec.ny; ++gj) {
const double worldY = spec.oy + gj * spec.dy;
for (int gi = 0; gi < spec.nx; ++gi) {
const double worldX = spec.ox + gi * spec.dx;
const int tIdx = nearestIndex(worldX, s.x0, s.dx, s.ntraces);
// X / Z 越界 → blank。
if (tIdx < 0 || sIdx < 0 || nChan == 0) {
out.vol.at(gi, gj, gk) = ScalarVolumeI16::kBlank;
continue;
}
// Y → 跨通道 1D 线性插值channelY 升序)。
double phys = 0.0;
bool blank = false;
if (worldY <= s.channelY.front()) {
// 边界外不外推:在 maxDist 内 clamp 到首通道,否则 blank。
if (s.channelY.front() - worldY > spec.maxDist) {
blank = true;
} else {
phys = s.at(0, tIdx, sIdx);
}
} else if (worldY >= s.channelY.back()) {
if (worldY - s.channelY.back() > spec.maxDist) {
blank = true;
} else {
phys = s.at(nChan - 1, tIdx, sIdx);
}
} else {
// 定位 worldY 落在哪两相邻通道间,线性插值。
int lo = 0;
while (lo + 1 < nChan && s.channelY[lo + 1] < worldY) ++lo;
const int hi = lo + 1;
const double yLo = s.channelY[lo];
const double yHi = s.channelY[hi];
const double span = yHi - yLo;
const double w = (span > 0.0) ? (worldY - yLo) / span : 0.0;
const double vLo = s.at(lo, tIdx, sIdx);
const double vHi = s.at(hi, tIdx, sIdx);
phys = vLo + (vHi - vLo) * w;
}
if (blank || std::isnan(phys)) {
out.vol.at(gi, gj, gk) = ScalarVolumeI16::kBlank;
} else {
out.vol.at(gi, gj, gk) = out.quant.toQ(phys);
}
}
}
}
return out;
}
} // namespace geopro::core

View File

@ -0,0 +1,27 @@
#ifndef GEOPRO_CORE_ALGO_GPRVOLUMEBUILDER_HPP
#define GEOPRO_CORE_ALGO_GPRVOLUMEBUILDER_HPP
#include <array>
#include "algo/IInterpolator.hpp" // GridSpec
#include "model/GprSurvey.hpp"
#include "model/ScalarVolumeI16.hpp" // ScalarVolumeI16, Quant
namespace geopro::core {
// 结构化建体产物int16 量化体 + 量化映射 + 几何origin/spacing+ 物理值域。
struct BuiltI16 {
ScalarVolumeI16 vol{0, 0, 0};
Quant quant;
std::array<double, 3> origin{{0, 0, 0}};
std::array<double, 3> spacing{{0, 0, 0}};
double vminPhys = 0, vmaxPhys = 0;
};
// 结构化建体X/Z 直接落格(取最近道/采样)+ 仅 Y 向跨通道 1D 线性插值。
// 超 X/Z 范围或 Y 越界且超 maxDist 的网格点 → ScalarVolumeI16::kBlank不外推
BuiltI16 buildGprVolume(const GprSurvey& s, const GridSpec& spec);
} // namespace geopro::core
#endif // GEOPRO_CORE_ALGO_GPRVOLUMEBUILDER_HPP

View File

@ -0,0 +1,34 @@
#ifndef GEOPRO_CORE_MODEL_GPRSURVEY_HPP
#define GEOPRO_CORE_MODEL_GPRSURVEY_HPP
#include <cstddef>
#include <vector>
namespace geopro::core {
// 规则化建体输入GPR 三维体)。雷达数据沿测线(X)、深度(Z)规则密采样,
// 仅横向/跨通道(Y)稀疏。结构化建体据此做 X/Z 落格 + 仅 Y 向 1D 线性插值。
//
// 放 core/modelgeopro::core以保持 geopro_core 自洽buildGprVolume 在 core
// 若 GprSurvey 落 io/gpr 则 core 需反向 include io故置于 core。
// 真实数据 → GprSurvey 的装配在后续 POC 台Task 9完成。
struct GprSurvey {
int ntraces = 0; // 沿测线道数X
int samples = 0; // 每道采样数Z 深度)
double x0 = 0, dx = 1; // 沿测线轴:第 t 道 X = x0 + t*dx
double z0 = 0, dz = 1; // 深度轴:第 s 采样 Z = z0 + s*dz
std::vector<double> channelY; // 各通道横向位置Y稀疏需按 Y 升序)
// 值channelY.size() × ntraces × samples布局 [(c*ntraces + t)*samples + s]。
std::vector<double> values;
// 64 位索引访问values[(size_t(c)*ntraces + t)*samples + s]。
double at(int c, int t, int s) const {
return values[(static_cast<std::size_t>(c) * ntraces + t) * samples + s];
}
};
} // namespace geopro::core
#endif // GEOPRO_CORE_MODEL_GPRSURVEY_HPP

View File

@ -33,6 +33,8 @@ target_sources(geopro_tests PRIVATE core/test_crs_transform.cpp)
target_sources(geopro_tests PRIVATE core/test_model_data.cpp) target_sources(geopro_tests PRIVATE core/test_model_data.cpp)
target_sources(geopro_tests PRIVATE core/test_geo_frame.cpp) target_sources(geopro_tests PRIVATE core/test_geo_frame.cpp)
target_sources(geopro_tests PRIVATE core/test_scalar_volume_i16.cpp) target_sources(geopro_tests PRIVATE core/test_scalar_volume_i16.cpp)
# GprVolumeBuilderX/Z + Y 1D 线 int16
target_sources(geopro_tests PRIVATE core/test_gpr_volume_builder.cpp)
target_link_libraries(geopro_tests PRIVATE geopro_core) target_link_libraries(geopro_tests PRIVATE geopro_core)
target_sources(geopro_tests PRIVATE data/test_parsers.cpp) target_sources(geopro_tests PRIVATE data/test_parsers.cpp)

View File

@ -0,0 +1,63 @@
#include "core/algo/GprVolumeBuilder.hpp"
#include "core/algo/IInterpolator.hpp" // GridSpec
#include <gtest/gtest.h>
using namespace geopro::core;
// GprSurvey 放 core/modelgeopro::core 命名空间),保持 core 自洽,避免 core->io 反向依赖。
TEST(GprVolumeBuilder, InterpolatesAcrossChannelsOnly) {
geopro::core::GprSurvey s;
s.ntraces = 1;
s.samples = 1;
s.x0 = 0;
s.dx = 1;
s.z0 = 0;
s.dz = 1;
s.channelY = {0.0, 1.0}; // 两通道Y=0 值0、Y=1 值100
s.values = {0.0, 100.0}; // [(c*1+0)*1+0] => c=0->0, c=1->100
GridSpec spec{};
spec.nx = 1;
spec.ny = 3;
spec.nz = 1;
spec.ox = 0;
spec.oy = 0;
spec.oz = 0;
spec.dx = 1;
spec.dy = 0.5;
spec.dz = 1;
spec.power = 2;
spec.maxDist = 9.9;
auto b = buildGprVolume(s, spec);
EXPECT_EQ(b.vol.nx(), 1);
EXPECT_EQ(b.vol.ny(), 3);
EXPECT_EQ(b.vol.nz(), 1);
EXPECT_NEAR(b.quant.toPhys(b.vol.at(0, 0, 0)), 0.0, 1.0); // Y=0 -> ch0
EXPECT_NEAR(b.quant.toPhys(b.vol.at(0, 1, 0)), 50.0, 1.5); // Y=0.5 -> 中点≈50
EXPECT_NEAR(b.quant.toPhys(b.vol.at(0, 2, 0)), 100.0, 1.0); // Y=1.0 -> ch1
}
TEST(GprVolumeBuilder, OutOfRangeBecomesBlank) {
geopro::core::GprSurvey s;
s.ntraces = 1;
s.samples = 1;
s.x0 = 0;
s.dx = 1;
s.z0 = 0;
s.dz = 1;
s.channelY = {0.0};
s.values = {5.0};
GridSpec spec{};
spec.nx = 3;
spec.ny = 1;
spec.nz = 1; // X 超出仅 1 道的范围
spec.ox = 0;
spec.oy = 0;
spec.oz = 0;
spec.dx = 10;
spec.dy = 1;
spec.dz = 1;
spec.power = 2;
spec.maxDist = 0.5;
auto b = buildGprVolume(s, spec);
EXPECT_EQ(b.vol.at(2, 0, 0), ScalarVolumeI16::kBlank); // 远端无覆盖->blank
}