diff --git a/src/core/CMakeLists.txt b/src/core/CMakeLists.txt index a0d2a72..f1fb57d 100644 --- a/src/core/CMakeLists.txt +++ b/src/core/CMakeLists.txt @@ -1,13 +1,15 @@ find_package(Eigen3 CONFIG REQUIRED) +find_package(PROJ CONFIG REQUIRED) add_library(geopro_core STATIC geo/LocalFrame.cpp + geo/CrsTransform.cpp model/ColorScale.cpp algo/IdwInterpolator.cpp ) target_include_directories(geopro_core PUBLIC ${CMAKE_CURRENT_SOURCE_DIR}) -target_link_libraries(geopro_core PUBLIC Eigen3::Eigen) +target_link_libraries(geopro_core PUBLIC Eigen3::Eigen PROJ::proj) target_compile_features(geopro_core PUBLIC cxx_std_17) # 铁律:core 为纯业务逻辑层,绝不链接 Qt / VTK。 diff --git a/src/core/geo/CrsTransform.cpp b/src/core/geo/CrsTransform.cpp new file mode 100644 index 0000000..c49f734 --- /dev/null +++ b/src/core/geo/CrsTransform.cpp @@ -0,0 +1,37 @@ +#include "geo/CrsTransform.hpp" +#include +#include +namespace geopro::core { + +struct CrsTransform::Impl { + PJ_CONTEXT* ctx = nullptr; + PJ* pj = nullptr; +}; + +CrsTransform::CrsTransform(const std::string& src, const std::string& dst) + : impl_(std::make_unique()) { + impl_->ctx = proj_context_create(); + impl_->pj = proj_create_crs_to_crs(impl_->ctx, src.c_str(), dst.c_str(), nullptr); + if (!impl_->pj) throw std::runtime_error("CrsTransform: failed to create " + src + "->" + dst); + PJ* norm = proj_normalize_for_visualization(impl_->ctx, impl_->pj); + if (norm) { proj_destroy(impl_->pj); impl_->pj = norm; } +} + +CrsTransform::~CrsTransform() { + if (impl_->pj) proj_destroy(impl_->pj); + if (impl_->ctx) proj_context_destroy(impl_->ctx); +} + +Xy CrsTransform::forward(double x, double y) const { + PJ_COORD c = proj_coord(x, y, 0, 0); + PJ_COORD r = proj_trans(impl_->pj, PJ_FWD, c); + return Xy{r.xy.x, r.xy.y}; +} + +Xy CrsTransform::inverse(double x, double y) const { + PJ_COORD c = proj_coord(x, y, 0, 0); + PJ_COORD r = proj_trans(impl_->pj, PJ_INV, c); + return Xy{r.xy.x, r.xy.y}; +} + +} // namespace geopro::core diff --git a/src/core/geo/CrsTransform.hpp b/src/core/geo/CrsTransform.hpp new file mode 100644 index 0000000..6ac7409 --- /dev/null +++ b/src/core/geo/CrsTransform.hpp @@ -0,0 +1,20 @@ +#pragma once +#include +#include +namespace geopro::core { + +struct Xy { double x, y; }; + +// PROJ 封装:源 CRS <-> 目标 CRS。用于世界系<->GIS<->经纬度、影像异源 CRS 重投影(设计 §5)。 +class CrsTransform { +public: + CrsTransform(const std::string& srcCrs, const std::string& dstCrs); + ~CrsTransform(); + Xy forward(double x, double y) const; // src -> dst + Xy inverse(double x, double y) const; // dst -> src +private: + struct Impl; + std::unique_ptr impl_; +}; + +} // namespace geopro::core diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt index 4d7024a..d296baf 100644 --- a/tests/CMakeLists.txt +++ b/tests/CMakeLists.txt @@ -7,12 +7,25 @@ add_executable(geopro_tests smoke_test.cpp) target_link_libraries(geopro_tests PRIVATE GTest::gtest GTest::gtest_main) include(GoogleTest) -gtest_discover_tests(geopro_tests) + +# PROJ 运行时需要数据目录(proj.db 等)。vcpkg 把数据装在 share/proj 下。 +# 通过测试环境变量 PROJ_DATA 注入,使 `ctest -R CrsTransform` 无需调用方手动 set。 +file(GLOB _proj_data_dirs + "${CMAKE_BINARY_DIR}/vcpkg_installed/*/share/proj" +) +if(_proj_data_dirs) + list(GET _proj_data_dirs 0 GEOPRO_PROJ_DATA) + gtest_discover_tests(geopro_tests + PROPERTIES ENVIRONMENT "PROJ_DATA=${GEOPRO_PROJ_DATA}") +else() + gtest_discover_tests(geopro_tests) +endif() target_sources(geopro_tests PRIVATE core/test_local_frame.cpp) target_sources(geopro_tests PRIVATE core/test_model.cpp) target_sources(geopro_tests PRIVATE core/test_color_scale.cpp) target_sources(geopro_tests PRIVATE core/test_idw.cpp) +target_sources(geopro_tests PRIVATE core/test_crs_transform.cpp) target_link_libraries(geopro_tests PRIVATE geopro_core) add_subdirectory(spike) # spike S3: banded contour 渲染验证 diff --git a/tests/core/test_crs_transform.cpp b/tests/core/test_crs_transform.cpp new file mode 100644 index 0000000..1c83454 --- /dev/null +++ b/tests/core/test_crs_transform.cpp @@ -0,0 +1,26 @@ +#include +#include "geo/CrsTransform.hpp" +using namespace geopro::core; + +// UTM 49N (CM 111°E, false easting 500000) <-> WGS84 经纬度往返。 +// 输入 easting≈516868 仅比中央经线偏东 ~16.9km,故经度落在 ~111.16°E(非 114°)。 +// normalize_for_visualization 后:forward 输出 (lon, lat),inverse 输入 (lon, lat)。 +TEST(CrsTransform, Utm49nToWgs84RoundTrip) { + CrsTransform t("EPSG:32649", "EPSG:4326"); + auto ll = t.forward(516868.0, 2494259.0); // (east,north) -> (lon,lat) + EXPECT_NEAR(ll.x, 111.16, 0.05); // lon + EXPECT_NEAR(ll.y, 22.55, 0.05); // lat + auto en = t.inverse(ll.x, ll.y); // (lon,lat) -> (east,north) + EXPECT_NEAR(en.x, 516868.0, 1.0); + EXPECT_NEAR(en.y, 2494259.0, 1.0); +} + +// WebMercator -> UTM 49N 异源重投影(设计 §5 影像重投影)。 +// tfw 原点经度真实约 114.16°E,距 49N 中央经线 111°E 偏东 ~3.16°, +// 故其 49N easting 远大于 false easting,落在 ~825km 处(合理且 < UTM 带宽上限)。 +TEST(CrsTransform, WebMercatorToUtm) { + CrsTransform t("EPSG:3857", "EPSG:32649"); + auto p = t.forward(12708343.88, 2577685.90); // tfw 原点 + EXPECT_GT(p.x, 800000.0); EXPECT_LT(p.x, 850000.0); // UTM 49N 合理 Easting + EXPECT_GT(p.y, 2400000.0); EXPECT_LT(p.y, 2600000.0); // 合理 Northing +} diff --git a/vcpkg.json b/vcpkg.json index 22cff44..77981b0 100644 --- a/vcpkg.json +++ b/vcpkg.json @@ -5,6 +5,7 @@ "dependencies": [ "eigen3", "gtest", - "nlohmann-json" + "nlohmann-json", + "proj" ] }