[postgis-users] touches() not returning correct values, maybe?

[Jeff D. Hamann]

I've been examining postgis for a little while (for production, teaching,
and research) and so far I'm very pleased with the results. I only wish
Arc/Info had SQL like this!

I've been loading shape files, generating adjacency lists, and then
exporting the results to a text file for a project and noticed the
touches() function isn't returning the correct results.

I've got two polygons (out of 53) that contain common points (from R),

> dp.stands$Shapes[[1]]$verts
          [,1]     [,2]
 [1,] 500141.9 505664.9
 [2,] 500237.9 505734.9
 [3,] 500435.0 505821.6
 [4,] 500435.0 504810.6
 [5,] 500167.6 504636.5
 [6,] 500143.3 504616.9
 [7,] 500143.1 504616.5
 [8,] 500142.8 504616.3
 [9,] 499962.4 504434.0
[10,] 499962.2 504433.7
[11,] 499950.9 504420.8
[12,] 499788.2 504212.4
[13,] 499786.5 504209.5
[14,] 499782.8 504205.1
[15,] 499585.7 503918.4
[16,] 499584.0 503915.2
[17,] 499582.0 503912.7
[18,] 499391.4 503606.6
[19,] 499386.0 503596.3
[20,] 499384.1 503593.6
[21,] 499253.8 503329.0
[22,] 499019.2 502965.0
[23,] 498866.0 502743.1
[24,] 498866.0 502743.1
[25,] 498775.5 502993.1
[26,] 498801.3 503256.2
[27,] 498845.7 503517.7
[28,] 498958.4 503821.9
[29,] 499112.3 504197.1
[30,] 499275.2 504463.3
[31,] 499564.9 504869.7
[32,] 499770.7 505108.3
[33,] 499917.6 505287.0
[34,] 500053.7 505517.0
[35,] 500141.9 505664.9
>

and

> dp.stands$Shapes[[39]]$verts
          [,1]     [,2]
 [1,] 497710.9 503306.0
 [2,] 497740.8 503426.9
 [3,] 497776.1 503659.6
 [4,] 497803.2 503884.1
 [5,] 497852.9 504042.5
 [6,] 497910.0 504186.1
 [7,] 497995.7 504276.2
 [8,] 498167.8 504313.2
 [9,] 498412.8 504385.5
[10,] 498551.0 504444.2
[11,] 498577.9 504481.8
[12,] 498966.3 503841.2
[13,] 498958.4 503821.9
[14,] 498845.7 503517.7
[15,] 498801.3 503256.2
[16,] 498775.5 502993.1
[17,] 498858.9 502762.7
[18,] 498225.1 502836.6
[19,] 497846.5 503108.7
[20,] 497710.9 503306.0
>

Using R's %in% operator,

> dp.stands$Shapes[[1]]$verts %in% dp.stands$Shapes[[39]]$verts
 [1] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
FALSE
[13] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
[25]  TRUE  TRUE  TRUE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
[37] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
[49] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE  TRUE
[61]  TRUE  TRUE  TRUE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
>

Showing that there are points that are common to the two sets.

When I used the touches() function to generate an list of adjacent polygons,

drop table adj_list;
select
  t1.gid as gid,
  t2.gid as adj_gid
into table
	adj_list
 from
  stands t1,
  stands t2
 where
  touches( t1.boundary, t2.boundary ) AND
  t1.boundary && t2.boundary
order by
	t1.gid,
	t2.gid;

The resulting list was,

1,16
1,27
1,38
2,4
2,18
...blah, blah, blah...
2,19
2,27
2,38

In which the first three (for polygon one) are correct but there's no
polygon 39. I've attached the pdf of the shapefile (hope it gets through)
and the accompanying adjacency list and sql file as an example for good
measure.

Is there something I should be doing to my query differently? Is this me?


Jeff.








-- 
Forest Informatics, Inc.
PO Box 1421
Corvallis, Oregon 97339-1421







-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: adj_list.txt
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20070102/d43c71d1/attachment.txt>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: dp_adj_prob.pdf
Type: application/pdf
Size: 43997 bytes
Desc: not available
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20070102/d43c71d1/attachment.pdf>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: stands.sql
URL: <http://lists.osgeo.org/pipermail/postgis-users/attachments/20070102/d43c71d1/attachment.ksh>