Movatterモバイル変換


[0]ホーム

URL:


Skip to content

Navigation Menu

Sign in
Appearance settings

Search code, repositories, users, issues, pull requests...

Provide feedback

We read every piece of feedback, and take your input very seriously.

Saved searches

Use saved searches to filter your results more quickly

Sign up
Appearance settings

Commit49475aa

Browse files
committed
Correct calculations of overlap and contains operations over polygons.
1 parent1f4b046 commit49475aa

File tree

5 files changed

+346
-49
lines changed

5 files changed

+346
-49
lines changed

‎src/backend/utils/adt/geo_ops.c

Lines changed: 203 additions & 49 deletions
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@
88
*
99
*
1010
* IDENTIFICATION
11-
* $PostgreSQL: pgsql/src/backend/utils/adt/geo_ops.c,v 1.102 2009/06/23 16:25:02 tgl Exp $
11+
* $PostgreSQL: pgsql/src/backend/utils/adt/geo_ops.c,v 1.103 2009/07/28 09:47:59 teodor Exp $
1212
*
1313
*-------------------------------------------------------------------------
1414
*/
@@ -66,6 +66,8 @@ static bool has_interpt_sl(LSEG *lseg, LINE *line);
6666
staticdoubledist_pl_internal(Point*pt,LINE*line);
6767
staticdoubledist_ps_internal(Point*pt,LSEG*lseg);
6868
staticPoint*line_interpt_internal(LINE*l1,LINE*l2);
69+
staticboollseg_inside_poly(Point*a,Point*b,POLYGON*poly,intstart);
70+
staticPoint*lseg_interpt_internal(LSEG*l1,LSEG*l2);
6971

7072

7173
/*
@@ -2352,15 +2354,9 @@ lseg_center(PG_FUNCTION_ARGS)
23522354
PG_RETURN_POINT_P(result);
23532355
}
23542356

2355-
2356-
/* lseg_interpt -
2357-
*Find the intersection point of two segments (if any).
2358-
*/
2359-
Datum
2360-
lseg_interpt(PG_FUNCTION_ARGS)
2357+
staticPoint*
2358+
lseg_interpt_internal(LSEG*l1,LSEG*l2)
23612359
{
2362-
LSEG*l1=PG_GETARG_LSEG_P(0);
2363-
LSEG*l2=PG_GETARG_LSEG_P(1);
23642360
Point*result;
23652361
LINEtmp1,
23662362
tmp2;
@@ -2372,15 +2368,18 @@ lseg_interpt(PG_FUNCTION_ARGS)
23722368
line_construct_pts(&tmp2,&l2->p[0],&l2->p[1]);
23732369
result=line_interpt_internal(&tmp1,&tmp2);
23742370
if (!PointerIsValid(result))
2375-
PG_RETURN_NULL();
2371+
returnNULL;
23762372

23772373
/*
23782374
* If the line intersection point isn't within l1 (or equivalently l2),
23792375
* there is no valid segment intersection point at all.
23802376
*/
23812377
if (!on_ps_internal(result,l1)||
23822378
!on_ps_internal(result,l2))
2383-
PG_RETURN_NULL();
2379+
{
2380+
pfree(result);
2381+
returnNULL;
2382+
}
23842383

23852384
/*
23862385
* If there is an intersection, then check explicitly for matching
@@ -2400,6 +2399,23 @@ lseg_interpt(PG_FUNCTION_ARGS)
24002399
result->y=l1->p[1].y;
24012400
}
24022401

2402+
returnresult;
2403+
}
2404+
2405+
/* lseg_interpt -
2406+
*Find the intersection point of two segments (if any).
2407+
*/
2408+
Datum
2409+
lseg_interpt(PG_FUNCTION_ARGS)
2410+
{
2411+
LSEG*l1=PG_GETARG_LSEG_P(0);
2412+
LSEG*l2=PG_GETARG_LSEG_P(1);
2413+
Point*result;
2414+
2415+
result=lseg_interpt_internal(l1,l2);
2416+
if (!PointerIsValid(result))
2417+
PG_RETURN_NULL();
2418+
24032419
PG_RETURN_POINT_P(result);
24042420
}
24052421

@@ -3742,10 +3758,7 @@ poly_same(PG_FUNCTION_ARGS)
37423758
}
37433759

37443760
/*-----------------------------------------------------------------
3745-
* Determine if polygon A overlaps polygon B by determining if
3746-
* their bounding boxes overlap.
3747-
*
3748-
* XXX ought to do a more correct check!
3761+
* Determine if polygon A overlaps polygon B
37493762
*-----------------------------------------------------------------*/
37503763
Datum
37513764
poly_overlap(PG_FUNCTION_ARGS)
@@ -3754,7 +3767,54 @@ poly_overlap(PG_FUNCTION_ARGS)
37543767
POLYGON*polyb=PG_GETARG_POLYGON_P(1);
37553768
boolresult;
37563769

3757-
result=box_ov(&polya->boundbox,&polyb->boundbox);
3770+
/* Quick check by bounding box */
3771+
result= (polya->npts>0&&polyb->npts>0&&
3772+
box_ov(&polya->boundbox,&polyb->boundbox)) ? true : false;
3773+
3774+
/*
3775+
* Brute-force algorithm - try to find intersected edges,
3776+
* if so then polygons are overlapped else check is one
3777+
* polygon inside other or not by testing single point
3778+
* of them.
3779+
*/
3780+
if (result)
3781+
{
3782+
intia,ib;
3783+
LSEGsa,sb;
3784+
3785+
/* Init first of polya's edge with last point */
3786+
sa.p[0]=polya->p[polya->npts-1];
3787+
result= false;
3788+
3789+
for(ia=0;ia<polya->npts&&result== false;ia++)
3790+
{
3791+
/* Second point of polya's edge is a current one */
3792+
sa.p[1]=polya->p[ia];
3793+
3794+
/* Init first of polyb's edge with last point */
3795+
sb.p[0]=polyb->p[polyb->npts-1];
3796+
3797+
for(ib=0;ib<polyb->npts&&result== false;ib++)
3798+
{
3799+
sb.p[1]=polyb->p[ib];
3800+
result=lseg_intersect_internal(&sa,&sb);
3801+
sb.p[0]=sb.p[1];
3802+
}
3803+
3804+
/*
3805+
* move current endpoint to the first point
3806+
* of next edge
3807+
*/
3808+
sa.p[0]=sa.p[1];
3809+
}
3810+
3811+
if (result==false)
3812+
{
3813+
result= (point_inside(polya->p,polyb->npts,polyb->p)
3814+
||
3815+
point_inside(polyb->p,polya->npts,polya->p) );
3816+
}
3817+
}
37583818

37593819
/*
37603820
* Avoid leaking memory for toasted inputs ... needed for rtree indexes
@@ -3765,6 +3825,119 @@ poly_overlap(PG_FUNCTION_ARGS)
37653825
PG_RETURN_BOOL(result);
37663826
}
37673827

3828+
/*
3829+
* Tests special kind of segment for in/out of polygon.
3830+
* Special kind means:
3831+
* - point a should be on segment s
3832+
* - segment (a,b) should not be contained by s
3833+
* Returns true if:
3834+
* - segment (a,b) is collinear to s and (a,b) is in polygon
3835+
* - segment (a,b) s not collinear to s. Note: that doesn't
3836+
* mean that segment is in polygon!
3837+
*/
3838+
3839+
staticbool
3840+
touched_lseg_inside_poly(Point*a,Point*b,LSEG*s,POLYGON*poly,intstart)
3841+
{
3842+
/* point a is on s, b is not */
3843+
LSEGt;
3844+
3845+
t.p[0]=*a;
3846+
t.p[1]=*b;
3847+
3848+
#definePOINTEQ(pt1,pt2) (FPeq((pt1)->x, (pt2)->x) && FPeq((pt1)->y, (pt2)->y))
3849+
if (POINTEQ(a,s->p) )
3850+
{
3851+
if (on_ps_internal(s->p+1,&t) )
3852+
returnlseg_inside_poly(b,s->p+1,poly,start);
3853+
}
3854+
elseif (POINTEQ(a,s->p+1))
3855+
{
3856+
if (on_ps_internal(s->p,&t) )
3857+
returnlseg_inside_poly(b,s->p,poly,start);
3858+
}
3859+
elseif (on_ps_internal(s->p,&t) )
3860+
{
3861+
returnlseg_inside_poly(b,s->p,poly,start);
3862+
}
3863+
elseif (on_ps_internal(s->p+1,&t) )
3864+
{
3865+
returnlseg_inside_poly(b,s->p+1,poly,start);
3866+
}
3867+
3868+
return true;/* may be not true, but that will check later */
3869+
}
3870+
3871+
/*
3872+
* Returns true if segment (a,b) is in polygon, option
3873+
* start is used for optimization - function checks
3874+
* polygon's edges started from start
3875+
*/
3876+
staticbool
3877+
lseg_inside_poly(Point*a,Point*b,POLYGON*poly,intstart)
3878+
{
3879+
LSEGs,
3880+
t;
3881+
inti;
3882+
boolres= true,
3883+
intersection= false;
3884+
3885+
t.p[0]=*a;
3886+
t.p[1]=*b;
3887+
s.p[0]=poly->p[(start==0) ? (poly->npts-1) : (start-1)];
3888+
3889+
for(i=start;i<poly->npts&&res== true;i++)
3890+
{
3891+
Point*interpt;
3892+
3893+
s.p[1]=poly->p[i];
3894+
3895+
if (on_ps_internal(t.p,&s) )
3896+
{
3897+
if (on_ps_internal(t.p+1,&s) )
3898+
return true;/* t is contained by s */
3899+
3900+
/* Y-cross */
3901+
res=touched_lseg_inside_poly(t.p,t.p+1,&s,poly,i+1);
3902+
}
3903+
elseif (on_ps_internal(t.p+1,&s) )
3904+
{
3905+
/* Y-cross */
3906+
res=touched_lseg_inside_poly(t.p+1,t.p,&s,poly,i+1);
3907+
}
3908+
elseif ( (interpt=lseg_interpt_internal(&t,&s))!=NULL )
3909+
{
3910+
/*
3911+
* segments are X-crossing, go to check each subsegment
3912+
*/
3913+
3914+
intersection= true;
3915+
res=lseg_inside_poly(t.p,interpt,poly,i+1);
3916+
if (res)
3917+
res=lseg_inside_poly(t.p+1,interpt,poly,i+1);
3918+
pfree(interpt);
3919+
}
3920+
3921+
s.p[0]=s.p[1];
3922+
}
3923+
3924+
if (res&& !intersection)
3925+
{
3926+
Pointp;
3927+
3928+
/*
3929+
* if X-intersection wasn't found then check central point
3930+
* of tested segment. In opposite case we already check all
3931+
* subsegments
3932+
*/
3933+
p.x= (t.p[0].x+t.p[1].x) /2.0;
3934+
p.y= (t.p[0].y+t.p[1].y) /2.0;
3935+
3936+
res=point_inside(&p,poly->npts,poly->p);
3937+
}
3938+
3939+
returnres;
3940+
}
37683941

37693942
/*-----------------------------------------------------------------
37703943
* Determine if polygon A contains polygon B.
@@ -3775,49 +3948,30 @@ poly_contain(PG_FUNCTION_ARGS)
37753948
POLYGON*polya=PG_GETARG_POLYGON_P(0);
37763949
POLYGON*polyb=PG_GETARG_POLYGON_P(1);
37773950
boolresult;
3778-
inti;
37793951

37803952
/*
37813953
* Quick check to see if bounding box is contained.
37823954
*/
3783-
if (DatumGetBool(DirectFunctionCall2(box_contain,
3784-
BoxPGetDatum(&polya->boundbox),
3785-
BoxPGetDatum(&polyb->boundbox))))
3955+
if (polya->npts>0&&polyb->npts>0&&
3956+
DatumGetBool(DirectFunctionCall2(box_contain,
3957+
BoxPGetDatum(&polya->boundbox),
3958+
BoxPGetDatum(&polyb->boundbox))))
37863959
{
3787-
result= true;/* assume true for now */
3788-
for (i=0;i<polyb->npts;i++)
3789-
{
3790-
if (point_inside(&(polyb->p[i]),polya->npts,&(polya->p[0]))==0)
3791-
{
3792-
#ifdefGEODEBUG
3793-
printf("poly_contain- point (%f,%f) not in polygon\n",polyb->p[i].x,polyb->p[i].y);
3794-
#endif
3795-
result= false;
3796-
break;
3797-
}
3798-
}
3799-
if (result)
3960+
inti;
3961+
LSEGs;
3962+
3963+
s.p[0]=polyb->p[polyb->npts-1];
3964+
result= true;
3965+
3966+
for(i=0;i<polyb->npts&&result== true;i++)
38003967
{
3801-
for (i=0;i<polya->npts;i++)
3802-
{
3803-
if (point_inside(&(polya->p[i]),polyb->npts,&(polyb->p[0]))==1)
3804-
{
3805-
#ifdefGEODEBUG
3806-
printf("poly_contain- point (%f,%f) in polygon\n",polya->p[i].x,polya->p[i].y);
3807-
#endif
3808-
result= false;
3809-
break;
3810-
}
3811-
}
3968+
s.p[1]=polyb->p[i];
3969+
result=lseg_inside_poly(s.p,s.p+1,polya,0);
3970+
s.p[0]=s.p[1];
38123971
}
38133972
}
38143973
else
38153974
{
3816-
#ifdefGEODEBUG
3817-
printf("poly_contain- bound box ((%f,%f),(%f,%f)) not inside ((%f,%f),(%f,%f))\n",
3818-
polyb->boundbox.low.x,polyb->boundbox.low.y,polyb->boundbox.high.x,polyb->boundbox.high.y,
3819-
polya->boundbox.low.x,polya->boundbox.low.y,polya->boundbox.high.x,polya->boundbox.high.y);
3820-
#endif
38213975
result= false;
38223976
}
38233977

‎src/test/regress/expected/create_index.out

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -53,6 +53,8 @@ CREATE INDEX gpolygonind ON polygon_tbl USING gist (f1);
5353
CREATE INDEX gcircleind ON circle_tbl USING gist (f1);
5454
CREATE TEMP TABLE gpolygon_tbl AS
5555
SELECT polygon(home_base) AS f1 FROM slow_emp4000;
56+
INSERT INTO gpolygon_tbl VALUES ( '(1000,0,0,1000)' );
57+
INSERT INTO gpolygon_tbl VALUES ( '(0,1000,1000,1000)' );
5658
CREATE TEMP TABLE gcircle_tbl AS
5759
SELECT circle(home_base) AS f1 FROM slow_emp4000;
5860
CREATE INDEX ggpolygonind ON gpolygon_tbl USING gist (f1);

0 commit comments

Comments
 (0)

[8]ページ先頭

©2009-2025 Movatter.jp