PostGIS  2.1.10dev-r@@SVN_REVISION@@
int pt_in_ring_3d ( const POINT3DZ p,
const POINTARRAY ring,
PLANE3D plane 
)

pt_in_ring_3d(): crossing number test for a point in a polygon input: p = a point, pa = vertex points of a ring V[n+1] with V[n]=V[0] plane=the plane that the vertex points are lying on returns: 0 = outside, 1 = inside

Our polygons have first and last point the same,

The difference in 3D variant is that we exclude the dimension that faces the plane least. That is the dimension with the highest number in pv

Definition at line 1027 of file measures3d.c.

References getPoint3dz_p(), LW_FALSE, LWDEBUGF, lwerror(), POINTARRAY::npoints, PLANE3D::pv, VECTOR3D::x, POINT3DZ::x, VECTOR3D::y, POINT3DZ::y, VECTOR3D::z, and POINT3DZ::z.

Referenced by lw_dist3d_pt_poly(), and lw_dist3d_ptarray_poly().

1028 {
1029 
1030  int cn = 0; /* the crossing number counter */
1031  int i;
1032  POINT3DZ v1, v2;
1033 
1034  POINT3DZ first, last;
1035 
1036  getPoint3dz_p(ring, 0, &first);
1037  getPoint3dz_p(ring, ring->npoints-1, &last);
1038  if ( memcmp(&first, &last, sizeof(POINT3DZ)) )
1039  {
1040  lwerror("pt_in_ring_3d: V[n] != V[0] (%g %g %g!= %g %g %g)",
1041  first.x, first.y, first.z, last.x, last.y, last.z);
1042  return LW_FALSE;
1043  }
1044 
1045  LWDEBUGF(2, "pt_in_ring_3d called with point: %g %g %g", p->x, p->y, p->z);
1046  /* printPA(ring); */
1047 
1048  /* loop through all edges of the polygon */
1049  getPoint3dz_p(ring, 0, &v1);
1050 
1051 
1052  if(fabs(plane->pv.z)>=fabs(plane->pv.x)&&fabs(plane->pv.z)>=fabs(plane->pv.y)) /*If the z vector of the normal vector to the plane is larger than x and y vector we project the ring to the xy-plane*/
1053  {
1054  for (i=0; i<ring->npoints-1; i++)
1055  {
1056  double vt;
1057  getPoint3dz_p(ring, i+1, &v2);
1058 
1059  /* edge from vertex i to vertex i+1 */
1060  if
1061  (
1062  /* an upward crossing */
1063  ((v1.y <= p->y) && (v2.y > p->y))
1064  /* a downward crossing */
1065  || ((v1.y > p->y) && (v2.y <= p->y))
1066  )
1067  {
1068 
1069  vt = (double)(p->y - v1.y) / (v2.y - v1.y);
1070 
1071  /* P.x <intersect */
1072  if (p->x < v1.x + vt * (v2.x - v1.x))
1073  {
1074  /* a valid crossing of y=p.y right of p.x */
1075  ++cn;
1076  }
1077  }
1078  v1 = v2;
1079  }
1080  }
1081  else if(fabs(plane->pv.y)>=fabs(plane->pv.x)&&fabs(plane->pv.y)>=fabs(plane->pv.z)) /*If the y vector of the normal vector to the plane is larger than x and z vector we project the ring to the xz-plane*/
1082  {
1083  for (i=0; i<ring->npoints-1; i++)
1084  {
1085  double vt;
1086  getPoint3dz_p(ring, i+1, &v2);
1087 
1088  /* edge from vertex i to vertex i+1 */
1089  if
1090  (
1091  /* an upward crossing */
1092  ((v1.z <= p->z) && (v2.z > p->z))
1093  /* a downward crossing */
1094  || ((v1.z > p->z) && (v2.z <= p->z))
1095  )
1096  {
1097 
1098  vt = (double)(p->z - v1.z) / (v2.z - v1.z);
1099 
1100  /* P.x <intersect */
1101  if (p->x < v1.x + vt * (v2.x - v1.x))
1102  {
1103  /* a valid crossing of y=p.y right of p.x */
1104  ++cn;
1105  }
1106  }
1107  v1 = v2;
1108  }
1109  }
1110  else /*Hopefully we only have the cases where x part of the normal vector is largest left*/
1111  {
1112  for (i=0; i<ring->npoints-1; i++)
1113  {
1114  double vt;
1115  getPoint3dz_p(ring, i+1, &v2);
1116 
1117  /* edge from vertex i to vertex i+1 */
1118  if
1119  (
1120  /* an upward crossing */
1121  ((v1.z <= p->z) && (v2.z > p->z))
1122  /* a downward crossing */
1123  || ((v1.z > p->z) && (v2.z <= p->z))
1124  )
1125  {
1126 
1127  vt = (double)(p->z - v1.z) / (v2.z - v1.z);
1128 
1129  /* P.x <intersect */
1130  if (p->y < v1.y + vt * (v2.y - v1.y))
1131  {
1132  /* a valid crossing of y=p.y right of p.x */
1133  ++cn;
1134  }
1135  }
1136  v1 = v2;
1137  }
1138  }
1139  LWDEBUGF(3, "pt_in_ring_3d returning %d", cn&1);
1140 
1141  return (cn&1); /* 0 if even (out), and 1 if odd (in) */
1142 }
double z
Definition: liblwgeom.h:290
double y
Definition: liblwgeom.h:290
double x
Definition: liblwgeom.h:290
double z
Definition: measures3d.h:36
int npoints
Definition: liblwgeom.h:327
VECTOR3D pv
Definition: measures3d.h:43
double y
Definition: measures3d.h:36
void lwerror(const char *fmt,...)
Write a notice out to the error handler.
Definition: lwutil.c:67
int getPoint3dz_p(const POINTARRAY *pa, int n, POINT3DZ *point)
Definition: lwgeom_api.c:305
#define LW_FALSE
Definition: liblwgeom.h:52
double x
Definition: measures3d.h:36
#define LWDEBUGF(level, msg,...)
Definition: lwgeom_log.h:55

Here is the call graph for this function:

Here is the caller graph for this function: