Re: Spatial join issues
Paul Ramsey via postgis-users <[email protected]> Wed, 29 Oct 2025 11:52:52 -0700
| Newsgroups | gmane.comp.gis.postgis |
|---|---|
| Message-ID | <CACowWR132zdWkhZNoME16VofBJ=pDbrr1v1hnM8D93UtNxW6WA@mail.gmail.com> |
--00000000000061fc4f064250a472 Content-Type: text/plain; charset="UTF-8" Content-Transfer-Encoding: quoted-printable It's quite robust, but of course it has a big magic number in the middle, in that you have to choose what area ratio value you will use as your "they are the same, as far as I am concerned" threshold. P. On Wed, Oct 29, 2025 at 10:58=E2=80=AFAM Shaozhong SHI <shishaozhong@gmail.= com> wrote: > Basically, is these a robust way to test approximate equal when spatial > fidelity exists. > > On Wed, 29 Oct 2025, 17:39 Shaozhong SHI, <[email protected]> wrote: > >> Very clever. I will spend time to understand. Approximate equal sounds >> very interesting. I just did intersection. Expecting a polygon but go= t >> partial lines. Suspect only part of lines of two seemingly the same >> polygons intersect. Thanks. >> >> On Wed, 29 Oct 2025, 17:25 Paul Ramsey, <[email protected]> >> wrote: >> >>> I would imagine that my little example is a bit error prone, but the >>> classic solution to "do these differently structured things with slight= ly >>> different coordinates cover basically the same space" is to inspect the >>> ratio "area(intersection(a, b)) / area(union(a, b))". You can see >>> intuitively how the more the same two polygons are, the closer that rat= io >>> will be to 1.0. >>> >>> WITH p AS ( >>> SELECT 'POLYGON((0 0, 10 0, 10 10, 0 10, 0 0))'::geometry AS p1, >>> 'POLYGON((0 0, 10 0, 10 5, 10 10, 0 10, 0 0))'::geometry AS p2= , >>> 0.001 AS shift >>> ), >>> shifted AS ( >>> SELECT p1, ST_Translate(p2, shift, shift) AS p2, shift >>> FROM p >>> ), >>> areas AS ( >>> SELECT ST_Area(ST_Union(p1, p2, shift/10)) AS area_union, >>> ST_Area(ST_Intersection(p1, p2, shift/10)) AS area_inter >>> FROM shifted >>> ) >>> SELECT ST_AsText(shifted.p1) AS p1_orig, >>> ST_AsText(shifted.p2) AS p2_orig, >>> area_inter/area_union AS ratio >>> FROM areas, shifted; >>> >>> On Wed, Oct 29, 2025 at 10:17=E2=80=AFAM Paul Ramsey <pramsey@cleverele= phant.ca> >>> wrote: >>> >>>> Greg isn't passing judgement on how hard your problem is, he's saying >>>> you haven't explained it particularly well. Pictures help. Taking a gu= ess >>>> at what you mean, here's some SQL that creates two polygons, with slig= htly >>>> different structure, and slightly different coordinates, that describe= the >>>> same general space in the universe, and then massages them until they = pass >>>> an equals test. >>>> >>>> WITH p AS ( >>>> SELECT 'POLYGON((0 0, 10 0, 10 10, 0 10, 0 0))'::geometry AS p1, >>>> 'POLYGON((0 0, 10 0, 10 5, 10 10, 0 10, 0 0))'::geometry AS p= 2 >>>> ), >>>> shifted AS ( >>>> SELECT p1, ST_Translate(p2, 0.0001, 0.0001) AS p2 >>>> FROM p >>>> ), >>>> rp AS ( >>>> SELECT ST_ReducePrecision(p1,0.1) AS p1, >>>> ST_ReducePrecision(p2,0.1) AS p2 >>>> FROM shifted >>>> ), >>>> snap AS ( >>>> SELECT ST_Snap(p1,p2,0.1) AS p1, >>>> ST_Snap(p2,p1,0.1) AS p2 >>>> FROM rp >>>> ) >>>> SELECT ST_AsText(shifted.p1) AS p1_orig, >>>> ST_AsText(shifted.p2) AS p2_orig, >>>> ST_AsText(snap.p1) AS p1_snap, >>>> ST_AsText(snap.p2) AS p2_snap, >>>> ST_Equals(snap.p1, snap.p2) >>>> FROM snap, shifted; >>>> >>>> -[ RECORD 1 >>>> ]---------------------------------------------------------------------= ----------------------------- >>>> p1_orig | POLYGON((0 0,10 0,10 10,0 10,0 0)) >>>> p2_orig | POLYGON((0.0001 0.0001,10.0001 0.0001,10.0001 >>>> 5.0001,10.0001 10.0001,0.0001 10.0001,0.0001 0.0001)) >>>> p1_snap | POLYGON((0 10,10 10,10 5,10 0,0 0,0 10)) >>>> p2_snap | POLYGON((0 10,10 10,10 5,10 0,0 0,0 10)) >>>> st_equals | t >>>> >>>> On Wed, Oct 29, 2025 at 10:00=E2=80=AFAM Shaozhong SHI <shishaozhong@g= mail.com> >>>> wrote: >>>> >>>>> This is very challenging. Take my words for it. Try on any polygons >>>>> you created and modified. >>>>> >>>>> On Wed, 29 Oct 2025 at 14:36, Greg Troxel <[email protected]> wrote: >>>>> >>>>>> Shaozhong SHI <[email protected]> writes: >>>>>> >>>>>> > Visually, there appears some matching polygons. Even if two >>>>>> geometries >>>>>> > represent the same shape visually, they might not be considered >>>>>> equal due >>>>>> > to tiny differences in precision or metadata. Have you encountere= d >>>>>> > problems of failure of spatial join? How did you overcome the >>>>>> problems? >>>>>> > Regards, David >>>>>> >>>>>> Could you post your example polygons, and the queries you are using? >>>>>> Your question is much too open ended. It even sounds like it might >>>>>> be a >>>>>> request for help with GIS homework, but it's hard to tell :-) >>>>>> >>>>>> --00000000000061fc4f064250a472 Content-Type: text/html; charset="UTF-8" Content-Transfer-Encoding: quoted-printable <div dir=3D"ltr">It's quite robust, but of course it has a big magic nu= mber in the middle, in that you have to choose what area ratio value you wi= ll use as your "they are the same, as far as I am concerned" thre= shold.<div><br></div><div>P.</div></div><br><div class=3D"gmail_quote gmail= _quote_container"><div dir=3D"ltr" class=3D"gmail_attr">On Wed, Oct 29, 202= 5 at 10:58=E2=80=AFAM Shaozhong SHI <<a href=3D"mailto:shishaozhong@gmai= l.com">[email protected]</a>> wrote:<br></div><blockquote class=3D"= gmail_quote" style=3D"margin:0px 0px 0px 0.8ex;border-left-width:1px;border= -left-style:solid;border-left-color:rgb(204,204,204);padding-left:1ex"><div= dir=3D"auto">Basically, is these a robust way to test approximate equal wh= en spatial fidelity exists.</div><br><div class=3D"gmail_quote"><div dir=3D= "ltr" class=3D"gmail_attr">On Wed, 29 Oct 2025, 17:39 Shaozhong SHI, <<a= href=3D"mailto:[email protected]" target=3D"_blank">shishaozhong@gmai= l.com</a>> wrote:<br></div><blockquote class=3D"gmail_quote" style=3D"ma= rgin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border= -left-color:rgb(204,204,204);padding-left:1ex"><div dir=3D"auto">Very cleve= r.=C2=A0 I will spend time to understand.=C2=A0 Approximate equal sounds ve= ry interesting.=C2=A0 =C2=A0I just did intersection.=C2=A0 Expecting a poly= gon but got partial lines.=C2=A0 Suspect only part of lines of two seemingl= y the same polygons intersect. Thanks.=C2=A0=C2=A0</div><br><div class=3D"g= mail_quote"><div dir=3D"ltr" class=3D"gmail_attr">On Wed, 29 Oct 2025, 17:2= 5 Paul Ramsey, <<a href=3D"mailto:[email protected]" rel=3D"nore= ferrer" target=3D"_blank">[email protected]</a>> wrote:<br></div= ><blockquote class=3D"gmail_quote" style=3D"margin:0px 0px 0px 0.8ex;border= -left-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);= padding-left:1ex"><div dir=3D"ltr">I would imagine that my little example i= s a bit error prone, but the classic solution to "do these differently= structured things with slightly different coordinates cover basically the = same space" is to inspect the ratio =C2=A0"area(intersection(a, b= )) / area(union(a, b))". You can see intuitively how the more the same= two polygons are, the closer that ratio will be to 1.0.=C2=A0<div><br></di= v><div>WITH p AS (<br>=C2=A0 SELECT 'POLYGON((0 0, 10 0, 10 10, 0 10, 0= 0))'::geometry AS p1,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0'POLYGO= N((0 0, 10 0, 10 5, 10 10, 0 10, 0 0))'::geometry AS p2,<br>=C2=A0 =C2= =A0 =C2=A0 =C2=A0 =C2=A00.001 AS shift<br>),<br>shifted AS (<br>=C2=A0 SELE= CT p1, ST_Translate(p2, shift, shift) AS p2, shift<br>=C2=A0 FROM p<br>),<b= r>areas AS (<br>=C2=A0 SELECT ST_Area(ST_Union(p1, p2, shift/10)) AS area_u= nion,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0ST_Area(ST_Intersection(p1, p2, = shift/10)) AS area_inter<br>=C2=A0 FROM shifted<br>)<br>SELECT ST_AsText(sh= ifted.p1) AS p1_orig, <br>=C2=A0 =C2=A0 =C2=A0 =C2=A0ST_AsText(shifted.p2) = AS p2_orig, <br>=C2=A0 =C2=A0 =C2=A0 =C2=A0area_inter/area_union AS ratio<b= r>FROM areas, shifted;<br></div></div><br><div class=3D"gmail_quote"><div d= ir=3D"ltr" class=3D"gmail_attr">On Wed, Oct 29, 2025 at 10:17=E2=80=AFAM Pa= ul Ramsey <<a href=3D"mailto:[email protected]" rel=3D"noreferre= r noreferrer" target=3D"_blank">[email protected]</a>> wrote:<br= ></div><blockquote class=3D"gmail_quote" style=3D"margin:0px 0px 0px 0.8ex;= border-left-width:1px;border-left-style:solid;border-left-color:rgb(204,204= ,204);padding-left:1ex"><div dir=3D"ltr">Greg isn't passing judgement o= n how hard your problem is, he's saying you haven't explained it pa= rticularly well. Pictures help. Taking a guess at what you mean, here's= some SQL that creates two polygons, with slightly different structure, and= slightly different coordinates, that describe the same general space in th= e universe, and then massages them until they pass an equals test.<div><br>= </div><div>WITH p AS (<br>=C2=A0 SELECT 'POLYGON((0 0, 10 0, 10 10, 0 1= 0, 0 0))'::geometry AS p1,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0'PO= LYGON((0 0, 10 0, 10 5, 10 10, 0 10, 0 0))'::geometry AS p2<br>),<br>sh= ifted AS (<br>=C2=A0 SELECT p1, ST_Translate(p2, 0.0001, 0.0001) AS p2<br>= =C2=A0 FROM p<br>),<br>rp AS (<br>=C2=A0 SELECT ST_ReducePrecision(p1,0.1) = AS p1,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0ST_ReducePrecision(p2,0.1) AS p= 2<br>=C2=A0 FROM shifted<br>),<br>snap AS (<br>=C2=A0 SELECT ST_Snap(p1,p2,= 0.1) AS p1, <br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0ST_Snap(p2,p1,0.1) AS p2<= br>=C2=A0 FROM rp<br>)<br>SELECT ST_AsText(shifted.p1) AS p1_orig, <br>=C2= =A0 =C2=A0 =C2=A0 =C2=A0ST_AsText(shifted.p2) AS p2_orig, <br>=C2=A0 =C2=A0= =C2=A0 =C2=A0ST_AsText(snap.p1) AS p1_snap, <br>=C2=A0 =C2=A0 =C2=A0 =C2= =A0ST_AsText(snap.p2) AS p2_snap,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0ST_Equals(s= nap.p1, snap.p2)<br>FROM snap, shifted;<br></div><div><br></div><div>-[ REC= ORD 1 ]--------------------------------------------------------------------= ------------------------------<br>p1_orig =C2=A0 | POLYGON((0 0,10 0,10 10,= 0 10,0 0))<br>p2_orig =C2=A0 | POLYGON((0.0001 0.0001,10.0001 0.0001,10.000= 1 5.0001,10.0001 10.0001,0.0001 10.0001,0.0001 0.0001))<br>p1_snap =C2=A0 |= POLYGON((0 10,10 10,10 5,10 0,0 0,0 10))<br>p2_snap =C2=A0 | POLYGON((0 10= ,10 10,10 5,10 0,0 0,0 10))<br>st_equals | t<br></div></div><br><div class= =3D"gmail_quote"><div dir=3D"ltr" class=3D"gmail_attr">On Wed, Oct 29, 2025= at 10:00=E2=80=AFAM Shaozhong SHI <<a href=3D"mailto:shishaozhong@gmail= .com" rel=3D"noreferrer noreferrer" target=3D"_blank">[email protected]= m</a>> wrote:<br></div><blockquote class=3D"gmail_quote" style=3D"margin= :0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;border-lef= t-color:rgb(204,204,204);padding-left:1ex"><div dir=3D"ltr">This is very ch= allenging.=C2=A0 Take my words for it.=C2=A0 Try on any polygons you create= d and modified.</div><br><div class=3D"gmail_quote"><div dir=3D"ltr" class= =3D"gmail_attr">On Wed, 29 Oct 2025 at 14:36, Greg Troxel <<a href=3D"ma= ilto:[email protected]" rel=3D"noreferrer noreferrer" target=3D"_blank">gdt@le= xort.com</a>> wrote:<br></div><blockquote class=3D"gmail_quote" style=3D= "margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-style:solid;bor= der-left-color:rgb(204,204,204);padding-left:1ex">Shaozhong SHI <<a href= =3D"mailto:[email protected]" rel=3D"noreferrer noreferrer" target=3D"= _blank">[email protected]</a>> writes:<br> <br> >=C2=A0 =C2=A0Visually, there appears some matching polygons.=C2=A0 Even= if two geometries<br> > represent the same shape visually, they might not be considered equal = due<br> > to tiny differences in precision or metadata.=C2=A0 Have you encounter= ed<br> > problems of failure of spatial join?=C2=A0 How did you overcome the pr= oblems?<br> > Regards, David<br> <br> Could you post your example polygons, and the queries you are using?<br> Your question is much too open ended.=C2=A0 It even sounds like it might be= a<br> request for help with GIS homework, but it's hard to tell :-)<br> <br> </blockquote></div> </blockquote></div> </blockquote></div> </blockquote></div> </blockquote></div> </blockquote></div> --00000000000061fc4f064250a472--