Changeset: 5398aea38dae for MonetDB
URL: http://dev.monetdb.org/hg/MonetDB?cmd=changeset;node=5398aea38dae
Modified Files:
geom/monetdb5/geomPoints.c
Branch: geo
Log Message:
contains: reading coordSeq in array
diffs (265 lines):
diff --git a/geom/monetdb5/geomPoints.c b/geom/monetdb5/geomPoints.c
--- a/geom/monetdb5/geomPoints.c
+++ b/geom/monetdb5/geomPoints.c
@@ -266,6 +266,8 @@ static str pnpoly_(int *out, const GEOSG
const GEOSCoordSequence *coordSeq;
unsigned int geometryPointsNum=0 ;
+ double *xPoints, *yPoints; //arrays with the points of the ring
+
/* get the coordinates of the points comprising the geometry */
if(!(coordSeq = GEOSGeom_getCoordSeq(geosGeometry)))
@@ -274,6 +276,10 @@ static str pnpoly_(int *out, const GEOSG
/* get the number of points in the geometry */
GEOSCoordSeq_getSize(coordSeq, &geometryPointsNum);
+ //allocate space for the coordinates
+ xPoints = (double*)GDKmalloc(sizeof(double*)*geometryPointsNum);
+ yPoints = (double*)GDKmalloc(sizeof(double*)*geometryPointsNum);
+
/*Get the BATs*/
if ((bpx = BATdescriptor(*point_x)) == NULL) {
throw(MAL, "geom.point", RUNTIME_OBJECT_MISSING);
@@ -305,22 +311,35 @@ static str pnpoly_(int *out, const GEOSG
// gettimeofday(&start, NULL);
cs = (bte*) Tloc(bo,BUNfirst(bo));
- for (i = 0; i < BATcount(bpx); i++) { //for each point in the x, y BATs
- int wn = 0;
- for (j = 0; j < geometryPointsNum-1; j++) { //check each point
in the geometry (the exteriorRing)
- double xCurrent=0.0, yCurrent=0.0, xNext=0.0, yNext=0.0;
- if(GEOSCoordSeq_getX(coordSeq, j, &xCurrent) == -1 ||
GEOSCoordSeq_getX(coordSeq, j+1, &xNext) == -1)
+ for (i = 0; i < BATcount(bpx); i++) {
+ int wn = 0;
+
+ if(i==0) {
+ //read the first point from the external ring
+ if(GEOSCoordSeq_getX(coordSeq, j, &xPoints[0]) == -1)
return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getX failed");
- if(GEOSCoordSeq_getY(coordSeq, j, &yCurrent) == -1 ||
GEOSCoordSeq_getY(coordSeq, j+1, &yNext) == -1)
+ if(GEOSCoordSeq_getY(coordSeq, j, &yPoints[0]) == -1)
return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getY failed");
- //the edge goes from small y to big y (upward
direction) and the point is somewhere in there
- if (yCurrent <= py[i] && yNext >= py[i]) {
- wn+=isLeft( xCurrent, yCurrent, xNext, yNext,
px[i], py[i]);
+ }
+ /*Check if the point is inside the polygon*/
+ for (j = 1; j < geometryPointsNum; j++) { //check each point in
the exterior ring)
+ //double xCurrent=0.0, yCurrent=0.0, xNext=0.0,
yNext=0.0;
+ if( i == 0) {
+ //it is the first iteration, I need to read the
points
+ if(GEOSCoordSeq_getX(coordSeq, j, &xPoints[j])
== -1)
+ return createException(MAL,
"batgeom.Contains", "GEOSCoordSeq_getX failed");
+ if(GEOSCoordSeq_getY(coordSeq, j, &yPoints[j])
== -1)
+ return createException(MAL,
"batgeom.Contains", "GEOSCoordSeq_getY failed");
+
+ }
+ //the edge goes from small y to big y (upward
direction) and the point is somewhere in there
+ if (yPoints[j-1] <= py[i] && yPoints[j] >= py[i]) {
+ wn+=isLeft(xPoints[j-1], yPoints[j-1],
xPoints[j], yPoints[j], px[i], py[i]);
}
//the edge goes from big y to small y (downward
direction) and the point is somewhere in there
- else if (yCurrent >= py[i] && yNext <= py[i]) {
- wn+=isRight( xCurrent, yCurrent, xNext, yNext,
px[i], py[i]);
- }
+ else if (yPoints[j-1] >= py[i] && yPoints[j] <= py[i]) {
+ wn+=isRight(xPoints[j-1], yPoints[j-1],
xPoints[j], yPoints[j], px[i], py[i]);
+ }
}
*cs++ = wn & 1;
}
@@ -338,6 +357,10 @@ static str pnpoly_(int *out, const GEOSG
BBPreleaseref(bpx->batCacheid);
BBPreleaseref(bpy->batCacheid);
BBPkeepref(*out = bo->batCacheid);
+
+ GDKfree(xPoints);
+ GDKfree(yPoints);
+
return MAL_SUCCEED;
}
@@ -351,7 +374,9 @@ static str pnpolyWithHoles_(int *out, GE
const GEOSGeometry *exteriorRingGeometry;
const GEOSCoordSequence *exteriorRingCoordSeq;
- unsigned int exteriorRingPointsNum=0 ;
+ double **xPoints, **yPoints; //arrays with the points of the rings so
that we do not read them every time
+ unsigned int *pointsNum; //array with the number of points in each ring
+ bool checked = false; //used to know when the internal rings have been
checked
/*Get the BATs*/
if ((bpx = BATdescriptor(*point_x)) == NULL) {
@@ -379,6 +404,10 @@ static str pnpolyWithHoles_(int *out, GE
throw(MAL, "geom.point", MAL_MALLOC_FAIL);
}
BATseqbase(bo, bpx->hseqbase);
+
+ xPoints = (double**)GDKmalloc(sizeof(double**)*(interiorRingsNum+1));
+ yPoints = (double**)GDKmalloc(sizeof(double**)*(interiorRingsNum+1));
+ pointsNum = (unsigned int*)GDKmalloc(sizeof(unsigned
int*)*(interiorRingsNum+1));
//get the exterior ring of the geometry
if(!(exteriorRingGeometry = GEOSGetExteriorRing(geosGeometry)))
@@ -389,7 +418,11 @@ static str pnpolyWithHoles_(int *out, GE
return createException(MAL, "batgeom.Contains",
"GEOSGeom_getCoordSeq failed");
/* get the number of points in the exterior ring */
- GEOSCoordSeq_getSize(exteriorRingCoordSeq, &exteriorRingPointsNum);
+ GEOSCoordSeq_getSize(exteriorRingCoordSeq, &pointsNum[0]);
+
+ //the exteriorRing is in position 0
+ xPoints[0] = (double*)GDKmalloc(sizeof(double*)*pointsNum[0]);
+ yPoints[0] = (double*)GDKmalloc(sizeof(double*)*pointsNum[0]);
/*Iterate over the Point BATs and determine if they are in Polygon
represented by vertex BATs*/
px = (dbl *) Tloc(bpx, BUNfirst(bpx));
@@ -398,20 +431,31 @@ static str pnpolyWithHoles_(int *out, GE
for (i = 0; i < BATcount(bpx); i++) {
int wn = 0;
- /*Check if the point is insode the polygon*/
- for (j = 0; j < exteriorRingPointsNum-1; j++) { //check each point in
the exterior ring)
- double xCurrent=0.0, yCurrent=0.0, xNext=0.0, yNext=0.0;
- if(GEOSCoordSeq_getX(exteriorRingCoordSeq, j, &xCurrent) == -1
|| GEOSCoordSeq_getX(exteriorRingCoordSeq, j+1, &xNext) == -1)
+ if(i==0 && pointsNum[0]>0) {
+ //read the first point from the external ring
+ if(GEOSCoordSeq_getX(exteriorRingCoordSeq, 0, &xPoints[0][0])
== -1)
return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getX failed");
- if(GEOSCoordSeq_getY(exteriorRingCoordSeq, j, &yCurrent) == -1
|| GEOSCoordSeq_getY(exteriorRingCoordSeq, j+1, &yNext) == -1)
+ if(GEOSCoordSeq_getY(exteriorRingCoordSeq, 0, &yPoints[0][0])
== -1)
return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getY failed");
+ }
+ /*Check if the point is inside the polygon*/
+ for (j = 1; j < pointsNum[0]; j++) { //check each point in the exterior
ring)
+ //double xCurrent=0.0, yCurrent=0.0, xNext=0.0, yNext=0.0;
+ if( i == 0) {
+ //it is the first iteration, I need to read the points
+ if(GEOSCoordSeq_getX(exteriorRingCoordSeq, j,
&xPoints[0][j]) == -1)
+ return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getX failed");
+ if(GEOSCoordSeq_getY(exteriorRingCoordSeq, j,
&yPoints[0][j]) == -1)
+ return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getY failed");
+
+ }
//the edge goes from small y to big y (upward direction) and
the point is somewhere in there
- if (yCurrent <= py[i] && yNext >= py[i]) {
- wn+=isLeft( xCurrent, yCurrent, xNext, yNext, px[i],
py[i]);
+ if (yPoints[0][j-1] <= py[i] && yPoints[0][j] >= py[i]) {
+ wn+=isLeft(xPoints[0][j-1], yPoints[0][j-1],
xPoints[0][j], yPoints[0][j], px[i], py[i]);
}
//the edge goes from big y to small y (downward direction) and
the point is somewhere in there
- else if (yCurrent >= py[i] && yNext <= py[i]) {
- wn+=isRight( xCurrent, yCurrent, xNext, yNext, px[i],
py[i]);
+ else if (yPoints[0][j-1] >= py[i] && yPoints[0][j] <= py[i]) {
+ wn+=isRight(xPoints[0][j-1], yPoints[0][j-1],
xPoints[0][j], yPoints[0][j], px[i], py[i]);
}
}
@@ -421,38 +465,51 @@ static str pnpolyWithHoles_(int *out, GE
if(!wn)
continue;
- //Inside the polygon, check the holes
+ //inside the polygon, check the holes
for (h = 0; h < interiorRingsNum; h++) {
- const GEOSGeometry *interiorRingGeometry;
const GEOSCoordSequence *interiorRingCoordSeq;
- unsigned int interiorRingPointsNum=0 ;
- //get the interior ring
- if(!(interiorRingGeometry = GEOSGetInteriorRingN(geosGeometry,
h)))
- return createException(MAL, "batgeom.Contains",
"GEOSGetInteriorRingN failed");
+ if(!checked) {
+ const GEOSGeometry *interiorRingGeometry;
+
+ //get the interior ring
+ if(!(interiorRingGeometry =
GEOSGetInteriorRingN(geosGeometry, h)))
+ return createException(MAL, "batgeom.Contains",
"GEOSGetInteriorRingN failed");
- /* get the coordinates of the points comprising the interior
ring */
- if(!(interiorRingCoordSeq =
GEOSGeom_getCoordSeq(interiorRingGeometry)))
- return createException(MAL, "batgeom.Contains",
"GEOSGeom_getCoordSeq failed");
+ // get the coordinates of the points comprising the
interior ring
+ if(!(interiorRingCoordSeq =
GEOSGeom_getCoordSeq(interiorRingGeometry)))
+ return createException(MAL, "batgeom.Contains",
"GEOSGeom_getCoordSeq failed");
- /* get the number of points in the interior ring */
- GEOSCoordSeq_getSize(interiorRingCoordSeq,
&interiorRingPointsNum);
-
+ // get the number of points in the interior ring
+ GEOSCoordSeq_getSize(interiorRingCoordSeq,
&pointsNum[h+1]);
+
+ xPoints[h+1] =
(double*)GDKmalloc(sizeof(double*)*pointsNum[h+1]);
+ yPoints[h+1] =
(double*)GDKmalloc(sizeof(double*)*pointsNum[h+1]);
+
+ //read the first point from the ring
+ if(GEOSCoordSeq_getX(interiorRingCoordSeq, 0,
&xPoints[h+1][0]) == -1)
+ return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getX failed");
+ if(GEOSCoordSeq_getY(interiorRingCoordSeq, 0,
&yPoints[h+1][0]) == -1)
+ return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getY failed");
+ }
+
wn = 0;
- for (j = 0; j < interiorRingPointsNum-1; j++) { //check each
point in the interior ring
- double xCurrent=0.0, yCurrent=0.0, xNext=0.0, yNext=0.0;
- if(GEOSCoordSeq_getX(interiorRingCoordSeq, j,
&xCurrent) == -1 || GEOSCoordSeq_getX(interiorRingCoordSeq, j+1, &xNext) == -1)
- return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getX failed");
- if(GEOSCoordSeq_getY(interiorRingCoordSeq, j,
&yCurrent) == -1 || GEOSCoordSeq_getY(interiorRingCoordSeq, j+1, &yNext) == -1)
- return createException(MAL, "batgeom.Contains",
"GEOSCoordSeq_getY failed");
- //the edge goes from small y to big y (upward
direction) and the point is somewhere in there
- if (yCurrent <= py[i] && yNext >= py[i]) {
- wn+=isLeft( xCurrent, yCurrent, xNext, yNext,
px[i], py[i]);
+ for (j = 1; j < pointsNum[h+1]; j++) { //check each point in
the interior ring
+ if(!checked) {
+ //it is the first iteration, I need to read the
points
+ if(GEOSCoordSeq_getX(interiorRingCoordSeq, j,
&xPoints[h+1][j]) == -1)
+ return createException(MAL,
"batgeom.Contains", "GEOSCoordSeq_getX failed");
+ if(GEOSCoordSeq_getY(interiorRingCoordSeq, j,
&yPoints[h+1][j]) == -1)
+ return createException(MAL,
"batgeom.Contains", "GEOSCoordSeq_getY failed");
+ }
+ //the edge goes from small y to big y (upward
direction) and the point is somewhere in there
+ if (yPoints[h+1][j-1] <= py[i] && yPoints[h+1][j] >=
py[i]) {
+ wn+=isLeft(xPoints[h+1][j-1],
yPoints[h+1][j-1], xPoints[h+1][j], yPoints[h+1][j], px[i], py[i]);
}
//the edge goes from big y to small y (downward
direction) and the point is somewhere in there
- else if (yCurrent >= py[i] && yNext <= py[i]) {
- wn+=isRight( xCurrent, yCurrent, xNext, yNext,
px[i], py[i]);
+ else if (yPoints[h+1][j-1] >= py[i] && yPoints[h+1][j]
<= py[i]) {
+ wn+=isRight(xPoints[h+1][j-1],
yPoints[h+1][j-1], xPoints[h+1][j], yPoints[h+1][j], px[i], py[i]);
}
}
@@ -462,6 +519,8 @@ static str pnpolyWithHoles_(int *out, GE
break;
}
}
+
+ checked = true;
}
BATsetcount(bo,BATcount(bpx));
@@ -469,6 +528,22 @@ static str pnpolyWithHoles_(int *out, GE
BBPreleaseref(bpx->batCacheid);
BBPreleaseref(bpy->batCacheid);
BBPkeepref(*out = bo->batCacheid);
+
+ GDKfree(xPoints[0]);
+ GDKfree(yPoints[0]);
+
+ //if there are no points this is not allocated and thus does not need
to be freed
+ if(checked) {
+ for(i = 1; i< interiorRingsNum+1; i++) {
+ GDKfree(xPoints[i]);
+ GDKfree(yPoints[i]);
+ }
+ }
+
+ GDKfree(xPoints);
+ GDKfree(yPoints);
+ GDKfree(pointsNum);
+
return MAL_SUCCEED;
}
_______________________________________________
checkin-list mailing list
[email protected]
https://www.monetdb.org/mailman/listinfo/checkin-list