PostGIS  2.5.0dev-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  * PostGIS is free software: you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation, either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * PostGIS is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with PostGIS. If not, see <http://www.gnu.org/licenses/>.
18  *
19  **********************************************************************
20  *
21  * Copyright (C) 2001-2003 Refractions Research Inc.
22  *
23  **********************************************************************/
24 
25 
26 #include "../postgis_config.h"
27 #include "liblwgeom_internal.h"
28 #include "lwgeom_log.h"
29 #include <string.h>
30 
31 
33 static void
35 {
36  pt->x *= M_PI/180.0;
37  pt->y *= M_PI/180.0;
38 }
39 
41 static void
43 {
44  pt->x *= 180.0/M_PI;
45  pt->y *= 180.0/M_PI;
46 }
47 
52 int
53 ptarray_transform(POINTARRAY *pa, projPJ inpj, projPJ outpj)
54 {
55  uint32_t i;
56  POINT4D p;
57 
58  for ( i = 0; i < pa->npoints; i++ )
59  {
60  getPoint4d_p(pa, i, &p);
61  if ( ! point4d_transform(&p, inpj, outpj) ) return LW_FAILURE;
62  ptarray_set_point4d(pa, i, &p);
63  }
64 
65  return LW_SUCCESS;
66 }
67 
68 
73 int
74 lwgeom_transform(LWGEOM *geom, projPJ inpj, projPJ outpj)
75 {
76  uint32_t i;
77 
78  /* No points to transform in an empty! */
79  if ( lwgeom_is_empty(geom) )
80  return LW_SUCCESS;
81 
82  switch(geom->type)
83  {
84  case POINTTYPE:
85  case LINETYPE:
86  case CIRCSTRINGTYPE:
87  case TRIANGLETYPE:
88  {
89  LWLINE *g = (LWLINE*)geom;
90  if ( ! ptarray_transform(g->points, inpj, outpj) ) return LW_FAILURE;
91  break;
92  }
93  case POLYGONTYPE:
94  {
95  LWPOLY *g = (LWPOLY*)geom;
96  for ( i = 0; i < g->nrings; i++ )
97  {
98  if ( ! ptarray_transform(g->rings[i], inpj, outpj) ) return LW_FAILURE;
99  }
100  break;
101  }
102  case MULTIPOINTTYPE:
103  case MULTILINETYPE:
104  case MULTIPOLYGONTYPE:
105  case COLLECTIONTYPE:
106  case COMPOUNDTYPE:
107  case CURVEPOLYTYPE:
108  case MULTICURVETYPE:
109  case MULTISURFACETYPE:
111  case TINTYPE:
112  {
113  LWCOLLECTION *g = (LWCOLLECTION*)geom;
114  for ( i = 0; i < g->ngeoms; i++ )
115  {
116  if ( ! lwgeom_transform(g->geoms[i], inpj, outpj) ) return LW_FAILURE;
117  }
118  break;
119  }
120  default:
121  {
122  lwerror("lwgeom_transform: Cannot handle type '%s'",
123  lwtype_name(geom->type));
124  return LW_FAILURE;
125  }
126  }
127  return LW_SUCCESS;
128 }
129 
130 int
131 point4d_transform(POINT4D *pt, projPJ srcpj, projPJ dstpj)
132 {
133  POINT3D orig_pt = {pt->x, pt->y, pt->z}; /* Copy for error report*/
134 
135  if (pj_is_latlong(srcpj)) to_rad(pt) ;
136 
137  LWDEBUGF(4, "transforming POINT(%f %f) from '%s' to '%s'",
138  orig_pt.x, orig_pt.y, pj_get_def(srcpj,0), pj_get_def(dstpj,0));
139 
140  if (pj_transform(srcpj, dstpj, 1, 0, &(pt->x), &(pt->y), &(pt->z)) != 0)
141  {
142  int pj_errno_val = *pj_get_errno_ref();
143  if (pj_errno_val == -38)
144  {
145  lwnotice("PostGIS was unable to transform the point because either no grid"
146  " shift files were found, or the point does not lie within the "
147  "range for which the grid shift is defined. Refer to the "
148  "ST_Transform() section of the PostGIS manual for details on how "
149  "to configure PostGIS to alter this behaviour.");
150  lwerror("transform: couldn't project point (%g %g %g): %s (%d)",
151  orig_pt.x, orig_pt.y, orig_pt.z,
152  pj_strerrno(pj_errno_val), pj_errno_val);
153  }
154  else
155  {
156  lwerror("transform: couldn't project point (%g %g %g): %s (%d)",
157  orig_pt.x, orig_pt.y, orig_pt.z,
158  pj_strerrno(pj_errno_val), pj_errno_val);
159  }
160  return LW_FAILURE;
161  }
162 
163  if (pj_is_latlong(dstpj)) to_dec(pt);
164  return LW_SUCCESS;
165 }
166 
167 projPJ
168 lwproj_from_string(const char *str1)
169 {
170  if (!str1 || str1[0] == '\0')
171  {
172  return NULL;
173  }
174  return pj_init_plus(str1);
175 }
176 
177 
178 
double x
Definition: liblwgeom.h:354
#define LINETYPE
Definition: liblwgeom.h:85
static void to_dec(POINT4D *pt)
convert radians to decimal degress
#define MULTICURVETYPE
Definition: liblwgeom.h:94
void lwnotice(const char *fmt,...)
Write a notice out to the notice handler.
Definition: lwutil.c:177
double y
Definition: liblwgeom.h:342
int lwgeom_transform(LWGEOM *geom, projPJ inpj, projPJ outpj)
Transform given SERIALIZED geometry from inpj projection to outpj projection.
#define POLYGONTYPE
Definition: liblwgeom.h:86
#define CURVEPOLYTYPE
Definition: liblwgeom.h:93
#define COMPOUNDTYPE
Definition: liblwgeom.h:92
#define MULTIPOINTTYPE
Definition: liblwgeom.h:87
#define LW_SUCCESS
Definition: liblwgeom.h:79
void ptarray_set_point4d(POINTARRAY *pa, uint32_t n, const POINT4D *p4d)
Definition: lwgeom_api.c:425
#define TRIANGLETYPE
Definition: liblwgeom.h:97
double x
Definition: liblwgeom.h:342
#define POLYHEDRALSURFACETYPE
Definition: liblwgeom.h:96
uint32_t ngeoms
Definition: liblwgeom.h:509
uint32_t nrings
Definition: liblwgeom.h:457
#define LW_FAILURE
Definition: liblwgeom.h:78
double z
Definition: liblwgeom.h:342
unsigned int uint32_t
Definition: uthash.h:78
projPJ lwproj_from_string(const char *str1)
Get a projection from a string representation.
const char * lwtype_name(uint8_t type)
Return the type name string associated with a type number (e.g.
Definition: lwutil.c:218
int point4d_transform(POINT4D *pt, projPJ srcpj, projPJ dstpj)
LWGEOM ** geoms
Definition: liblwgeom.h:511
#define TINTYPE
Definition: liblwgeom.h:98
POINTARRAY ** rings
Definition: liblwgeom.h:459
double z
Definition: liblwgeom.h:354
int getPoint4d_p(const POINTARRAY *pa, uint32_t n, POINT4D *point)
Definition: lwgeom_api.c:113
#define MULTIPOLYGONTYPE
Definition: liblwgeom.h:89
#define MULTISURFACETYPE
Definition: liblwgeom.h:95
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition: liblwgeom.h:84
uint8_t type
Definition: liblwgeom.h:398
#define CIRCSTRINGTYPE
Definition: liblwgeom.h:91
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:1393
double y
Definition: liblwgeom.h:354
#define MULTILINETYPE
Definition: liblwgeom.h:88
#define LWDEBUGF(level, msg,...)
Definition: lwgeom_log.h:88
void lwerror(const char *fmt,...)
Write a notice out to the error handler.
Definition: lwutil.c:190
#define COLLECTIONTYPE
Definition: liblwgeom.h:90
int ptarray_transform(POINTARRAY *pa, projPJ inpj, projPJ outpj)
Transform given POINTARRAY from inpj projection to outpj projection.
POINTARRAY * points
Definition: liblwgeom.h:424
uint32_t npoints
Definition: liblwgeom.h:373