sij/mapnik
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

wip: support custom non-crs proj transforms again

Browse source 
closes: #4306
This commit is contained in:
Johannes Schauer Marin Rodrigues 2022-04-15 07:21:48 +02:00
parent 1ba1278b42
commit 10801d4b48
No known key found for this signature in database
GPG key ID: F2CBA5C78FBD83E1
2 changed files with 100 additions and 23 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

118
src/proj_transform.cpp
View file

@ -35,6 +35,7 @@
#ifdef MAPNIK_USE_PROJ
// proj
#include <proj.h>
#include <proj/io.hpp>
#endif
// stl
@ -92,6 +93,29 @@ auto envelope_points(box2d<T> const& env, std::size_t num_points) -> geometry::m
} // namespace
// from src/4D_api.cpp
/** Adds a " +type=crs" suffix to a PROJ string (if it is a PROJ string) */
std::string pj_add_type_crs_if_needed(const std::string &str) {
std::string ret(str);
if ((str.rfind("proj=", 0) == 0 || str.rfind("+proj=", 0) == 0 ||
str.rfind("+init=", 0) == 0 || str.rfind("+title=", 0) == 0) &&
str.find("type=crs") == std::string::npos) {
ret += " +type=crs";
}
return ret;
}
void throw_proj_exception(PJ_CONTEXT *ctx, std::string msg)
{
throw std::runtime_error(msg
#if MAPNIK_PROJ_VERSION >= 80000
+ " because of "
+ std::string(proj_context_errno_string(ctx, proj_context_errno(ctx))));
#else
);
#endif
}
proj_transform::proj_transform(projection const& source, projection const& dest)
: is_source_longlat_(false)
, is_dest_longlat_(false)
@ -125,30 +149,80 @@ proj_transform::proj_transform(projection const& source, projection const& dest)
#ifdef MAPNIK_USE_PROJ
ctx_ = proj_context_create();
proj_log_level(ctx_, PJ_LOG_ERROR);
transform_ = proj_create_crs_to_crs(ctx_, source.params().c_str(), dest.params().c_str(), nullptr);
if (transform_ == nullptr)
{
throw std::runtime_error(std::string("Cannot initialize proj_transform (crs_to_crs) for given projections: '") +
source.params() + "'->'" + dest.params() +
#if MAPNIK_PROJ_VERSION >= 80000
"' because of " + std::string(proj_context_errno_string(ctx_, proj_context_errno(ctx_))));
#else
"'");
#endif
// we replicate what proj_create_crs_to_crs() does and then call
// proj_create_crs_to_crs_from_pj because we need a PJ object for
// proj_is_crs and we want to avoid running proj_create twice as
// as much (once ourselves and once in proj_create_crs_to_crs)
std::string crs_source = pj_add_type_crs_if_needed(source.params());
std::string crs_dest = pj_add_type_crs_if_needed(dest.params());
PJ* src;
PJ* dst;
try {
src = proj_create(ctx_, crs_source.c_str());
} catch( const std::exception& ) {
src = nullptr;
}
PJ* transform_gis = proj_normalize_for_visualization(ctx_, transform_);
if (transform_gis == nullptr)
{
throw std::runtime_error(std::string("Cannot initialize proj_transform (normalize) for given projections: '") +
source.params() + "'->'" + dest.params() +
#if MAPNIK_PROJ_VERSION >= 80000
"' because of " + std::string(proj_context_errno_string(ctx_, proj_context_errno(ctx_))));
#else
"'");
#endif
if (!src) {
throw_proj_exception(ctx_,
std::string("Cannot instantiate source projection for '") + source.params() + "'");
}
proj_destroy(transform_);
transform_ = transform_gis;
try {
dst = proj_create(ctx_, crs_dest.c_str());
} catch( const std::exception& ) {
dst = nullptr;
}
if (!dst) {
throw_proj_exception(ctx_,
std::string("Cannot instantiate dest projection for '") + dest.params() + "'");
}
if (proj_is_crs(src) && proj_is_crs(dst)) {
// if both source and destination are a CRS, we can use
// proj_create_crs_to_crs
transform_ = proj_create_crs_to_crs_from_pj(ctx_, src, dst, nullptr, nullptr);
if (transform_ == nullptr)
{
throw_proj_exception(ctx_,
std::string("Cannot initialize proj_transform (crs_to_crs) for given projections: '") +
source.params() + "'->'" + dest.params() + "'");
}
PJ* transform_gis = proj_normalize_for_visualization(ctx_, transform_);
if (transform_gis == nullptr)
{
throw_proj_exception(ctx_,
std::string("Cannot initialize proj_transform (normalize) for given projections: '") +
source.params() + "'->'" + dest.params() + "'");
}
proj_destroy(transform_);
transform_ = transform_gis;
} else {
// if one of the projections is not a CRS, then this was
// probably a user-supplied transform and we have to create
// the transform ourselves
auto formatter = osgeo::proj::io::PROJStringFormatter::create();
formatter->startInversion();
try {
formatter->ingestPROJString(source.params());
} catch (const osgeo::proj::io::ParsingException &e) {
throw_proj_exception(ctx_, std::string("Failed to ingest source: '") + source.params() + "'");
}
formatter->stopInversion();
try {
formatter->ingestPROJString(dest.params());
} catch (const osgeo::proj::io::ParsingException &e) {
throw_proj_exception(ctx_, std::string("Failed to ingest dest: '") + dest.params() + "'");
}
try {
transform_ = proj_create(ctx_, formatter->toString().c_str());
} catch( const std::exception& ) {
transform_ = nullptr;
}
if (!transform_) {
throw_proj_exception(ctx_,
std::string("Cannot instantiate custom projection pipeline: '") + formatter->toString() + "'");
}
}
proj_destroy(src);
proj_destroy(dst);
#else
throw std::runtime_error(
std::string(

5
test/unit/projection/proj_transform.cpp
View file

@ -88,13 +88,16 @@ TEST_CASE("projection transform")
mapnik::projection proj_4937("epsg:4937");
//"Webmercator" WGS 84 / Pseudo-Mercator -- Spherical Mercator, Google Maps, OpenStreetMap, Bing, ArcGIS, ESRI
mapnik::projection proj_3857("epsg:3857");
//axisswap, which is not a CRS
mapnik::projection proj_axisswap("+proj=axisswap +order=2,1");
mapnik::proj_transform tr1(proj_4236, proj_27700);
mapnik::proj_transform tr2(proj_4236, proj_8857);
mapnik::proj_transform tr3(proj_4236, proj_4236);
mapnik::proj_transform tr4(proj_4236, proj_4937);
mapnik::proj_transform tr5(proj_4236, proj_3857);
std::initializer_list<std::reference_wrapper<mapnik::proj_transform>> transforms = {tr1, tr2, tr3, tr4, tr5};
mapnik::proj_transform tr6(proj_4236, proj_axisswap);
std::initializer_list<std::reference_wrapper<mapnik::proj_transform>> transforms = {tr1, tr2, tr3, tr4, tr5, tr6};
std::initializer_list<std::tuple<double, double>> coords = {
{-4.0278869, 57.8796955}, // Dórnach, Highland