[GRASS-SVN] r72448 - grass/trunk/lib/proj
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue Mar 20 15:24:58 PDT 2018
Author: mmetz
Date: 2018-03-20 15:24:58 -0700 (Tue, 20 Mar 2018)
New Revision: 72448
Modified:
grass/trunk/lib/proj/do_proj.c
grass/trunk/lib/proj/do_proj_ll.c
Log:
libproj: simplify PROJ 5 usage
Modified: grass/trunk/lib/proj/do_proj.c
===================================================================
--- grass/trunk/lib/proj/do_proj.c 2018-03-20 22:16:17 UTC (rev 72447)
+++ grass/trunk/lib/proj/do_proj.c 2018-03-20 22:24:58 UTC (rev 72448)
@@ -53,6 +53,28 @@
return pjdef;
}
+
+/* TODO: add to gprojects.h */
+/* Create a transformation object */
+int GPJ_prepare_pjinfo(const struct pj_info *info_in,
+ const struct pj_info *info_out,
+ struct pj_info *info_new)
+{
+ char *projdef, *projdefin, *projdefout;
+
+ projdefin = gpj_get_def(info_in->pj);
+ projdefout = gpj_get_def(info_out->pj);
+ projdef = NULL;
+ G_asprintf(&projdef, "+proj=pipeline +step %s +inv +step %s", projdefin, projdefout);
+ info_new->pj = proj_create(PJ_DEFAULT_CTX, projdef);
+ if (info_new->pj == NULL)
+ G_fatal_error(_("proj_create() failed"));
+ G_free(projdefin);
+ G_free(projdefout);
+ G_free(projdef);
+
+ return 1;
+}
#endif
/* TODO: rename pj_ to GPJ_ to avoid symbol clash with PROJ lib */
@@ -79,20 +101,10 @@
{
int ok;
#ifdef HAVE_PROJ_H
- PJ *P;
- char *projdef, *projdefin, *projdefout;
+ struct pj_info info_new;
PJ_COORD c;
- projdefin = gpj_get_def(info_in->pj);
- projdefout = gpj_get_def(info_out->pj);
- projdef = NULL;
- G_asprintf(&projdef, "+proj=pipeline +step %s +inv +step %s", projdefin, projdefout);
- P = proj_create(PJ_DEFAULT_CTX, projdef);
- if (P == NULL)
- G_fatal_error(_("proj_create() failed"));
- G_free(projdefin);
- G_free(projdefout);
- G_free(projdef);
+ GPJ_prepare_pjinfo(info_in, info_out, &info_new);
METERS_in = info_in->meters;
METERS_out = info_out->meters;
@@ -103,8 +115,8 @@
c.lpzt.phi = (*y) / RAD_TO_DEG;
c.lpzt.z = 0;
c.lpzt.t = 0;
- c = proj_trans(P, PJ_FWD, c);
- ok = proj_errno(P);
+ c = proj_trans(info_new.pj, PJ_FWD, c);
+ ok = proj_errno(info_new.pj);
if (strncmp(info_out->proj, "ll", 2) == 0) {
/* convert to degrees */
@@ -123,8 +135,8 @@
c.xyzt.y = *y * METERS_in;
c.xyzt.z = 0;
c.xyzt.t = 0;
- c = proj_trans(P, PJ_FWD, c);
- ok = proj_errno(P);
+ c = proj_trans(info_new.pj, PJ_FWD, c);
+ ok = proj_errno(info_new.pj);
if (strncmp(info_out->proj, "ll", 2) == 0) {
/* convert to degrees */
@@ -137,7 +149,7 @@
*y = c.xy.y / METERS_out;
}
}
- proj_destroy(P);
+ proj_destroy(info_new.pj);
if (ok < 0) {
G_warning(_("proj_trans() failed: %d"), ok);
@@ -218,20 +230,10 @@
int has_h = 1;
#ifdef HAVE_PROJ_H
int i;
- PJ *P;
- char *projdef, *projdefin, *projdefout;
+ struct pj_info info_new;
PJ_COORD c;
- projdefin = gpj_get_def(info_in->pj);
- projdefout = gpj_get_def(info_out->pj);
- projdef = NULL;
- G_asprintf(&projdef, "+proj=pipeline +step %s +inv +step %s", projdefin, projdefout);
- P = proj_create(PJ_DEFAULT_CTX, projdef);
- if (P == NULL)
- G_fatal_error(_("proj_create() failed"));
- G_free(projdefin);
- G_free(projdefout);
- G_free(projdef);
+ GPJ_prepare_pjinfo(info_in, info_out, &info_new);
METERS_in = info_in->meters;
METERS_out = info_out->meters;
@@ -252,8 +254,8 @@
c.lpzt.lam = x[i] / RAD_TO_DEG;
c.lpzt.phi = y[i] / RAD_TO_DEG;
c.lpzt.z = h[i];
- c = proj_trans(P, PJ_FWD, c);
- if ((ok = proj_errno(P)) < 0)
+ c = proj_trans(info_new.pj, PJ_FWD, c);
+ if ((ok = proj_errno(info_new.pj)) < 0)
break;
/* convert to degrees */
x[i] = c.lp.lam * RAD_TO_DEG;
@@ -266,8 +268,8 @@
c.lpzt.lam = x[i] / RAD_TO_DEG;
c.lpzt.phi = y[i] / RAD_TO_DEG;
c.lpzt.z = h[i];
- c = proj_trans(P, PJ_FWD, c);
- if ((ok = proj_errno(P)) < 0)
+ c = proj_trans(info_new.pj, PJ_FWD, c);
+ if ((ok = proj_errno(info_new.pj)) < 0)
break;
/* convert to map units */
x[i] = c.xy.x / METERS_out;
@@ -283,8 +285,8 @@
c.xyzt.x = x[i] * METERS_in;
c.xyzt.y = y[i] * METERS_in;
c.xyzt.z = h[i];
- c = proj_trans(P, PJ_FWD, c);
- if ((ok = proj_errno(P)) < 0)
+ c = proj_trans(info_new.pj, PJ_FWD, c);
+ if ((ok = proj_errno(info_new.pj)) < 0)
break;
/* convert to degrees */
x[i] = c.lp.lam * RAD_TO_DEG;
@@ -297,8 +299,8 @@
c.xyzt.x = x[i] * METERS_in;
c.xyzt.y = y[i] * METERS_in;
c.xyzt.z = h[i];
- c = proj_trans(P, PJ_FWD, c);
- if ((ok = proj_errno(P)) < 0)
+ c = proj_trans(info_new.pj, PJ_FWD, c);
+ if ((ok = proj_errno(info_new.pj)) < 0)
break;
/* convert to map units */
x[i] = c.xy.x / METERS_out;
@@ -308,7 +310,7 @@
}
if (!has_h)
G_free(h);
- proj_destroy(P);
+ proj_destroy(info_new.pj);
if (ok < 0) {
G_warning(_("proj_trans() failed: %d"), ok);
Modified: grass/trunk/lib/proj/do_proj_ll.c
===================================================================
--- grass/trunk/lib/proj/do_proj_ll.c 2018-03-20 22:16:17 UTC (rev 72447)
+++ grass/trunk/lib/proj/do_proj_ll.c 2018-03-20 22:24:58 UTC (rev 72448)
@@ -24,17 +24,6 @@
static double METERS_in = 1.0;
-#ifdef HAVE_PROJ_H
-static char *gpj_get_def(PJ *P)
-{
- char *pjdef;
- PJ_PROJ_INFO pj_proj_info = proj_pj_info(P);
-
- pjdef = G_store(pj_proj_info.definition);
-
- return pjdef;
-}
-
/**
* \brief Re-project a point between a co-ordinate system and its
* latitude / longitude equivalent, using the same datum,
@@ -43,8 +32,8 @@
* This function takes a pointer to a pj_info structure as argument,
* and projects a point between the co-ordinate system and its
* latitude / longitude equivalent. With forward projection, the point
- * is projected from the co-ordinate system to its latitude / longitude
- * equivalent, while backward projection does the reverse.
+ * is projected from the equivalent latitude / longitude system to the
+ * projected co-ordinate system, while backward projection does the reverse.
* The easting and northing of the point are contained in two pointers
* passed to the function; these will be overwritten by the co-ordinates
* of the re-projected point.
@@ -58,11 +47,11 @@
**/
int GPJ_do_proj_ll(double *x, double *y,
- const struct pj_info *info_in, int direction)
+ const struct pj_info *info_in, int direction)
{
int ok;
- PJ *P;
- char *projdef;
+
+#ifdef HAVE_PROJ_H
PJ_COORD c;
ok = 0;
@@ -71,12 +60,6 @@
return ok;
}
- projdef = gpj_get_def(info_in->pj);
- P = proj_create(PJ_DEFAULT_CTX, projdef);
- if (P == NULL)
- G_fatal_error(_("proj_create() failed"));
- G_free(projdef);
-
METERS_in = info_in->meters;
if (direction == PJ_FWD) {
@@ -88,8 +71,8 @@
c.lpzt.z = 0;
c.lpzt.t = 0;
- c = proj_trans(P, PJ_FWD, c);
- ok = proj_errno(P);
+ c = proj_trans(info_in->pj, PJ_FWD, c);
+ ok = proj_errno(info_in->pj);
/* convert to map units */
*x = c.xy.x / METERS_in;
@@ -104,18 +87,34 @@
c.xyzt.z = 0;
c.xyzt.t = 0;
- c = proj_trans(P, PJ_INV, c);
- ok = proj_errno(P);
+ c = proj_trans(info_in->pj, PJ_INV, c);
+ ok = proj_errno(info_in->pj);
/* convert to degrees */
*x = c.lp.lam * RAD_TO_DEG;
*y = c.lp.phi * RAD_TO_DEG;
}
- proj_destroy(P);
if (ok < 0) {
G_warning(_("proj_trans() failed: %d"), ok);
}
+#else
+ struct pj_info info_out;
+
+ if ((ok = GPJ_get_equivalent_latlong(&info_out, info_in)) != 1)
+ return ok;
+
+ if (direction == PJ_FWD) {
+ /* from ll to projected */
+ ok = pj_do_proj(x, y, &info_out, info_in);
+ }
+ else {
+ /* from projected to ll */
+ ok = pj_do_proj(x, y, info_in, &info_out);
+ }
+ pj_free(info_out.pj);
+#endif
+
return ok;
}
@@ -127,8 +126,8 @@
* This function takes a pointer to a pj_info structure as argument,
* and projects an array of points between the co-ordinate system and its
* latitude / longitude equivalent. With forward projection, the points
- * are projected from the co-ordinate system to its latitude / longitude
- * equivalent, while backward projection does the reverse.
+ * are projected from the equivalent latitude / longitude system to the
+ * projected co-ordinate system, while backward projection does the reverse.
* The easting and northing of the points are contained in two pointers passed
* to the function; these will be overwritten by the co-ordinates of the
* re-projected point.
@@ -149,10 +148,10 @@
const struct pj_info *info_in, int direction)
{
int ok;
+
+#ifdef HAVE_PROJ_H
int has_h = 1;
int i;
- PJ *P;
- char *projdef;
PJ_COORD c;
ok = 0;
@@ -161,12 +160,6 @@
return ok;
}
- projdef = gpj_get_def(info_in->pj);
- P = proj_create(PJ_DEFAULT_CTX, projdef);
- if (P == NULL)
- G_fatal_error(_("proj_create() failed"));
- G_free(projdef);
-
METERS_in = info_in->meters;
if (h == NULL) {
@@ -187,8 +180,8 @@
c.lpzt.z = h[i];
c.lpzt.t = 0;
- c = proj_trans(P, PJ_FWD, c);
- if ((ok = proj_errno(P)) < 0)
+ c = proj_trans(info_in->pj, PJ_FWD, c);
+ if ((ok = proj_errno(info_in->pj)) < 0)
break;
/* convert to map units */
@@ -206,8 +199,8 @@
c.xyzt.z = h[i];
c.xyzt.t = 0;
- c = proj_trans(P, PJ_INV, c);
- if ((ok = proj_errno(P)) < 0)
+ c = proj_trans(info_in->pj, PJ_INV, c);
+ if ((ok = proj_errno(info_in->pj)) < 0)
break;
/* convert to degrees */
@@ -217,13 +210,27 @@
}
if (!has_h)
G_free(h);
- proj_destroy(P);
if (ok < 0) {
G_warning(_("proj_trans() failed: %d"), ok);
}
+#else
+ struct pj_info info_out;
+
+ if ((ok = GPJ_get_equivalent_latlong(&info_out, info_in)) != 1)
+ return ok;
+
+ if (direction == PJ_FWD) {
+ /* from ll to projected */
+ ok = pj_do_transform(count, x, y, h, &info_out, info_in);
+ }
+ else {
+ /* from projected to ll */
+ ok = pj_do_transform(count, x, y, h, info_in, &info_out);
+ }
+ pj_free(info_out.pj);
+#endif
return ok;
}
-#endif /* HAVE_PROJ_H */
More information about the grass-commit
mailing list