C library for Geodesics: Information for PROJ user (original) (raw)
The library is included as part of PROJ starting with version 4.9.0 (released in 2015), where it is used as the computational backend for geod(1) and to implement some projections.
If PROJ is installed on your system, that you can compiled your application with cmake. A skeleton CMakeLists.txt
is
cmake_minimum_required (VERSION 3.13.0)
project (proj-example C)
find_package (PROJ REQUIRED)
add_executable (proj-example proj-example.c)
target_link_libraries (proj-example ${PROJ_LIBRARIES})
This compiles and links the sample code [proj-example.c](proj-example%5F8c.html "An example of computing geodesics with the PROJ library.")
.
#include <stdio.h>
#include <proj.h>
int main (void) {
double a, invf;
const int n = 30;
double lat[n+1], lon[n+1], azi[n+1];
int i;
{
PJ_CONTEXT* C = proj_context_create();
char* argv[3] = {"type=crs", "proj=longlat", "ellps=WGS84"};
PJ* P = proj_create_argv(C, 3, argv);
if (P == 0) {
fprintf(stderr, "Failed to create transformation object.\n");
return 1;
}
proj_ellipsoid_get_parameters(C, proj_get_ellipsoid(C, P),
&a, 0, 0, &invf);
proj_destroy(P);
proj_context_destroy(C);
}
geod_init(&g, a, invf != 0 ? 1/invf : 0);
printf("latitude longitude azimuth\n");
for (i = 0; i <= n; ++i) {
i * l.a13 / n, lat + i, lon + i, azi + i, 0, 0, 0, 0, 0);
printf("%.5f %.5f %.5f\n", lat[i], lon[i], azi[i]);
}
return 0;
}
API for the geodesic routines in C.
void GEOD_DLL geod_init(struct geod_geodesic *g, double a, double f)
double GEOD_DLL geod_genposition(const struct geod_geodesicline *l, unsigned flags, double s12_a12, double *plat2, double *plon2, double *pazi2, double *ps12, double *pm12, double *pM12, double *pM21, double *pS12)
void GEOD_DLL geod_inverseline(struct geod_geodesicline *l, const struct geod_geodesic *g, double lat1, double lon1, double lat2, double lon2, unsigned caps)
In this example, PROJ provides the following
- access to the header file
[geodesic.h](geodesic%5F8h.html "API for the geodesic routines in C.")
; - access to the geodesic code bundled into the PROJ library;
- a mechanism for retrieving the parameters (equatorial radius and flattening) defining the ellipsoid.
To compile and run this example, do
cmake -B BUILD -S . cd BUILD make ./proj-example
which will print out 31 way points from Los Angeles to Christchurch starting with
latitude longitude azimuth 33.94000 -118.41000 -136.41544 31.49945 -121.08843 -137.86373 29.00508 -123.62942 -139.14437 ...
Other points to note for PROJ programmers:
- The shape of the ellipsoid is given by the flattening f. f = 0 denotes a sphere. f < 0 denotes a prolate ellipsoid.
- The size of the ellipsoid is given by the equatorial radius a. This is usually given in meters. However, any units may be used and these are the units used for all distance-like quantites computed by the geodesic routines.
- All angle like quantities are expressed in degrees, not radians.