Re: Spatial join issues

Shaozhong SHI <[email protected]> Wed, 29 Oct 2025 17:58:07 +0000
Newsgroups gmane.comp.gis.postgis
Message-ID <CA+i5JwY6qQBw8bYDreuJveCed2APX+6cS2f9zm70CKLEGO-RvA@mail.gmail.com>
--000000000000a0a51406424fe04f
Content-Type: text/plain; charset="UTF-8"
Content-Transfer-Encoding: quoted-printable

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 got
> 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 slightl=
y
>> 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 rati=
o
>> 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@cleverelep=
hant.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 gue=
ss
>>> at what you mean, here's some SQL that creates two polygons, with sligh=
tly
>>> different structure, and slightly different coordinates, that describe =
the
>>> same general space in the universe, and then massages them until they p=
ass
>>> 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 p2
>>> ),
>>> 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.000=
1
>>> 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@gm=
ail.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 encountered
>>>>> > 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 b=
e
>>>>> a
>>>>> request for help with GIS homework, but it's hard to tell :-)
>>>>>
>>>>>

--000000000000a0a51406424fe04f
Content-Type: text/html; charset="UTF-8"
Content-Transfer-Encoding: quoted-printable

<div dir=3D"auto">Basically, is these a robust way to test approximate equa=
l when spatial fidelity exists.</div><br><div class=3D"gmail_quote gmail_qu=
ote_container"><div dir=3D"ltr" class=3D"gmail_attr">On Wed, 29 Oct 2025, 1=
7:39 Shaozhong SHI, &lt;<a href=3D"mailto:[email protected]">shishaozh=
[email protected]</a>&gt; wrote:<br></div><blockquote class=3D"gmail_quote" sty=
le=3D"margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div d=
ir=3D"auto">Very clever.=C2=A0 I will spend time to understand.=C2=A0 Appro=
ximate equal sounds very interesting.=C2=A0 =C2=A0I just did intersection.=
=C2=A0 Expecting a polygon but got partial lines.=C2=A0 Suspect only part o=
f lines of two seemingly the same polygons intersect. Thanks.=C2=A0=C2=A0</=
div><br><div class=3D"gmail_quote"><div dir=3D"ltr" class=3D"gmail_attr">On=
 Wed, 29 Oct 2025, 17:25 Paul Ramsey, &lt;<a href=3D"mailto:pramsey@clevere=
lephant.ca" target=3D"_blank" rel=3D"noreferrer">[email protected]<=
/a>&gt; wrote:<br></div><blockquote class=3D"gmail_quote" style=3D"margin:0=
 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><div dir=3D"ltr">I w=
ould imagine that my little example is a bit error prone, but the classic s=
olution to &quot;do these differently structured things with slightly diffe=
rent coordinates cover basically the same space&quot; is to inspect the rat=
io =C2=A0&quot;area(intersection(a, b)) / area(union(a, b))&quot;. You can =
see intuitively how the more the same two polygons are, the closer that rat=
io will be to 1.0.=C2=A0<div><br></div><div>WITH p AS (<br>=C2=A0 SELECT &#=
39;POLYGON((0 0, 10 0, 10 10, 0 10, 0 0))&#39;::geometry AS p1,<br>=C2=A0 =
=C2=A0 =C2=A0 =C2=A0 =C2=A0&#39;POLYGON((0 0, 10 0, 10 5, 10 10, 0 10, 0 0)=
)&#39;::geometry AS p2,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A00.001 AS shift=
<br>),<br>shifted AS (<br>=C2=A0 SELECT p1, ST_Translate(p2, shift, shift) =
AS p2, shift<br>=C2=A0 FROM p<br>),<br>areas AS (<br>=C2=A0 SELECT ST_Area(=
ST_Union(p1, p2, shift/10)) AS area_union,<br>=C2=A0 =C2=A0 =C2=A0 =C2=A0 =
=C2=A0ST_Area(ST_Intersection(p1, p2, shift/10)) AS area_inter<br>=C2=A0 FR=
OM shifted<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=A0area_inter/area_union AS ratio<br>FROM areas, shifted;<br></div></di=
v><br><div class=3D"gmail_quote"><div dir=3D"ltr" class=3D"gmail_attr">On W=
ed, Oct 29, 2025 at 10:17=E2=80=AFAM Paul Ramsey &lt;<a href=3D"mailto:pram=
[email protected]" rel=3D"noreferrer noreferrer" target=3D"_blank">pram=
[email protected]</a>&gt; wrote:<br></div><blockquote class=3D"gmail_qu=
ote" style=3D"margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-st=
yle:solid;border-left-color:rgb(204,204,204);padding-left:1ex"><div dir=3D"=
ltr">Greg isn&#39;t passing judgement on how hard your problem is, he&#39;s=
 saying you haven&#39;t explained it particularly well. Pictures help. Taki=
ng a guess at what you mean, here&#39;s some SQL that creates two polygons,=
 with slightly different structure, and slightly different coordinates, tha=
t describe the same general space in the universe, and then massages them u=
ntil they pass an equals test.<div><br></div><div>WITH p AS (<br>=C2=A0 SEL=
ECT &#39;POLYGON((0 0, 10 0, 10 10, 0 10, 0 0))&#39;::geometry AS p1,<br>=
=C2=A0 =C2=A0 =C2=A0 =C2=A0 =C2=A0&#39;POLYGON((0 0, 10 0, 10 5, 10 10, 0 1=
0, 0 0))&#39;::geometry AS p2<br>),<br>shifted 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 p2<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) A=
S 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(snap.p1, snap.p2)<br>FROM snap, shift=
ed;<br></div><div><br></div><div>-[ RECORD 1 ]-----------------------------=
---------------------------------------------------------------------<br>p1=
_orig =C2=A0 | POLYGON((0 0,10 0,10 10,0 10,0 0))<br>p2_orig =C2=A0 | POLYG=
ON((0.0001 0.0001,10.0001 0.0001,10.0001 5.0001,10.0001 10.0001,0.0001 10.0=
001,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_eq=
uals | t<br></div></div><br><div class=3D"gmail_quote"><div dir=3D"ltr" cla=
ss=3D"gmail_attr">On Wed, Oct 29, 2025 at 10:00=E2=80=AFAM Shaozhong SHI &l=
t;<a href=3D"mailto:[email protected]" rel=3D"noreferrer noreferrer" t=
arget=3D"_blank">[email protected]</a>&gt; 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">This is very challenging.=C2=A0 Take my words for it=
.=C2=A0 Try on any polygons you created 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 &lt;<a href=3D"mailto:[email protected]" rel=3D"noreferr=
er noreferrer" target=3D"_blank">[email protected]</a>&gt; wrote:<br></div><bl=
ockquote class=3D"gmail_quote" style=3D"margin:0px 0px 0px 0.8ex;border-lef=
t-width:1px;border-left-style:solid;border-left-color:rgb(204,204,204);padd=
ing-left:1ex">Shaozhong SHI &lt;<a href=3D"mailto:[email protected]" r=
el=3D"noreferrer noreferrer" target=3D"_blank">[email protected]</a>&g=
t; writes:<br>
<br>
&gt;=C2=A0 =C2=A0Visually, there appears some matching polygons.=C2=A0 Even=
 if two geometries<br>
&gt; represent the same shape visually, they might not be considered equal =
due<br>
&gt; to tiny differences in precision or metadata.=C2=A0 Have you encounter=
ed<br>
&gt; problems of failure of spatial join?=C2=A0 How did you overcome the pr=
oblems?<br>
&gt; 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&#39;s hard to tell :-)<br>
<br>
</blockquote></div>
</blockquote></div>
</blockquote></div>
</blockquote></div>
</blockquote></div>

--000000000000a0a51406424fe04f--