PostGIS  2.2.8dev-r@@SVN_REVISION@@

◆ gserialized_distance_nd()

Datum gserialized_distance_nd ( PG_FUNCTION_ARGS  )

Definition at line 668 of file gserialized_gist_nd.c.

References distance(), gidx_distance_leaf_centroid(), gserialized_within(), LINETYPE, LW_SUCCESS, lwgeom_as_lwline(), lwgeom_closest_line(), lwgeom_closest_line_3d(), lwgeom_free(), lwgeom_from_gserialized(), lwgeom_get_type(), lwgeom_has_m(), lwgeom_has_z(), lwgeom_interpolate_point(), lwgeom_length(), lwgeom_length_2d(), lwline_get_lwpoint(), lwpoint_free(), lwpoint_getPoint4d_p(), POINT4D::m, PG_FUNCTION_INFO_V1(), and POINTTYPE.

Referenced by gserialized_expand().

669 {
670  char b1mem[GIDX_MAX_SIZE];
671  GIDX *b1 = (GIDX*)b1mem;
672  char b2mem[GIDX_MAX_SIZE];
673  GIDX *b2 = (GIDX*)b2mem;
674 
675 #if POSTGIS_PGSQL_VERSION < 95
676 
677  /* Centroid-to-centroid distance */
678  Datum gs1 = PG_GETARG_DATUM(0);
679  Datum gs2 = PG_GETARG_DATUM(1);
680  double box_distance = FLT_MAX;
681 
682  /* Must be able to build box for each argument (ie, not empty geometry). */
683  if ( (gserialized_datum_get_gidx_p(gs1, b1) == LW_SUCCESS) &&
684  (gserialized_datum_get_gidx_p(gs2, b2) == LW_SUCCESS) )
685  {
686  box_distance = gidx_distance_leaf_centroid(b1, b2);
687  POSTGIS_DEBUGF(3, "got boxes %s and %s", gidx_to_string(b1), gidx_to_string(b2));
688  }
689  PG_RETURN_FLOAT8(box_distance);
690 
691 #else /* POSTGIS_PGSQL_VERSION >= 96 */
692 
693  /* Feature-to-feature distance */
694  GSERIALIZED *geom1 = PG_GETARG_GSERIALIZED_P(0);
695  GSERIALIZED *geom2 = PG_GETARG_GSERIALIZED_P(1);
696  LWGEOM *lw1 = lwgeom_from_gserialized(geom1);
697  LWGEOM *lw2 = lwgeom_from_gserialized(geom2);
698  LWGEOM *closest;
699  double distance;
700 
701 
702  /* Find an exact shortest line w/ the dimensions we support */
703  if ( lwgeom_has_z(lw1) && lwgeom_has_z(lw2) )
704  {
705  closest = lwgeom_closest_line_3d(lw1, lw2);
706  distance = lwgeom_length(closest);
707  }
708  else
709  {
710  closest = lwgeom_closest_line(lw1, lw2);
711  distance = lwgeom_length_2d(closest);
712  }
713 
714  /* Un-sqrt the distance so we can add extra terms */
715  distance = distance*distance;
716 
717  /* Can only add the M term if both objects have M */
718  if ( lwgeom_has_m(lw1) && lwgeom_has_m(lw2) )
719  {
720  double m1, m2;
721  int usebox = false;
722 
723  if ( lwgeom_get_type(lw1) == POINTTYPE )
724  {
725  POINT4D p;
726  lwpoint_getPoint4d_p((LWPOINT*)lw1, &p);
727  m1 = p.m;
728  }
729  else if ( lwgeom_get_type(lw1) == LINETYPE )
730  {
731  LWPOINT *lwp1 = lwline_get_lwpoint(lwgeom_as_lwline(closest), 0);
732  m1 = lwgeom_interpolate_point(lw1, lwp1);
733  lwpoint_free(lwp1);
734  }
735  else
736  {
737  usebox = true;
738  }
739 
740  if ( lwgeom_get_type(lw2) == POINTTYPE )
741  {
742  POINT4D p;
743  lwpoint_getPoint4d_p((LWPOINT*)lw2, &p);
744  m2 = p.m;
745  }
746  else if ( lwgeom_get_type(lw2) == LINETYPE )
747  {
748  LWPOINT *lwp2 = lwline_get_lwpoint(lwgeom_as_lwline(closest), 1);
749  m2 = lwgeom_interpolate_point(lw2, lwp2);
750  lwpoint_free(lwp2);
751  }
752  else
753  {
754  usebox = true;
755  }
756 
757  if ( usebox )
758  {
759  double d;
760  gserialized_get_gidx_p(geom1, b1);
761  gserialized_get_gidx_p(geom2, b2);
762  d = gidx_distance_m(b1, b2);
763  distance += d*d;
764  }
765  else
766  {
767  distance += (m2-m1)*(m2-m1);
768  }
769  }
770 
771  lwgeom_free(closest);
772 
773  PG_FREE_IF_COPY(geom1, 0);
774  PG_FREE_IF_COPY(geom2, 1);
775  PG_RETURN_FLOAT8(sqrt(distance));
776 #endif /* POSTGIS_PGSQL_VERSION >= 96 */
777 }
#define LINETYPE
Definition: liblwgeom.h:71
double m
Definition: liblwgeom.h:336
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
uint32_t lwgeom_get_type(const LWGEOM *geom)
Return LWTYPE number.
Definition: lwgeom.c:829
void lwpoint_free(LWPOINT *pt)
Definition: lwpoint.c:182
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1050
#define LW_SUCCESS
Definition: liblwgeom.h:65
LWGEOM * lwgeom_closest_line_3d(const LWGEOM *lw1, const LWGEOM *lw2)
Definition: measures3d.c:64
double lwgeom_length_2d(const LWGEOM *geom)
Definition: lwgeom.c:1668
int lwgeom_has_z(const LWGEOM *geom)
Return LW_TRUE if geometry has Z ordinates.
Definition: lwgeom.c:836
double lwgeom_interpolate_point(const LWGEOM *lwin, const LWPOINT *lwpt)
Find the measure value at the location on the line closest to the point.
static double gidx_distance_leaf_centroid(const GIDX *a, const GIDX *b)
Calculate the centroid->centroid distance between the boxes.
int lwpoint_getPoint4d_p(const LWPOINT *point, POINT4D *out)
Definition: lwpoint.c:44
Datum distance(PG_FUNCTION_ARGS)
LWLINE * lwgeom_as_lwline(const LWGEOM *lwgeom)
Definition: lwgeom.c:89
double lwgeom_length(const LWGEOM *geom)
Definition: lwgeom.c:1646
#define POINTTYPE
LWTYPE numbers, used internally by PostGIS.
Definition: liblwgeom.h:70
LWPOINT * lwline_get_lwpoint(const LWLINE *line, int where)
Returns freshly allocated LWPOINT that corresponds to the index where.
Definition: lwline.c:295
LWGEOM * lwgeom_closest_line(const LWGEOM *lw1, const LWGEOM *lw2)
Definition: measures.c:28
int lwgeom_has_m(const LWGEOM *geom)
Return LW_TRUE if geometry has M ordinates.
Definition: lwgeom.c:843
Here is the call graph for this function:
Here is the caller graph for this function: