0

Given a zoo time series with columns "lat", "lon" and "value":

z <- zoo(...)
# anyone knows how to create a simple example series here?

I need to convert the lat and lon degrees to UTM in order to perform further operations like distance calculation. After a bit research I came up with a working snippet that converts everything to a data.frame:

d <- as.data.frame(z)
coordinates(d) <- ~ lon+lat
proj4string(d) <- CRS("+proj=longlat +datum=WGS84")
d.final <- spTransform(d, CRS(paste("+proj=utm +zone=",zone," ellps=WGS84",sep="")))

This removes the initial lat and lon columns and adds two new UTM-like lat and lon columns. Is there any simple way to do this directly with the zoo object?

Furthermore, this snippet clearly lacks information from the variable zone. From Determining UTM zone (to convert) from longitude/latitude we learn how to get the zone number from the longitude:

zone.number <- (floor((z$lon + 180)/6) %% 60) + 1

But I still do not now how to calculate the zone letter, eg. "N". How to obtain that? In general, I can expect the zone for all my time series to be constant.

This solution modified from an answer was what was used:

# coming from a zoo time series z with columns z$lat and z$lon
library(rgdal)
# assume that z stays within a zone
zone <- (floor((z$lon[1] + 180)/6) %% 60) + 1
# convert
utm <- project(merge(z$lon, z$lat), paste0("+proj=utm +zone=", zone))
# assign UTM values to new columns in z
z$utmx <- utm[,1]
z$utmy <- utm[,2]
Community
  • 1
  • 1
Martin
  • 1,395
  • 1
  • 11
  • 33

2 Answers2

1

Try this:

library(zoo)
library(rgdal)
z <- zoo(cbind(lat = 12:14, lon = 10:12, value = 1:3)) # test data

zone <- (floor((z$lon[1] + 180)/6) %% 60) + 1
cbind(utm = project(z[, 2:1], paste0("+proj=utm +zone=", zone)), value = z[, 3])

giving:

     utm.1   utm.2 value
1 608864.2 1326751     1
2 716922.6 1437988     2
3 824104.0 1549780     3

where utm.1 and utm.2 are x and y (or easting and northing) respectively. Note that the assumptions as to positive/negative are as in

http://home.hiwaay.net/~taylorc/toolbox/geography/geoutm.html

and that the above R code output double checks to that page. It also double checks to the output in the link in the comments if you look at the standard UTM output on that page.

Update: corrections and modified example.

G. Grothendieck
  • 254,981
  • 17
  • 203
  • 341
  • How can the conversion be independent of the zone letter? I expect the curvature shape of the UTM geometry to be related to latitude much more than longitude. – Martin Sep 17 '14 at 07:52
  • It double checks for me. Please provide an example that it does not double check to. – G. Grothendieck Sep 17 '14 at 12:15
  • I just found out that `project` needs "lon" first, than "lat", which is contradicting and I can not see documented clearly. So this simple one-liner converts latlon to UTM in R: `project(t(c(lon,lat)), paste0("+proj=utm +zone=",zone)`. – Martin Sep 17 '14 at 14:16
  • That's why its `z[, 2:1]` in the answer and not `z[, 1:2]`. – G. Grothendieck Sep 17 '14 at 15:54
-1

This post is basically a restatement of another answer and has been moved to the bottom of the question.

G. Grothendieck
  • 254,981
  • 17
  • 203
  • 341
Martin
  • 1,395
  • 1
  • 11
  • 33
  • (1) It creates an UTM object from a zoo time series, with which one can do anything: replace colums, append, use elsewhere, ... (2) It uses an arbitrary zoo series where column numbers and order is not known. – Martin Sep 17 '14 at 17:45
  • People have a large amount of columns in their time serieses, which change so often, that you never know at which position your column is. That is why column names exist. The "value" column is just an example. In this answer here, there could be hundreds of columns with different values, and you can replace or add an UTM column without caring about the other columns. That is a clear advantage to your solution. Besides, I upvoted your answer hours ago and stressed that credits should go to you. I definitely do not understand your emotions here. – Martin Sep 17 '14 at 18:03