aboutsummaryrefslogtreecommitdiff
path: root/src/apps/cs2cs.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'src/apps/cs2cs.cpp')
-rw-r--r--src/apps/cs2cs.cpp640
1 files changed, 640 insertions, 0 deletions
diff --git a/src/apps/cs2cs.cpp b/src/apps/cs2cs.cpp
new file mode 100644
index 00000000..f63bedcc
--- /dev/null
+++ b/src/apps/cs2cs.cpp
@@ -0,0 +1,640 @@
+/******************************************************************************
+ * Project: PROJ.4
+ * Purpose: Mainline program sort of like ``proj'' for converting between
+ * two coordinate systems.
+ * Author: Frank Warmerdam, warmerda@home.com
+ *
+ ******************************************************************************
+ * Copyright (c) 2000, Frank Warmerdam
+ *
+ * Permission is hereby granted, free of charge, to any person obtaining a
+ * copy of this software and associated documentation files (the "Software"),
+ * to deal in the Software without restriction, including without limitation
+ * the rights to use, copy, modify, merge, publish, distribute, sublicense,
+ * and/or sell copies of the Software, and to permit persons to whom the
+ * Software is furnished to do so, subject to the following conditions:
+ *
+ * The above copyright notice and this permission notice shall be included
+ * in all copies or substantial portions of the Software.
+ *
+ * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
+ * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+ * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
+ * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+ * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
+ * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
+ * DEALINGS IN THE SOFTWARE.
+ *****************************************************************************/
+
+#define FROM_PROJ_CPP
+
+#include <ctype.h>
+#include <locale.h>
+#include <math.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+
+#include <cassert>
+#include <string>
+
+#include <proj/internal/internal.hpp>
+
+// PROJ include order is sensitive
+// clang-format off
+#include "proj.h"
+#include "projects.h"
+#include "emess.h"
+// clang-format on
+
+#define MAX_LINE 1000
+
+static PJ *transformation = nullptr;
+
+static bool srcIsGeog = false;
+static double srcToRadians = 0.0;
+
+static bool destIsGeog = false;
+static double destToRadians = 0.0;
+static bool destIsLatLong = false;
+
+static int reversein = 0, /* != 0 reverse input arguments */
+ reverseout = 0, /* != 0 reverse output arguments */
+ echoin = 0, /* echo input data to output line */
+ tag = '#'; /* beginning of line tag character */
+
+static const char *oform =
+ nullptr; /* output format for x-y or decimal degrees */
+static char oform_buffer[16]; /* buffer for oform when using -d */
+static const char *oterr = "*\t*"; /* output line for unprojectable input */
+static const char *usage =
+ "%s\nusage: %s [ -dDeEfIlrstvwW [args] ] [ +opts[=arg] ]\n"
+ " [+to [+opts[=arg] [ files ]\n";
+
+static double (*informat)(const char *,
+ char **); /* input data deformatter function */
+
+/************************************************************************/
+/* process() */
+/* */
+/* File processing function. */
+/************************************************************************/
+static void process(FILE *fid)
+
+{
+ char line[MAX_LINE + 3], *s, pline[40];
+ projUV data;
+
+ for (;;) {
+ double z;
+
+ ++emess_dat.File_line;
+ if (!(s = fgets(line, MAX_LINE, fid)))
+ break;
+ if (!strchr(s, '\n')) { /* overlong line */
+ int c;
+ (void)strcat(s, "\n");
+ /* gobble up to newline */
+ while ((c = fgetc(fid)) != EOF && c != '\n')
+ ;
+ }
+ if (*s == tag) {
+ fputs(line, stdout);
+ continue;
+ }
+
+ if (reversein) {
+ data.v = (*informat)(s, &s);
+ data.u = (*informat)(s, &s);
+ } else {
+ data.u = (*informat)(s, &s);
+ data.v = (*informat)(s, &s);
+ }
+
+ z = strtod(s, &s);
+
+ if (data.v == HUGE_VAL)
+ data.u = HUGE_VAL;
+
+ if (!*s && (s > line))
+ --s; /* assumed we gobbled \n */
+
+ if (echoin) {
+ char t;
+ t = *s;
+ *s = '\0';
+ (void)fputs(line, stdout);
+ *s = t;
+ putchar('\t');
+ }
+
+ if (data.u != HUGE_VAL) {
+
+ if (srcIsGeog) {
+ /* dmstor gives values to radians. Convert now to the SRS unit
+ */
+ data.u /= srcToRadians;
+ data.v /= srcToRadians;
+ }
+
+ PJ_COORD coord;
+ coord.xyzt.x = data.u;
+ coord.xyzt.y = data.v;
+ coord.xyzt.z = z;
+ coord.xyzt.t = HUGE_VAL;
+ coord = proj_trans(transformation, PJ_FWD, coord);
+ data.u = coord.xyz.x;
+ data.v = coord.xyz.y;
+ z = coord.xyz.z;
+ }
+
+ if (data.u == HUGE_VAL) /* error output */
+ fputs(oterr, stdout);
+
+ else if (destIsGeog && !oform) { /*ascii DMS output */
+
+ // rtodms() expect radians: convert from the output SRS unit
+ data.u *= destToRadians;
+ data.v *= destToRadians;
+
+ if (destIsLatLong) {
+ if (reverseout) {
+ fputs(rtodms(pline, data.v, 'E', 'W'), stdout);
+ putchar('\t');
+ fputs(rtodms(pline, data.u, 'N', 'S'), stdout);
+ } else {
+ fputs(rtodms(pline, data.u, 'N', 'S'), stdout);
+ putchar('\t');
+ fputs(rtodms(pline, data.v, 'E', 'W'), stdout);
+ }
+ } else if (reverseout) {
+ fputs(rtodms(pline, data.v, 'N', 'S'), stdout);
+ putchar('\t');
+ fputs(rtodms(pline, data.u, 'E', 'W'), stdout);
+ } else {
+ fputs(rtodms(pline, data.u, 'E', 'W'), stdout);
+ putchar('\t');
+ fputs(rtodms(pline, data.v, 'N', 'S'), stdout);
+ }
+
+ } else { /* x-y or decimal degree ascii output */
+ if (destIsGeog) {
+ data.v *= destToRadians * RAD_TO_DEG;
+ data.u *= destToRadians * RAD_TO_DEG;
+ }
+ if (reverseout) {
+ printf(oform, data.v);
+ putchar('\t');
+ printf(oform, data.u);
+ } else {
+ printf(oform, data.u);
+ putchar('\t');
+ printf(oform, data.v);
+ }
+ }
+
+ putchar(' ');
+ if (oform != nullptr)
+ printf(oform, z);
+ else
+ printf("%.3f", z);
+ if (s)
+ printf("%s", s);
+ else
+ printf("\n");
+ }
+}
+
+/************************************************************************/
+/* instanciate_crs() */
+/************************************************************************/
+
+static PJ_OBJ *instanciate_crs(const std::string &definition,
+ const char *const *optionsImportCRS,
+ bool &isGeog, double &toRadians,
+ bool &isLatFirst) {
+ PJ_OBJ *crs = proj_obj_create_from_user_input(nullptr, definition.c_str(),
+ optionsImportCRS);
+ if (!crs) {
+ return nullptr;
+ }
+
+ isGeog = false;
+ toRadians = 0.0;
+ isLatFirst = false;
+
+ auto type = proj_obj_get_type(crs);
+ if (type == PJ_OBJ_TYPE_BOUND_CRS) {
+ auto base = proj_obj_get_source_crs(nullptr, crs);
+ proj_obj_destroy(crs);
+ crs = base;
+ type = proj_obj_get_type(crs);
+ }
+ if (type == PJ_OBJ_TYPE_GEOGRAPHIC_2D_CRS ||
+ type == PJ_OBJ_TYPE_GEOGRAPHIC_3D_CRS) {
+ auto cs = proj_obj_crs_get_coordinate_system(nullptr, crs);
+ assert(cs);
+
+ isGeog = true;
+ const char *axisName = "";
+ proj_obj_cs_get_axis_info(nullptr, cs, 0,
+ &axisName, // name,
+ nullptr, // abbrev
+ nullptr, // direction
+ &toRadians,
+ nullptr, // unit name
+ nullptr, // unit authority
+ nullptr // unit code
+ );
+ isLatFirst =
+ NS_PROJ::internal::ci_find(std::string(axisName), "latitude") !=
+ std::string::npos;
+
+ proj_obj_destroy(cs);
+ }
+
+ return crs;
+}
+
+/************************************************************************/
+/* get_geog_crs_proj_string_from_proj_crs() */
+/************************************************************************/
+
+static std::string get_geog_crs_proj_string_from_proj_crs(PJ_OBJ *src,
+ double &toRadians,
+ bool &isLatFirst) {
+ auto srcType = proj_obj_get_type(src);
+ if (srcType == PJ_OBJ_TYPE_BOUND_CRS) {
+ auto base = proj_obj_get_source_crs(nullptr, src);
+ assert(base);
+ proj_obj_destroy(src);
+ src = base;
+ srcType = proj_obj_get_type(src);
+ }
+ if (srcType != PJ_OBJ_TYPE_PROJECTED_CRS) {
+ return std::string();
+ }
+
+ auto base = proj_obj_get_source_crs(nullptr, src);
+ assert(base);
+ auto baseType = proj_obj_get_type(base);
+ if (baseType != PJ_OBJ_TYPE_GEOGRAPHIC_2D_CRS &&
+ baseType != PJ_OBJ_TYPE_GEOGRAPHIC_3D_CRS) {
+ proj_obj_destroy(base);
+ return std::string();
+ }
+
+ auto cs = proj_obj_crs_get_coordinate_system(nullptr, base);
+ assert(cs);
+
+ const char *axisName = "";
+ proj_obj_cs_get_axis_info(nullptr, cs, 0,
+ &axisName, // name,
+ nullptr, // abbrev
+ nullptr, // direction
+ &toRadians,
+ nullptr, // unit name
+ nullptr, // unit authority
+ nullptr // unit code
+ );
+ isLatFirst = NS_PROJ::internal::ci_find(std::string(axisName),
+ "latitude") != std::string::npos;
+
+ proj_obj_destroy(cs);
+
+ auto retCStr = proj_obj_as_proj_string(nullptr, base, PJ_PROJ_5, nullptr);
+ std::string ret(retCStr ? retCStr : "");
+ proj_obj_destroy(base);
+ return ret;
+}
+
+/************************************************************************/
+/* main() */
+/************************************************************************/
+
+int main(int argc, char **argv) {
+ char *arg;
+ char **eargv = argv;
+ std::string fromStr;
+ std::string toStr;
+ FILE *fid;
+ int eargc = 0, mon = 0;
+ int have_to_flag = 0, inverse = 0;
+ int use_env_locale = 0;
+
+ /* This is just to check that pj_init() is locale-safe */
+ /* Used by nad/testvarious */
+ if (getenv("PROJ_USE_ENV_LOCALE") != nullptr)
+ use_env_locale = 1;
+
+ /* Enable compatibility mode for init=epsg:XXXX by default */
+ if (getenv("PROJ_USE_PROJ4_INIT_RULES") == nullptr) {
+ proj_context_use_proj4_init_rules(nullptr, true);
+ }
+
+ if ((emess_dat.Prog_name = strrchr(*argv, DIR_CHAR)) != nullptr)
+ ++emess_dat.Prog_name;
+ else
+ emess_dat.Prog_name = *argv;
+ inverse = !strncmp(emess_dat.Prog_name, "inv", 3);
+ if (argc <= 1) {
+ (void)fprintf(stderr, usage, pj_get_release(), emess_dat.Prog_name);
+ exit(0);
+ }
+
+ // First pass to check if we have "cs2cs [-bla]* <SRC> <DEST>" syntax
+ int countNonOptionArg = 0;
+ for (int i = 1; i < argc; i++) {
+ if (argv[i][0] == '-') {
+ if (argv[i][1] == '\0') {
+ countNonOptionArg++;
+ }
+ } else {
+ if (strcmp(argv[i], "+to") == 0) {
+ countNonOptionArg = -1;
+ break;
+ }
+ countNonOptionArg++;
+ }
+ }
+ const bool isSrcDestSyntax = (countNonOptionArg == 2);
+
+ /* process run line arguments */
+ while (--argc > 0) { /* collect run line arguments */
+ if (**++argv == '-') {
+ for (arg = *argv;;) {
+ switch (*++arg) {
+ case '\0': /* position of "stdin" */
+ if (arg[-1] == '-')
+ eargv[eargc++] = const_cast<char *>("-");
+ break;
+ case 'v': /* monitor dump of initialization */
+ mon = 1;
+ continue;
+ case 'I': /* alt. method to spec inverse */
+ inverse = 1;
+ continue;
+ case 'E': /* echo ascii input to ascii output */
+ echoin = 1;
+ continue;
+ case 't': /* set col. one char */
+ if (arg[1])
+ tag = *++arg;
+ else
+ emess(1, "missing -t col. 1 tag");
+ continue;
+ case 'l': /* list projections, ellipses or units */
+ if (!arg[1] || arg[1] == 'p' || arg[1] == 'P') {
+ /* list projections */
+ const struct PJ_LIST *lp;
+ int do_long = arg[1] == 'P', c;
+ const char *str;
+
+ for (lp = proj_list_operations(); lp->id; ++lp) {
+ (void)printf("%s : ", lp->id);
+ if (do_long) /* possibly multiline description */
+ (void)puts(*lp->descr);
+ else { /* first line, only */
+ str = *lp->descr;
+ while ((c = *str++) && c != '\n')
+ putchar(c);
+ putchar('\n');
+ }
+ }
+ } else if (arg[1] == '=') { /* list projection 'descr' */
+ const struct PJ_LIST *lp;
+
+ arg += 2;
+ for (lp = proj_list_operations(); lp->id; ++lp)
+ if (!strcmp(lp->id, arg)) {
+ (void)printf("%9s : %s\n", lp->id, *lp->descr);
+ break;
+ }
+ } else if (arg[1] == 'e') { /* list ellipses */
+ const struct PJ_ELLPS *le;
+
+ for (le = proj_list_ellps(); le->id; ++le)
+ (void)printf("%9s %-16s %-16s %s\n", le->id,
+ le->major, le->ell, le->name);
+ } else if (arg[1] == 'u') { /* list units */
+ const struct PJ_UNITS *lu;
+
+ for (lu = proj_list_units(); lu->id; ++lu)
+ (void)printf("%12s %-20s %s\n", lu->id,
+ lu->to_meter, lu->name);
+ } else if (arg[1] == 'd') { /* list datums */
+ const struct PJ_DATUMS *ld;
+
+ printf("__datum_id__ __ellipse___ "
+ "__definition/"
+ "comments______________________________\n");
+ for (ld = pj_get_datums_ref(); ld->id; ++ld) {
+ printf("%12s %-12s %-30s\n", ld->id, ld->ellipse_id,
+ ld->defn);
+ if (ld->comments != nullptr &&
+ strlen(ld->comments) > 0)
+ printf("%25s %s\n", " ", ld->comments);
+ }
+ } else if (arg[1] == 'm') { /* list prime meridians */
+ const struct PJ_PRIME_MERIDIANS *lpm;
+
+ for (lpm = proj_list_prime_meridians(); lpm->id; ++lpm)
+ (void)printf("%12s %-30s\n", lpm->id, lpm->defn);
+ } else
+ emess(1, "invalid list option: l%c", arg[1]);
+ exit(0);
+ /* cppcheck-suppress duplicateBreak */
+ continue; /* artificial */
+ case 'e': /* error line alternative */
+ if (--argc <= 0)
+ noargument:
+ emess(1, "missing argument for -%c", *arg);
+ oterr = *++argv;
+ continue;
+ case 'W': /* specify seconds precision */
+ case 'w': /* -W for constant field width */
+ {
+ char c = arg[1];
+ if (c != 0 && isdigit(c)) {
+ set_rtodms(c - '0', *arg == 'W');
+ ++arg;
+ } else
+ emess(1, "-W argument missing or non-digit");
+ continue;
+ }
+ case 'f': /* alternate output format degrees or xy */
+ if (--argc <= 0)
+ goto noargument;
+ oform = *++argv;
+ continue;
+ case 'r': /* reverse input */
+ reversein = 1;
+ continue;
+ case 's': /* reverse output */
+ reverseout = 1;
+ continue;
+ case 'D': /* set debug level */
+ if (--argc <= 0)
+ goto noargument;
+ pj_ctx_set_debug(pj_get_default_ctx(), atoi(*++argv));
+ continue;
+ case 'd':
+ if (--argc <= 0)
+ goto noargument;
+ sprintf(oform_buffer, "%%.%df", atoi(*++argv));
+ oform = oform_buffer;
+ break;
+ default:
+ emess(1, "invalid option: -%c", *arg);
+ break;
+ }
+ break;
+ }
+ } else if (isSrcDestSyntax) {
+ if (fromStr.empty())
+ fromStr = *argv;
+ else
+ toStr = *argv;
+ } else if (strcmp(*argv, "+to") == 0) {
+ have_to_flag = 1;
+
+ } else if (**argv == '+') { /* + argument */
+ if (have_to_flag) {
+ if (!toStr.empty())
+ toStr += ' ';
+ toStr += *argv;
+ } else {
+ if (!fromStr.empty())
+ fromStr += ' ';
+ fromStr += *argv;
+ }
+ } else if (!have_to_flag) {
+ fromStr = *argv;
+ } else if (toStr.empty()) {
+ toStr = *argv;
+ } else /* assumed to be input file name(s) */
+ eargv[eargc++] = *argv;
+ }
+ if (eargc == 0) /* if no specific files force sysin */
+ eargv[eargc++] = const_cast<char *>("-");
+
+ /*
+ * If the user has requested inverse, then just reverse the
+ * coordinate systems.
+ */
+ if (inverse) {
+ std::swap(fromStr, toStr);
+ }
+
+ if (use_env_locale) {
+ /* Set locale from environment */
+ setlocale(LC_ALL, "");
+ }
+
+ if (fromStr.empty() && toStr.empty()) {
+ emess(3, "missing source and target coordinate systems");
+ }
+
+ const char *const optionsProj4Mode[] = {"USE_PROJ4_INIT_RULES=YES",
+ nullptr};
+ const char *const *optionsImportCRS =
+ proj_context_get_use_proj4_init_rules(nullptr, TRUE) ? optionsProj4Mode
+ : nullptr;
+
+ PJ_OBJ *src = nullptr;
+ if (!fromStr.empty()) {
+ bool ignored;
+ src = instanciate_crs(fromStr, optionsImportCRS, srcIsGeog,
+ srcToRadians, ignored);
+ if (!src) {
+ emess(3, "cannot instanciate source coordinate system");
+ }
+ }
+
+ PJ_OBJ *dst = nullptr;
+ if (!toStr.empty()) {
+ dst = instanciate_crs(toStr, optionsImportCRS, destIsGeog,
+ destToRadians, destIsLatLong);
+ if (!dst) {
+ emess(3, "cannot instanciate target coordinate system");
+ }
+ }
+
+ if (toStr.empty()) {
+ assert(src);
+ toStr = get_geog_crs_proj_string_from_proj_crs(src, destToRadians,
+ destIsLatLong);
+ if (toStr.empty()) {
+ emess(3,
+ "missing target CRS and source CRS is not a projected CRS");
+ }
+ destIsGeog = true;
+ } else if (fromStr.empty()) {
+ assert(dst);
+ bool ignored;
+ fromStr =
+ get_geog_crs_proj_string_from_proj_crs(dst, srcToRadians, ignored);
+ if (fromStr.empty()) {
+ emess(3,
+ "missing source CRS and target CRS is not a projected CRS");
+ }
+ srcIsGeog = true;
+ }
+
+ proj_obj_destroy(src);
+ proj_obj_destroy(dst);
+
+ transformation = proj_create_crs_to_crs(nullptr, fromStr.c_str(),
+ toStr.c_str(), nullptr);
+ if (!transformation) {
+ emess(3, "cannot initialize transformation\ncause: %s",
+ pj_strerrno(pj_errno));
+ }
+
+ if (use_env_locale) {
+ /* Restore C locale to avoid issues in parsing/outputting numbers*/
+ setlocale(LC_ALL, "C");
+ }
+
+ if (mon) {
+ printf("%c ---- From Coordinate System ----\n", tag);
+ printf("%s\n", fromStr.c_str());
+ printf("%c ---- To Coordinate System ----\n", tag);
+ printf("%s\n", toStr.c_str());
+ }
+
+ /* set input formatting control */
+ if (!srcIsGeog)
+ informat = strtod;
+ else {
+ informat = dmstor;
+ }
+
+ if (!destIsGeog && !oform)
+ oform = "%.2f";
+
+ /* process input file list */
+ for (; eargc--; ++eargv) {
+ if (**eargv == '-') {
+ fid = stdin;
+ emess_dat.File_name = const_cast<char *>("<stdin>");
+
+ } else {
+ if ((fid = fopen(*eargv, "rt")) == nullptr) {
+ emess(-2, *eargv, "input file");
+ continue;
+ }
+ emess_dat.File_name = *eargv;
+ }
+ emess_dat.File_line = 0;
+ process(fid);
+ fclose(fid);
+ emess_dat.File_name = nullptr;
+ }
+
+ proj_destroy(transformation);
+
+ pj_deallocate_grids();
+
+ exit(0); /* normal completion */
+}