Quick start for C++ API usageΒΆ

This is a short introduction to the PROJ C++ API. In the following section we create a simple program that illustrate how to transform points between two different coordinate systems.

Note: the essential osgeo::proj::operation::CoordinateOperation::coordinateTransformer() method was only added in PROJ 9.3. For earlier versions, coordinate operations must be exported as a PROJ string with the osgeo::proj::operation::CoordinateOperation::exportToPROJ() method, and passed to proj_create() to instance a PJ* object to use with proj_trans() (cf the Quick start for C API usage). The use of the C API proj_create_crs_to_crs() might still be easier to let PROJ automatically select the "best" coordinate operation when several ones are possible, as the corresponding automation is not currently available in the C++ API.

Before the PROJ API can be used it is necessary to include the various header files:

#include <cassert>
#include <cmath>   // for HUGE_VAL
#include <iomanip> // for std::setprecision()

#include <iostream>

#include "proj/coordinateoperation.hpp"
#include "proj/crs.hpp"
#include "proj/io.hpp"
#include "proj/util.hpp" // for nn_dynamic_pointer_cast

For convenience, we also declare using a few namespaces:

using namespace NS_PROJ::crs;
using namespace NS_PROJ::io;
using namespace NS_PROJ::operation;
using namespace NS_PROJ::util;

We start by creating a database context (osgeo::proj::io::DatabaseContext) with the default settings to find the PROJ database.

auto dbContext = DatabaseContext::create();

We then instantiate a generic authority factory (osgeo::proj::io::AuthorityFactory), that is not tied to a particular authority, to be able to get transformations registered by different authorities. This can only be used for a osgeo::proj::operation::CoordinateOperationContext, and not to instantiate objects of the database which are all tied to a non-generic authority.

    auto authFactory = AuthorityFactory::create(dbContext, std::string());

We create a coordinate operation context, that can be customized to ammend the way coordinate operations are computed. Here we ask for default settings, as we have a coordinate operation that just involves a "simple" map projection in the same datum.

auto coord_op_ctxt =
    CoordinateOperationContext::create(authFactory, nullptr, 0.0);

We instantiate a authority factory for EPSG related objects.

auto authFactoryEPSG = AuthorityFactory::create(dbContext, "EPSG");

We instantiate the source CRS from its code: 4326, for WGS 84 latitude/longitude.

auto sourceCRS = authFactoryEPSG->createCoordinateReferenceSystem("4326");

We instantiate the source CRS from its PROJ.4 string (it would be possible to instantiate it from its 32631 code, similarly to above), and cast the generic osgeo::proj::util::BaseObject to the osgeo::proj::crs::CRS class required later.

auto targetCRS =
    NN_CHECK_THROW(nn_dynamic_pointer_cast<CRS>(createFromUserInput(
        "+proj=utm +zone=31 +datum=WGS84 +type=crs", dbContext)));

Warning

The use of PROJ strings to describe a CRS is not recommended. One of the main weaknesses of PROJ strings is their inability to describe a geodetic datum, other than the few ones hardcoded in the +datum parameter.

We ask for the list of operations available to transform from the source to the target CRS with the osgeo::proj::operation::CoordinateOperationFactory::createOperations() method.

auto list = CoordinateOperationFactory::create()->createOperations(
    sourceCRS, targetCRS, coord_op_ctxt);

We check that we got a non-empty list of operations. The list is sorted from the most relevant to the less relevant one. Cf Filtering and sorting of coordinate operations for more details on the sorting of those operations. For a transformation between a projected CRS and its base CRS, like we do here, there will be only one operation.

assert(!list.empty());

We create an execution context (must only be used by one thread at a time) with the proj_context_create() function.

PJ_CONTEXT *ctx = proj_context_create();

We create a coordinate transformer (osgeo::proj::operation::CoordinateTransformer) from the first operation of the list:

auto transformer = list[0]->coordinateTransformer(ctx);

We can now transform a point with the osgeo::proj::operation::CoordinateTransformer::transform() method. Note that the the expected input values should be passed in the order and the unit of the successive axis of the input CRS. Similarly the values returned in the v[] array of the output PJ_COORD are in the order and the unit of the successive axis of the output CRS. For coordinate operations involving a time-dependent operation, coord.v[3] is the decimal year of the coordinate epoch of the input (or HUGE_VAL to indicate none).

PJ_COORD c = {{
    49.0,    // latitude in degree
    2.0,     // longitude in degree
    0.0,     // z ordinate. unused
    HUGE_VAL // time ordinate. unused
}};
c = transformer->transform(c);

and output the result:

std::cout << std::fixed << std::setprecision(3);
std::cout << "Easting: " << c.v[0] << std::endl;  // should be 426857.988
std::cout << "Northing: " << c.v[1] << std::endl; // should be 5427937.523

We need to clean up the PJ_CONTEXT handle before exiting with the proj_context_destroy() function.

proj_context_destroy(ctx);

A complete compilable version of the example code can be seen below:

 1#include <cassert>
 2#include <cmath>   // for HUGE_VAL
 3#include <iomanip> // for std::setprecision()
 4
 5#include <iostream>
 6
 7#include "proj/coordinateoperation.hpp"
 8#include "proj/crs.hpp"
 9#include "proj/io.hpp"
10#include "proj/util.hpp" // for nn_dynamic_pointer_cast
11
12using namespace NS_PROJ::crs;
13using namespace NS_PROJ::io;
14using namespace NS_PROJ::operation;
15using namespace NS_PROJ::util;
16
17int main(void) {
18    auto dbContext = DatabaseContext::create();
19
20    // Instantiate a generic authority factory, that is not tied to a particular
21    // authority, to be able to get transformations registered by different
22    // authorities. This can only be used for CoordinateOperationContext.
23    auto authFactory = AuthorityFactory::create(dbContext, std::string());
24
25    // Create a coordinate operation context, that can be customized to ammend
26    // the way coordinate operations are computed. Here we ask for default
27    // settings.
28    auto coord_op_ctxt =
29        CoordinateOperationContext::create(authFactory, nullptr, 0.0);
30
31    // Instantiate a authority factory for EPSG related objects.
32    auto authFactoryEPSG = AuthorityFactory::create(dbContext, "EPSG");
33
34    // Instantiate source CRS from EPSG code
35    auto sourceCRS = authFactoryEPSG->createCoordinateReferenceSystem("4326");
36
37    // Instantiate target CRS from PROJ.4 string (commented out, the equivalent
38    // from the EPSG code)
39    // auto targetCRS =
40    // authFactoryEPSG->createCoordinateReferenceSystem("32631");
41    auto targetCRS =
42        NN_CHECK_THROW(nn_dynamic_pointer_cast<CRS>(createFromUserInput(
43            "+proj=utm +zone=31 +datum=WGS84 +type=crs", dbContext)));
44
45    // List operations available to transform from EPSG:4326
46    // (WGS 84 latitude/longitude) to EPSG:32631 (WGS 84 / UTM zone 31N).
47    auto list = CoordinateOperationFactory::create()->createOperations(
48        sourceCRS, targetCRS, coord_op_ctxt);
49
50    // Check that we got a non-empty list of operations
51    // The list is sorted from the most relevant to the less relevant one.
52    // Cf
53    // https://proj.org/operations/operations_computation.html#filtering-and-sorting-of-coordinate-operations
54    // for more details on the sorting of those operations.
55    // For a transformation between a projected CRS and its base CRS, like
56    // we do here, there will be only one operation.
57    assert(!list.empty());
58
59    // Create an execution context (must only be used by one thread at a time)
60    PJ_CONTEXT *ctx = proj_context_create();
61
62    // Create a coordinate transformer from the first operation of the list
63    auto transformer = list[0]->coordinateTransformer(ctx);
64
65    // Perform the coordinate transformation.
66    PJ_COORD c = {{
67        49.0,    // latitude in degree
68        2.0,     // longitude in degree
69        0.0,     // z ordinate. unused
70        HUGE_VAL // time ordinate. unused
71    }};
72    c = transformer->transform(c);
73
74    // Display result
75    std::cout << std::fixed << std::setprecision(3);
76    std::cout << "Easting: " << c.v[0] << std::endl;  // should be 426857.988
77    std::cout << "Northing: " << c.v[1] << std::endl; // should be 5427937.523
78
79    // Destroy execution context
80    proj_context_destroy(ctx);
81
82    return 0;
83}