PostGIS  3.7.0dev-r@@SVN_REVISION@@

◆ ptarray_grid_in_place()

void ptarray_grid_in_place ( POINTARRAY pa,
gridspec grid 
)

Snap to grid.

Definition at line 2234 of file ptarray.c.

2235 {
2236  uint32_t j = 0;
2237  POINT4D *p, *p_out = NULL;
2238  double x, y, z = 0, m = 0;
2239  uint32_t ndims = FLAGS_NDIMS(pa->flags);
2240  uint32_t has_z = FLAGS_GET_Z(pa->flags);
2241  uint32_t has_m = FLAGS_GET_M(pa->flags);
2242 
2243  for (uint32_t i = 0; i < pa->npoints; i++)
2244  {
2245  /* Look straight into the abyss */
2246  p = (POINT4D *)(getPoint_internal(pa, i));
2247  x = p->x;
2248  y = p->y;
2249  if (ndims > 2)
2250  z = p->z;
2251  if (ndims > 3)
2252  m = p->m;
2253 
2254  /*
2255  * See https://github.com/libgeos/geos/pull/956
2256  * We use scale for rounding when gridsize is < 1 and
2257  * gridsize for rounding when scale < 1.
2258  */
2259  if (grid->xsize > 0) {
2260  if (grid->xsize < 1)
2261  x = rint((x - grid->ipx) * grid->xscale) / grid->xscale + grid->ipx;
2262  else
2263  x = rint((x - grid->ipx) / grid->xsize) * grid->xsize + grid->ipx;
2264  }
2265 
2266  if (grid->ysize > 0) {
2267  if (grid->ysize < 1)
2268  y = rint((y - grid->ipy) * grid->yscale) / grid->yscale + grid->ipy;
2269  else
2270  y = rint((y - grid->ipy) / grid->ysize) * grid->ysize + grid->ipy;
2271  }
2272 
2273  /* Read and round this point */
2274  /* Z is always in third position */
2275  if (has_z && grid->zsize > 0)
2276  z = rint((z - grid->ipz) / grid->zsize) * grid->zsize + grid->ipz;
2277 
2278  /* M might be in 3rd or 4th position */
2279  if (has_m && grid->msize > 0)
2280  {
2281  /* In POINT ZM, M is in 4th position, in POINT M, M is in 3rd position which is Z in POINT4D */
2282  if (has_z)
2283  m = rint((m - grid->ipm) / grid->msize) * grid->msize + grid->ipm;
2284  else
2285  z = rint((z - grid->ipm) / grid->msize) * grid->msize + grid->ipm;
2286  }
2287 
2288  /* Skip duplicates */
2289  if (p_out &&
2290  p_out->x == x &&
2291  p_out->y == y &&
2292  (ndims > 2 ? p_out->z == z : 1) &&
2293  (ndims > 3 ? p_out->m == m : 1))
2294  {
2295  continue;
2296  }
2297 
2298  /* Write rounded values into the next available point */
2299  p_out = (POINT4D *)(getPoint_internal(pa, j++));
2300  p_out->x = x;
2301  p_out->y = y;
2302  if (ndims > 2)
2303  p_out->z = z;
2304  if (ndims > 3)
2305  p_out->m = m;
2306  }
2307 
2308  /* Update output ptarray length */
2309  pa->npoints = j;
2310  return;
2311 }
#define FLAGS_GET_Z(flags)
Definition: liblwgeom.h:165
#define FLAGS_NDIMS(flags)
Definition: liblwgeom.h:179
#define FLAGS_GET_M(flags)
Definition: liblwgeom.h:166
static uint8_t * getPoint_internal(const POINTARRAY *pa, uint32_t n)
Definition: lwinline.h:75
double m
Definition: liblwgeom.h:414
double x
Definition: liblwgeom.h:414
double z
Definition: liblwgeom.h:414
double y
Definition: liblwgeom.h:414
lwflags_t flags
Definition: liblwgeom.h:431
uint32_t npoints
Definition: liblwgeom.h:427
double xscale
Definition: liblwgeom.h:1409
double ipm
Definition: liblwgeom.h:1404
double zsize
Definition: liblwgeom.h:1407
double ysize
Definition: liblwgeom.h:1406
double xsize
Definition: liblwgeom.h:1405
double yscale
Definition: liblwgeom.h:1410
double ipx
Definition: liblwgeom.h:1401
double msize
Definition: liblwgeom.h:1408
double ipy
Definition: liblwgeom.h:1402
double ipz
Definition: liblwgeom.h:1403

References POINTARRAY::flags, FLAGS_GET_M, FLAGS_GET_Z, FLAGS_NDIMS, getPoint_internal(), gridspec_t::ipm, gridspec_t::ipx, gridspec_t::ipy, gridspec_t::ipz, POINT4D::m, gridspec_t::msize, POINTARRAY::npoints, POINT4D::x, pixval::x, gridspec_t::xscale, gridspec_t::xsize, POINT4D::y, pixval::y, gridspec_t::yscale, gridspec_t::ysize, POINT4D::z, and gridspec_t::zsize.

Referenced by lwgeom_grid_in_place().

Here is the call graph for this function:
Here is the caller graph for this function: