c++boostgeospatialgeoboost-geometry

using boost::geometry to convert from Latitude and longitude to UTM

boost::geometry seems to have support to convert from lat/lon to UTM. Unfortunately I couldn't find any example on how exactly to do that. Does someone have an example they are willing to share?

Solution

  • I just spent an embarrassing amount of time looking for this (multiple hours).

    It seems like everything spherical (SRS) is just ... not documented at all. Perhaps there's something broken with the documentation generator.

    Regardless, several tedious searches through test code further I stumbled on a working answer.

    Here's examples

    With even more stumbling I found the corresponding EPSG code on https://epsg.io/?q=UTM+31N: enter image description here

    Live On Compiler Explorer

    #include <iostream>
#include <boost/geometry.hpp>
#include <boost/geometry/core/coordinate_system.hpp>
#include <boost/geometry/geometries/point_xy.hpp>

#include <boost/geometry/srs/epsg.hpp>
#ifdef COLIRU
#include <boost/geometry/srs/projection.hpp>
#endif

namespace bg = boost::geometry;
namespace bm = bg::model::d2;
namespace srs = bg::srs;

using LongLat = bm::point_xy<double, bg::cs::geographic<bg::degree> >;
using UTM     = bm::point_xy<double/*, srs::static_epsg<3043>*/>;

constexpr LongLat Amsterdam() { return { 4.897, 52.371}; }
constexpr LongLat Barcelona() { return { 2.1734, 41.3851 }; }

void report(LongLat src, auto const& prj) {
    UTM r {};
    prj.forward(src, r);
    std::cout << std::fixed << bg::wkt(src) << " -> " << bg::wkt(r) << "\n";
}

int main() {
#ifndef COLIRU
    // dynamic projection factory too heavy on Coliru
    srs::projection<> zone31 = srs::proj4("+proj=utm +ellps=GRS80 +lat_1=0.5 +lat_2=2 +n=0.5 +zone=31");
    report(Amsterdam(), zone31);
    report(Barcelona(), zone31);
#endif

    srs::projection<srs::static_epsg<3043> > epsg3043;
    report(Amsterdam(), epsg3043);
    report(Barcelona(), epsg3043);
}

    Prints

    POINT(4.897000 52.371000) -> POINT(629144.771310 5803996.656944)
POINT(2.173400 41.385100) -> POINT(430887.564331 4581837.853239)
POINT(4.897000 52.371000) -> POINT(629144.771310 5803996.656944)
POINT(2.173400 41.385100) -> POINT(430887.564331 4581837.853239)

    Disclaimer:

    I'm not actually convinced that the UTM point type I defined makes much sense. I think any point type that can receive 2 coordinates will do. At least the srs::... as the coordinate system doesn't seem to do anything - bg::convert, bg::assign and friends don't magically know how to apply the projection anyways. I think the "doxygen_z_article09" link you included is just badly out of date/not how things were eventually implemented.