No intersection found between polygons I know intersect!

Question

I'm using two geographies: one is made up of small areas (about 2200), the other large areas (about 70). Both cover all of Australia, and the small areas are not necessarily nested inside large areas. For the large areas, I was given a kml file, and a mapinfo file (both come from the Australian Bureau of Statistics). I converted both into .shp format, for use in R (using rgeos).

For each of the small areas, I would like the proportion that fits in each large area. This can be done with

gArea(gIntersection(small area, large area))/gArea(small area)

Which I execute inside a foreach loop that returns a 2200 by 70 matrix of weights.

My problem is that for some of the small regions, this method does not find any intersection (as in it returns zero for all larger areas). I have run it on shape files constructed from both the kml file and the mapinfo files, and areas with no intersection are different!

What would be the most likely cause of this sort of problem?

---

Edit-

Thanks guys for your help. I have put together a little zip file containing my two shape files and the R script I'm using. I'm pretty fresh to this (3 days) so reading up on projections. Thanks again. (The file is about 60 meg)

https://dl.dropbox.com/u/63100926/Example.zip

---

Edit 2:

Thanks to your advice, I now know that my small regions are in GRS80, and my big regions are in WGS84. So thanks! I am having a bit of difficulty finding out how to change them. Any pointers on where to look? Cheers.

Answer 1

It sounds as though projection is the problem.
The Kml is no doubt in the google required wgs projection
While the mapinfo file could probably not be in that projection.
For getting accurate measurements you should re-project to a "local" unitized projection.
such as meters.

Answer 2

Update: after downloading your example files, the main issue seems to be geometry errors - a small subset of LFS_Regions (Tassie only) gives 105 errors in QGIS, for example. This is probably where the inconsistencies come from. Try repairing geometries in GRASS, perhaps! Some notes on this from Faunalia here.

---

You'll need to put your shapefiles into the same projection - and, if you're concerned about area calculations, use a projection which is appropriate; GDA94 / Geoscience Australia Lambert (EPSG:3112) might work for you. Using an example from UTM 36S:

library(rgdal) # for spatial transformation
library(maptools) # for working with shapefiles
library(rgeos) # spTransform method 

Load shapefiles using readOGR (to preserve projection information which readShape* won't do!), and look at the summaries:

f <- readOGR(dsn = getwd(), layer = "gma_east_utm36")
g <- readOGR(dsn=getwd(), layer = "gma_overlap")
summary(f)

Object of class SpatialPolygonsDataFrame
Coordinates:
        min       max
x  189715.7  544014.4
y 8448195.7 8859829.1
Is projected: TRUE 
proj4string :
[+proj=utm +zone=36 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0]
Data attributes ... 

summary(g)
Object of class SpatialPolygonsDataFrame
Coordinates:
        min       max
x  32.45725  32.82836
y -12.58917 -10.92426
Is projected: FALSE 
proj4string :
[+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0]
Data attributes ...

f is projected (UTM / metres) and g is not. Project g:

g.utm36 <- spTransform(hrr.shp, CRS("+init=epsg:32736"))

Now you can create a new object with the overlap:

fg <- gIntersection(f, g.utm36)

Since your 'large areas' are in WGS84 (by which I assume you mean unprojected lat/lon?) you may want to transform both your shapes into the GDA94 projection mentioned above and do your calculations there.