Australian Bureau of Statistics (ABS) provide census information for all Australia. However, ABS has its own division of a state/suburb/town called Statistical Areas Level 1-4 (SA), which is different from the usual government defined suburb-like divisions. ABS also provides ESRI shapefiles for its SA1 divisions.

Problem at hand is that I have number of addresses and I used google maps to find Lat/Longs for them. Now I want to find out in R-Language as in which SA1 division each address belongs to?

I am following the approach mentioned at this Bear-in-Park problem page. Main question that I am struggling with is that I believe Google map's lat/long are in different CRS than ABS ones. I don't know which CRS system Google map uses but some information about ABS shapefile is as follows,

"+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs"

So how to assign a CRS to Google-Map downloaded Lat/Longs and then how to transform them in same CRS as ABS and then how to determine which address goes to which sa1?


As you can see in this answer: Which CRS to use for Google Maps?, Google Maps uses EPSG 4326. In proj4 format is +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs. So, in R the procedure for any projection is:


wgs84 <- '+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs '

pts <- cbind(c(-29.333,-29.444,-29.555), c(148.555,148.666,148.777)) # dummy point

pts_ <- SpatialPoints(pts, proj4string = CRS(wgs84))

pts_repro <- spTransform(pts_,CRS("+proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +no_defs"))

But, the output coordinates are almost the same:

## SpatialPoints:
##      coords.x1 coords.x2
## [1,]   -29.333   148.555
## [2,]   -29.444   148.666
## [3,]   -29.555   148.777
## Coordinate Reference System (CRS) arguments: +proj=longlat +ellps=GRS80 +towgs84=0,0,0,0,0,0,0
## +no_defs 

coordinates(pts_repro) == coordinates(pts_)
##      coords.x1 coords.x2
## [1,]      TRUE     FALSE
## [2,]     FALSE      TRUE
## [3,]      TRUE      TRUE

Both ellipsoids are pretty similar.

If this doesn't work, maybe coordinates are in web mercator (if you don't post examples coordinate we can help you more)


I offer an ASGS R package that offers this function.


latlon2SA(-33.8688, 151.21, to = "SA2", yr = "2016", return = "v")
#> [1] Sydney - Haymarket - The Rocks
#> 2310 Levels: Abbotsford Aberfoyle Park Acacia Gardens ACT - South West ... Zillmere

It's far too large for CRAN, so to install:

# Installation of dependencies:
install.packages(c("rgdal", "data.table", "sp", "spdep",
                   "dplyr", "leaflet", "htmltools", "magrittr"))

# Download the tar ball (in case the install fails)
# N.B. tarball is about 1 GB
tempf <- tempfile(fileext = ".tar.gz")
download.file(url = "https://dl.dropbox.com/s/zmggqb1wmmv7mqe/ASGS_0.4.0.tar.gz",
              destfile = tempf)

install.packages(tempf, repos = NULL, type = "source")

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.