data projection in R using package SP

Question

I am a bit stuck trying to project my GPS data in R.

I am going through the tutorial in the vignette for the R package 'T-LoCoH', and using those instructions to try to projecting my data from WGS84 to UTM37 (the zone for Ethiopia).

My dataset is called 'd', and the co-ordinates are in columns called 'Longitude' and 'Latitude'.

require(sp)

head(d)
  Latitude Longitude
1  6.36933  39.84300
2  6.37050  39.84417
3  6.41733  39.83800
4  6.41750  39.83800
5  6.41717  39.83750
6  6.41767  39.83733

d2 <- SpatialPoints(d, proj4string=CRS("+proj=longlat, ellps=WGS84"))
Error in CRS("+proj=longlat \nellps=WGS84") : unknown projection id

Can you let me know what I'm doing wrong please?

Also, if I have more columns in my dataset, how do I specify that I just want to project the data in the 2 columns containing Latitude and Longitude?

I would be very grateful for your help.

Thanks, Lindy

you will need to use spTransform() – MLavoie May 02 '16 at 17:01 — MLavoie, May 02 '16 at 17:01
@MLavoie yeah, once you get the data with the right CRS! – Spacedman May 02 '16 at 17:52 — Spacedman, May 02 '16 at 17:52

score 3 · Answer 1 · answered May 02 '16 at 17:46

Instead of:

> d2 <- SpatialPoints(d, proj4string=CRS("+proj=longlat, ellps=WGS84"))
Error in CRS("+proj=longlat, ellps=WGS84") : unknown projection id

You need:

> d2 <- SpatialPoints(d, proj4string=CRS("+proj=longlat +ellps=WGS84"))

But better would be to use the epsg code for GPS coordinates:

> d2 <- SpatialPoints(d, proj4string=CRS("+init=epsg:4326"))

Which is shorthand for this:

Coordinate Reference System (CRS) arguments: +init=epsg:4326
+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0

If you just want to pass the two columns you can do:

> d2 <- SpatialPoints(d[,c("Longitude","Latitude")], proj4string=CRS("+init=epsg:4326"))

which makes me think you should probably have the X and Y coordinates in the other order. You could use any other method for selecting data frame columns, so if your Lat-long are in columns 7 and 23, then d[,c(23,7)] will do it.

score 0 · Answer 2 · answered May 02 '16 at 17:24

You can use the PBSmapping library to do this for you.

require(PBSmapping)

# re-create your data
d = as.data.frame(list(Latitude = c(6.36933,6.37050,6.41733,
                                6.41750,6.41717,6.41767),
                   Longitude = c(39.84300,39.84417,39.83800,
                                 39.83800,39.83750,39.83733)))

# convUL function requires columns to be called "X" and "Y"
dataLL=as.data.frame(list(X=d$Longitude,Y=d$Latitude))

# add a projection attribute (see details in ?convUL)
attr(dataLL,"projection")="LL"

# create a new data frame called dataUTM
dataUTM=convUL(dataLL,km=F)

# output from running previous command
convUL: For the UTM conversion, automatically detected zone 37.
convUL: Converting coordinates within the northern hemisphere.

# extract UTM zone
attr(dataUTM,"zone")

data projection in R using package SP

2 Answers2

Linked