Skip to content

Commit

Permalink
added preliminary support for +axis (#18)
Browse files Browse the repository at this point in the history
git-svn-id: http://svn.osgeo.org/metacrs/proj/trunk@1825 4e78687f-474d-0410-85f9-8d5e500ac6b2
  • Loading branch information
warmerdam committed Mar 1, 2010
1 parent d3753ed commit 3405bf7
Show file tree
Hide file tree
Showing 5 changed files with 174 additions and 4 deletions.
3 changes: 3 additions & 0 deletions ChangeLog
@@ -1,5 +1,8 @@
2010-02-28 Frank Warmerdam <warmerdam@pobox.com>

* src/pj_init.c, src/pj_transform.c: added support for +axis setting
to control axis orientation (#18).

* nad/epsg: Regenerated from EPSG 7.4.1 with the big datum selection
upgrade.

Expand Down
28 changes: 28 additions & 0 deletions src/pj_init.c
Expand Up @@ -270,6 +270,7 @@ pj_init(int argc, char **argv) {
PIN->is_geocent = 0;
PIN->is_long_wrap_set = 0;
PIN->long_wrap_center = 0.0;
strcpy( PIN->axis, "enu" );

/* set datum parameters */
if (pj_datum_set(start, PIN)) goto bum_call;
Expand Down Expand Up @@ -304,6 +305,33 @@ pj_init(int argc, char **argv) {
PIN->over = pj_param(start, "bover").i;

/* longitude center for wrapping */
PIN->is_long_wrap_set = pj_param(start, "tlon_wrap").i;
if (PIN->is_long_wrap_set)
PIN->long_wrap_center = pj_param(start, "rlon_wrap").f;

/* axis orientation */
if( (pj_param(start,"saxis").s) != NULL )
{
static const char *axis_legal = "ewnsud";
const char *axis_arg = pj_param(start,"saxis").s;
if( strlen(axis_arg) != 3 )
{
pj_errno = PJD_ERR_AXIS;
goto bum_call;
}

if( strchr( axis_legal, axis_arg[0] ) == NULL
|| strchr( axis_legal, axis_arg[1] ) == NULL
|| (axis_arg[2] && strchr( axis_legal, axis_arg[1] ) == NULL))
{
pj_errno = PJD_ERR_AXIS;
goto bum_call;
}

/* it would be nice to validate we don't have on axis repeated */
strcpy( PIN->axis, axis_arg );
}

PIN->is_long_wrap_set = pj_param(start, "tlon_wrap").i;
if (PIN->is_long_wrap_set)
PIN->long_wrap_center = pj_param(start, "rlon_wrap").f;
Expand Down
1 change: 1 addition & 0 deletions src/pj_strerrno.c
Expand Up @@ -52,6 +52,7 @@ pj_err_list[] = {
"unparseable coordinate system definition", /* -44 */
"geocentric transformation missing z or ellps", /* -45 */
"unknown prime meridian conversion id", /* -46 */
"illegal axis orientation combination", /* -47 */
};
char *
pj_strerrno(int err)
Expand Down
140 changes: 138 additions & 2 deletions src/pj_transform.c
Expand Up @@ -36,6 +36,10 @@

PJ_CVSID("$Id$");

static int pj_adjust_axis( const char *axis, int denormalize_flag,
long point_count, int point_offset,
double *x, double *y, double *z );

#ifndef SRS_WGS84_SEMIMAJOR
#define SRS_WGS84_SEMIMAJOR 6378137.0
#endif
Expand Down Expand Up @@ -63,13 +67,13 @@ PJ_CVSID("$Id$");
** list or something, but while experimenting with it this should be fine.
*/

static const int transient_error[45] = {
static const int transient_error[50] = {
/* 0 1 2 3 4 5 6 7 8 9 */
/* 0 to 9 */ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
/* 10 to 19 */ 0, 0, 0, 0, 1, 1, 0, 1, 1, 1,
/* 20 to 29 */ 1, 0, 0, 0, 0, 0, 0, 1, 0, 0,
/* 30 to 39 */ 0, 0, 0, 0, 0, 0, 0, 0, 1, 0,
/* 40 to 44 */ 0, 0, 0, 0, 0 };
/* 40 to 49 */ 0, 0, 0, 0, 0, 0, 0, 0, 0, 0 };

/************************************************************************/
/* pj_transform() */
Expand All @@ -92,6 +96,20 @@ int pj_transform( PJ *srcdefn, PJ *dstdefn, long point_count, int point_offset,
if( point_offset == 0 )
point_offset = 1;

/* -------------------------------------------------------------------- */
/* Transform unusual input coordinate axis orientation to */
/* standard form if needed. */
/* -------------------------------------------------------------------- */
if( strcmp(srcdefn->axis,"enu") != 0 )
{
int err;

err = pj_adjust_axis( srcdefn->axis, 0, point_count, point_offset,
x, y, z );
if( err != 0 )
return err;
}

/* -------------------------------------------------------------------- */
/* Transform geocentric source coordinates to lat/long. */
/* -------------------------------------------------------------------- */
Expand Down Expand Up @@ -282,6 +300,20 @@ int pj_transform( PJ *srcdefn, PJ *dstdefn, long point_count, int point_offset,
}
}

/* -------------------------------------------------------------------- */
/* Transform normalized axes into unusual output coordinate axis */
/* orientation if needed. */
/* -------------------------------------------------------------------- */
if( strcmp(dstdefn->axis,"enu") != 0 )
{
int err;

err = pj_adjust_axis( dstdefn->axis, 1, point_count, point_offset,
x, y, z );
if( err != 0 )
return err;
}

return 0;
}

Expand Down Expand Up @@ -644,3 +676,107 @@ int pj_datum_transform( PJ *srcdefn, PJ *dstdefn,
return 0;
}

/************************************************************************/
/* pj_adjust_axis() */
/* */
/* Normalize or de-normalized the x/y/z axes. The normal form */
/* is "enu" (easting, northing, up). */
/************************************************************************/
static int pj_adjust_axis( const char *axis, int denormalize_flag,
long point_count, int point_offset,
double *x, double *y, double *z )

{
double x_in, y_in, z_in = 0.0;
int i, i_axis;

if( !denormalize_flag )
{
for( i = 0; i < point_count; i++ )
{
x_in = x[point_offset*i];
y_in = y[point_offset*i];
if( z )
z_in = z[point_offset*i];

for( i_axis = 0; i_axis < 3; i_axis++ )
{
double value;

if( i_axis == 0 )
value = x_in;
else if( i_axis == 1 )
value = y_in;
else
value = z_in;

switch( axis[i_axis] )
{
case 'e':
x[point_offset*i] = value; break;
case 'w':
x[point_offset*i] = -value; break;
case 'n':
y[point_offset*i] = value; break;
case 's':
y[point_offset*i] = -value; break;
case 'u':
if( z ) z[point_offset*i] = value; break;
case 'd':
if( z ) z[point_offset*i] = -value; break;
default:
pj_errno = PJD_ERR_AXIS;
return PJD_ERR_AXIS;
}
} /* i_axis */
} /* i (point) */
}

else /* denormalize */
{
for( i = 0; i < point_count; i++ )
{
x_in = x[point_offset*i];
y_in = y[point_offset*i];
if( z )
z_in = z[point_offset*i];

for( i_axis = 0; i_axis < 3; i_axis++ )
{
double *target;

if( i_axis == 2 && z == NULL )
continue;

if( i_axis == 0 )
target = x;
else if( i_axis == 1 )
target = y;
else
target = z;

switch( axis[i_axis] )
{
case 'e':
target[point_offset*i] = x_in; break;
case 'w':
target[point_offset*i] = -x_in; break;
case 'n':
target[point_offset*i] = y_in; break;
case 's':
target[point_offset*i] = -y_in; break;
case 'u':
target[point_offset*i] = z_in; break;
case 'd':
target[point_offset*i] = -z_in; break;
default:
pj_errno = PJD_ERR_AXIS;
return PJD_ERR_AXIS;
}
} /* i_axis */
} /* i (point) */
}

return 0;
}

6 changes: 4 additions & 2 deletions src/projects.h
Expand Up @@ -131,8 +131,9 @@ extern double hypot(double, double);
#define PJD_GRIDSHIFT 3
#define PJD_WGS84 4 /* WGS84 (or anything considered equivelent) */

/* datum system errors */
#define PJD_ERR_GEOCENTRIC -45
/* library errors */
#define PJD_ERR_GEOCENTRIC -45
#define PJD_ERR_AXIS -47

#define USE_PROJUV

Expand Down Expand Up @@ -235,6 +236,7 @@ typedef struct PJconsts {
double from_greenwich; /* prime meridian offset (in radians) */
double long_wrap_center; /* 0.0 for -180 to 180, actually in radians*/
int is_long_wrap_set;
char axis[4];

#ifdef PROJ_PARMS__
PROJ_PARMS__
Expand Down

0 comments on commit 3405bf7

Please sign in to comment.