2

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?

Esan
  • 83
  • 1
  • 6

2 Answers2

1

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:

library(sp)

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:

ptos_repro
## 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)

aldo_tapia
  • 13,560
  • 5
  • 30
  • 56
0

I offer an ASGS R package that offers this function. 

library(ASGS)

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")
PolyGeo
  • 65,136
  • 29
  • 109
  • 338
Hugh
  • 101
  • 1