PostGIS  2.1.10dev-r@@SVN_REVISION@@
liblwgeom/lwgeom_transform.c
Go to the documentation of this file.
1 /**********************************************************************
2  *
3  * PostGIS - Spatial Types for PostgreSQL
4  * http://postgis.net
5  *
6  * Copyright (C) 2001-2003 Refractions Research Inc.
7  *
8  * This is free software; you can redistribute and/or modify it under
9  * the terms of the GNU General Public Licence. See the COPYING file.
10  *
11  **********************************************************************/
12 
13 #include "../postgis_config.h"
14 #include "liblwgeom.h"
15 #include "lwgeom_log.h"
16 #include <string.h>
17 
18 
20 static void
22 {
23  pt->x *= M_PI/180.0;
24  pt->y *= M_PI/180.0;
25 }
26 
28 static void
30 {
31  pt->x *= 180.0/M_PI;
32  pt->y *= 180.0/M_PI;
33 }
34 
39 int
40 ptarray_transform(POINTARRAY *pa, projPJ inpj, projPJ outpj)
41 {
42  int i;
43  POINT4D p;
44 
45  for ( i = 0; i < pa->npoints; i++ )
46  {
47  getPoint4d_p(pa, i, &p);
48  if ( ! point4d_transform(&p, inpj, outpj) ) return LW_FAILURE;
49  ptarray_set_point4d(pa, i, &p);
50  }
51 
52  return LW_SUCCESS;
53 }
54 
55 
60 int
61 lwgeom_transform(LWGEOM *geom, projPJ inpj, projPJ outpj)
62 {
63  int i;
64 
65  /* No points to transform in an empty! */
66  if ( lwgeom_is_empty(geom) )
67  return LW_SUCCESS;
68 
69  switch(geom->type)
70  {
71  case POINTTYPE:
72  case LINETYPE:
73  case CIRCSTRINGTYPE:
74  case TRIANGLETYPE:
75  {
76  LWLINE *g = (LWLINE*)geom;
77  if ( ! ptarray_transform(g->points, inpj, outpj) ) return LW_FAILURE;
78  break;
79  }
80  case POLYGONTYPE:
81  {
82  LWPOLY *g = (LWPOLY*)geom;
83  for ( i = 0; i < g->nrings; i++ )
84  {
85  if ( ! ptarray_transform(g->rings[i], inpj, outpj) ) return LW_FAILURE;
86  }
87  break;
88  }
89  case MULTIPOINTTYPE:
90  case MULTILINETYPE:
91  case MULTIPOLYGONTYPE:
92  case COLLECTIONTYPE:
93  case COMPOUNDTYPE:
94  case CURVEPOLYTYPE:
95  case MULTICURVETYPE:
96  case MULTISURFACETYPE:
98  case TINTYPE:
99  {
100  LWCOLLECTION *g = (LWCOLLECTION*)geom;
101  for ( i = 0; i < g->ngeoms; i++ )
102  {
103  if ( ! lwgeom_transform(g->geoms[i], inpj, outpj) ) return LW_FAILURE;
104  }
105  break;
106  }
107  default:
108  {
109  lwerror("lwgeom_transform: Cannot handle type '%s'",
110  lwtype_name(geom->type));
111  return LW_FAILURE;
112  }
113  }
114  return LW_SUCCESS;
115 }
116 
117 int
118 point4d_transform(POINT4D *pt, projPJ srcpj, projPJ dstpj)
119 {
120  int* pj_errno_ref;
121  POINT4D orig_pt;
122 
123  /* Make a copy of the input point so we can report the original should an error occur */
124  orig_pt.x = pt->x;
125  orig_pt.y = pt->y;
126  orig_pt.z = pt->z;
127 
128  if (pj_is_latlong(srcpj)) to_rad(pt) ;
129 
130  LWDEBUGF(4, "transforming POINT(%f %f) from '%s' to '%s'", orig_pt.x, orig_pt.y, pj_get_def(srcpj,0), pj_get_def(dstpj,0));
131 
132  /* Perform the transform */
133  pj_transform(srcpj, dstpj, 1, 0, &(pt->x), &(pt->y), &(pt->z));
134 
135  /* For NAD grid-shift errors, display an error message with an additional hint */
136  pj_errno_ref = pj_get_errno_ref();
137 
138  if (*pj_errno_ref != 0)
139  {
140  if (*pj_errno_ref == -38)
141  {
142  lwnotice("PostGIS was unable to transform the point because either no grid shift files were found, or the point does not lie within the range for which the grid shift is defined. Refer to the ST_Transform() section of the PostGIS manual for details on how to configure PostGIS to alter this behaviour.");
143  lwerror("transform: couldn't project point (%g %g %g): %s (%d)",
144  orig_pt.x, orig_pt.y, orig_pt.z, pj_strerrno(*pj_errno_ref), *pj_errno_ref);
145  return 0;
146  }
147  else
148  {
149  lwerror("transform: couldn't project point (%g %g %g): %s (%d)",
150  orig_pt.x, orig_pt.y, orig_pt.z, pj_strerrno(*pj_errno_ref), *pj_errno_ref);
151  return 0;
152  }
153  }
154 
155  if (pj_is_latlong(dstpj)) to_dec(pt);
156  return 1;
157 }
158 
159 projPJ
160 lwproj_from_string(const char *str1)
161 {
162  int t;
163  char *params[1024]; /* one for each parameter */
164  char *loc;
165  char *str;
166  size_t slen;
167  projPJ result;
168 
169 
170  if (str1 == NULL) return NULL;
171 
172  slen = strlen(str1);
173 
174  if (slen == 0) return NULL;
175 
176  str = lwalloc(slen+1);
177  strcpy(str, str1);
178 
179  /*
180  * first we split the string into a bunch of smaller strings,
181  * based on the " " separator
182  */
183 
184  params[0] = str; /* 1st param, we'll null terminate at the " " soon */
185 
186  loc = str;
187  t = 1;
188  while ((loc != NULL) && (*loc != 0) )
189  {
190  loc = strchr(loc, ' ');
191  if (loc != NULL)
192  {
193  *loc = 0; /* null terminate */
194  params[t] = loc+1;
195  loc++; /* next char */
196  t++; /*next param */
197  }
198  }
199 
200  if (!(result=pj_init(t, params)))
201  {
202  lwfree(str);
203  return NULL;
204  }
205  lwfree(str);
206  return result;
207 }
208 
209 
210 
void ptarray_set_point4d(POINTARRAY *pa, int n, const POINT4D *p4d)
Definition: lwgeom_api.c:501
double x
Definition: liblwgeom.h:308
#define LINETYPE
Definition: liblwgeom.h:61
static void to_dec(POINT4D *pt)
convert radians to decimal degress
#define MULTICURVETYPE
Definition: liblwgeom.h:70
void lwfree(void *mem)
Definition: lwutil.c:190
int npoints
Definition: liblwgeom.h:327
int lwgeom_transform(LWGEOM *geom, projPJ inpj, projPJ outpj)
Transform given SERIALIZED geometry from inpj projection to outpj projection.
#define POLYGONTYPE
Definition: liblwgeom.h:62
#define CURVEPOLYTYPE
Definition: liblwgeom.h:69
#define COMPOUNDTYPE
Definition: liblwgeom.h:68
#define MULTIPOINTTYPE
Definition: liblwgeom.h:63
#define LW_SUCCESS
Definition: liblwgeom.h:55
#define TRIANGLETYPE
Definition: liblwgeom.h:73
#define POLYHEDRALSURFACETYPE
Definition: liblwgeom.h:72
char ** result
Definition: liblwgeom.h:218
#define LW_FAILURE
Definition: liblwgeom.h:54
void lwerror(const char *fmt,...)
Write a notice out to the error handler.
Definition: lwutil.c:67
projPJ lwproj_from_string(const char *str1)
Get a projection from a string representation.
void lwnotice(const char *fmt,...)
Write a notice out to the notice handler.
Definition: lwutil.c:54
const char * lwtype_name(uint8_t type)
Return the type name string associated with a type number (e.g.
Definition: lwutil.c:164
int point4d_transform(POINT4D *pt, projPJ srcpj, projPJ dstpj)
LWGEOM ** geoms
Definition: liblwgeom.h:465
#define TINTYPE
Definition: liblwgeom.h:74
POINTARRAY ** rings
Definition: liblwgeom.h:413
int nrings
Definition: liblwgeom.h:411
double z
Definition: liblwgeom.h:308
#define MULTIPOLYGONTYPE
Definition: liblwgeom.h:65
#define MULTISURFACETYPE
Definition: liblwgeom.h:71
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition: liblwgeom.h:60
uint8_t type
Definition: liblwgeom.h:352
#define CIRCSTRINGTYPE
Definition: liblwgeom.h:67
void * lwalloc(size_t size)
Definition: lwutil.c:175
static void to_rad(POINT4D *pt)
convert decimal degress to radians
int lwgeom_is_empty(const LWGEOM *geom)
Return true or false depending on whether a geometry is an "empty" geometry (no vertices members) ...
Definition: lwgeom.c:1229
double y
Definition: liblwgeom.h:308
#define MULTILINETYPE
Definition: liblwgeom.h:64
#define LWDEBUGF(level, msg,...)
Definition: lwgeom_log.h:55
int getPoint4d_p(const POINTARRAY *pa, int n, POINT4D *point)
Definition: lwgeom_api.c:217
#define COLLECTIONTYPE
Definition: liblwgeom.h:66
This library is the generic geometry handling section of PostGIS.
int ptarray_transform(POINTARRAY *pa, projPJ inpj, projPJ outpj)
Transform given POINTARRAY from inpj projection to outpj projection.
POINTARRAY * points
Definition: liblwgeom.h:378