/src/postgis/liblwgeom/lwinline.h

Source
/**********************************************************************
 *
 * PostGIS - Spatial Types for PostgreSQL
 * http://postgis.net
 *
 * PostGIS is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 2 of the License, or
 * (at your option) any later version.
 *
 * PostGIS is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with PostGIS.  If not, see <http://www.gnu.org/licenses/>.
 *
 **********************************************************************
 *
 * Copyright 2018 Darafei Praliaskouski <me@komzpa.net>
 * Copyright 2017-2018 Daniel Baston <dbaston@gmail.com>
 * Copyright 2011 Sandro Santilli <strk@kbt.io>
 * Copyright 2011 Paul Ramsey <pramsey@cleverelephant.ca>
 * Copyright 2007-2008 Mark Cave-Ayland
 * Copyright 2001-2006 Refractions Research Inc.
 *
 **********************************************************************/

#include <assert.h>

inline static double
distance2d_sqr_pt_pt(const POINT2D *p1, const POINT2D *p2)
{
  double hside = p2->x - p1->x;
  double vside = p2->y - p1->y;

  return hside * hside + vside * vside;
}

inline static double
distance3d_sqr_pt_pt(const POINT3D *p1, const POINT3D *p2)
{
  double hside = p2->x - p1->x;
  double vside = p2->y - p1->y;
  double zside = p2->z - p1->z;

  return hside * hside + vside * vside + zside * zside;
}

/*
 * Size of point represeneted in the POINTARRAY
 * 16 for 2d, 24 for 3d, 32 for 4d
 */
static inline size_t
ptarray_point_size(const POINTARRAY *pa)
{
  return sizeof(double) * FLAGS_NDIMS(pa->flags);
}

/*
 * Get a pointer to Nth point of a POINTARRAY
 * You'll need to cast it to appropriate dimensioned point.
 * Note that if you cast to a higher dimensional point you'll
 * possibly corrupt the POINTARRAY.
 *
 * Casting to returned pointer to POINT2D* should be safe,
 * as gserialized format always keeps the POINTARRAY pointer
 * aligned to double boundary.
 *
 * WARNING: Don't cast this to a POINT!
 * it would not be reliable due to memory alignment constraints
 */
static inline uint8_t *
getPoint_internal(const POINTARRAY *pa, uint32_t n)
{
  size_t size;
  uint8_t *ptr;

  assert(pa);
  assert(n <= pa->npoints);
  assert(n <= pa->maxpoints);

  size = ptarray_point_size(pa);
  ptr = pa->serialized_pointlist + size * n;

  return ptr;
}

/**
 * Returns a POINT2D pointer into the POINTARRAY serialized_ptlist,
 * suitable for reading from. This is very high performance
 * and declared const because you aren't allowed to muck with the
 * values, only read them.
 */
static inline const POINT2D *
getPoint2d_cp(const POINTARRAY *pa, uint32_t n)
{
  return (const POINT2D *)getPoint_internal(pa, n);
}

/**
 * Returns a POINT3D pointer into the POINTARRAY serialized_ptlist,
 * suitable for reading from. This is very high performance
 * and declared const because you aren't allowed to muck with the
 * values, only read them.
 */
static inline const POINT3D *
getPoint3d_cp(const POINTARRAY *pa, uint32_t n)
{
  return (const POINT3D *)getPoint_internal(pa, n);
}

/**
 * Returns a POINT4D pointer into the POINTARRAY serialized_ptlist,
 * suitable for reading from. This is very high performance
 * and declared const because you aren't allowed to muck with the
 * values, only read them.
 */
static inline const POINT4D *
getPoint4d_cp(const POINTARRAY *pa, uint32_t n)
{
  return (const POINT4D *)getPoint_internal(pa, n);
}

static inline LWPOINT *
lwgeom_as_lwpoint(const LWGEOM *lwgeom)
{
  if (!lwgeom)
    return NULL;
  if (lwgeom->type == POINTTYPE)
    return (LWPOINT *)lwgeom;
  else
    return NULL;
}

/**
 * Return LWTYPE number
 */
static inline uint32_t
lwgeom_get_type(const LWGEOM *geom)
{
  if (!geom)
    return 0;
  return geom->type;
}

static inline int
lwpoint_is_empty(const LWPOINT *point)
{
  return !point->point || point->point->npoints < 1;
}

static inline int
lwline_is_empty(const LWLINE *line)
{
  return !line->points || line->points->npoints < 1;
}

static inline int
lwcircstring_is_empty(const LWCIRCSTRING *circ)
{
  return !circ->points || circ->points->npoints < 1;
}

static inline int
lwpoly_is_empty(const LWPOLY *poly)
{
  return poly->nrings < 1 || !poly->rings || !poly->rings[0] || poly->rings[0]->npoints < 1;
}

static inline int
lwtriangle_is_empty(const LWTRIANGLE *triangle)
{
  return !triangle->points || triangle->points->npoints < 1;
}

static inline int lwgeom_is_empty(const LWGEOM *geom);

/**
 * Determine whether a collection geometry contains no non-empty geometries.
 *
 * Returns true if the collection has zero geometries, its geometry array is NULL,
 * or every contained geometry is empty according to lwgeom_is_empty(). Returns
 * false as soon as any contained geometry is found to be non-empty.
 *
 * @param col Collection to check; must be a valid LWCOLLECTION pointer.
 * @return LW_TRUE if the collection is empty, LW_FALSE otherwise.
 */
static inline int
lwcollection_is_empty(const LWCOLLECTION *col)
{
  uint32_t i;
  if (col->ngeoms == 0 || !col->geoms)
    return LW_TRUE;
  for (i = 0; i < col->ngeoms; i++)
  {
    if (!lwgeom_is_empty(col->geoms[i]))
      return LW_FALSE;
  }
  return LW_TRUE;
}

/**
 * Determine whether a NURBS curve is empty.
 *
 * A NURBS curve is considered empty if the curve pointer is NULL, its
 * associated points array is NULL, or the points array reports zero control
 * points.
 *
 * @param curve NURBS curve to test.
 * @return LW_TRUE if the curve is empty, otherwise LW_FALSE.
 */
static inline int
lwnurbscurve_is_empty(const LWNURBSCURVE *curve)
{
    if (!curve)
        return LW_TRUE;

    if (!curve->points)
        return LW_TRUE;

    if (curve->points->npoints == 0)
        return LW_TRUE;

    return LW_FALSE;
}

/**
 * Determine whether a geometry contains no vertices.
 *
 * Returns LW_TRUE if the provided LWGEOM has no coordinate vertices according to
 * its concrete type (points, lines, rings, polygons, collections, NURBS curves, etc.).
 * The check is delegated to type-specific emptiness helpers; unrecognized geometry
 * types are treated as non-empty (returns LW_FALSE).
 *
 * The caller must pass a non-NULL geom pointer whose ->type field is valid.
 *
 * @param geom Geometry to test for emptiness.
 * @return LW_TRUE if the geometry is empty, otherwise LW_FALSE.
 */
static inline int
lwgeom_is_empty(const LWGEOM *geom)
{
  switch (geom->type)
  {
  case POINTTYPE:
    return lwpoint_is_empty((LWPOINT *)geom);
    break;
  case LINETYPE:
    return lwline_is_empty((LWLINE *)geom);
    break;
  case CIRCSTRINGTYPE:
    return lwcircstring_is_empty((LWCIRCSTRING *)geom);
    break;
  case POLYGONTYPE:
    return lwpoly_is_empty((LWPOLY *)geom);
    break;
  case TRIANGLETYPE:
    return lwtriangle_is_empty((LWTRIANGLE *)geom);
    break;
  case MULTIPOINTTYPE:
  case MULTILINETYPE:
  case MULTIPOLYGONTYPE:
  case COMPOUNDTYPE:
  case CURVEPOLYTYPE:
  case MULTICURVETYPE:
  case MULTISURFACETYPE:
  case POLYHEDRALSURFACETYPE:
  case TINTYPE:
  case COLLECTIONTYPE:
    return lwcollection_is_empty((LWCOLLECTION *)geom);
    break;
  case NURBSCURVETYPE:
    return lwnurbscurve_is_empty((LWNURBSCURVE*)geom);
    break;
  default:
    return LW_FALSE;
    break;
  }
}

inline static uint64_t
uint64_interleave_2(uint64_t x, uint64_t y)
{
  x = (x | (x << 16)) & 0x0000FFFF0000FFFFULL;
  x = (x | (x << 8)) & 0x00FF00FF00FF00FFULL;
  x = (x | (x << 4)) & 0x0F0F0F0F0F0F0F0FULL;
  x = (x | (x << 2)) & 0x3333333333333333ULL;
  x = (x | (x << 1)) & 0x5555555555555555ULL;

  y = (y | (y << 16)) & 0x0000FFFF0000FFFFULL;
  y = (y | (y << 8)) & 0x00FF00FF00FF00FFULL;
  y = (y | (y << 4)) & 0x0F0F0F0F0F0F0F0FULL;
  y = (y | (y << 2)) & 0x3333333333333333ULL;
  y = (y | (y << 1)) & 0x5555555555555555ULL;

  return x | (y << 1);
}

/* Based on https://github.com/rawrunprotected/hilbert_curves Public Domain code */
inline static uint64_t
uint32_hilbert(uint32_t px, uint32_t py)
{
  uint64_t x = px;
  uint64_t y = py;

  uint64_t A, B, C, D;
  uint64_t a, b, c, d;
  uint64_t i0, i1;

  // Initial prefix scan round, prime with x and y
  {
    a = x ^ y;
    b = 0xFFFFFFFFULL ^ a;
    c = 0xFFFFFFFFULL ^ (x | y);
    d = x & (y ^ 0xFFFFFFFFULL);

    A = a | (b >> 1);
    B = (a >> 1) ^ a;
    C = ((c >> 1) ^ (b & (d >> 1))) ^ c;
    D = ((a & (c >> 1)) ^ (d >> 1)) ^ d;
  }

  {
    a = A;
    b = B;
    c = C;
    d = D;

    A = ((a & (a >> 2)) ^ (b & (b >> 2)));
    B = ((a & (b >> 2)) ^ (b & ((a ^ b) >> 2)));
    C ^= ((a & (c >> 2)) ^ (b & (d >> 2)));
    D ^= ((b & (c >> 2)) ^ ((a ^ b) & (d >> 2)));
  }

  {
    a = A;
    b = B;
    c = C;
    d = D;

    A = ((a & (a >> 4)) ^ (b & (b >> 4)));
    B = ((a & (b >> 4)) ^ (b & ((a ^ b) >> 4)));
    C ^= ((a & (c >> 4)) ^ (b & (d >> 4)));
    D ^= ((b & (c >> 4)) ^ ((a ^ b) & (d >> 4)));
  }

  {
    a = A;
    b = B;
    c = C;
    d = D;

    A = ((a & (a >> 8)) ^ (b & (b >> 8)));
    B = ((a & (b >> 8)) ^ (b & ((a ^ b) >> 8)));
    C ^= ((a & (c >> 8)) ^ (b & (d >> 8)));
    D ^= ((b & (c >> 8)) ^ ((a ^ b) & (d >> 8)));
  }

  {
    a = A;
    b = B;
    c = C;
    d = D;

    C ^= ((a & (c >> 16)) ^ (b & (d >> 16)));
    D ^= ((b & (c >> 16)) ^ ((a ^ b) & (d >> 16)));
  }

  // Undo transformation prefix scan
  a = C ^ (C >> 1);
  b = D ^ (D >> 1);

  // Recover index bits
  i0 = x ^ y;
  i1 = b | (0xFFFFFFFFULL ^ (i0 | a));

  return uint64_interleave_2(i0, i1);
}

/*
 * This macro is based on PG_FREE_IF_COPY, except that it accepts two pointers.
 * See PG_FREE_IF_COPY comment in src/include/fmgr.h in postgres source code
 * for more details.
 */
#define POSTGIS_FREE_IF_COPY_P(ptrsrc, ptrori) \
  do \
  { \
    if ((Pointer)(ptrsrc) != (Pointer)(ptrori)) \
      pfree(ptrsrc); \
  } while (0)

Coverage Report

Created: 2026-09-01 06:53