[R-sig-Geo] LL to UTM

[Gabor Grothendieck]

Great!  Your LL2UTM function was exactly what I wanted.

My data sets are each about a city block in size and I mostly work
in UTM on this project but I do get a few points in lat long. Its small
enough that I have been using that web site I mentioned but it gets
tedious after a while and this should be a big help.

I did notice that the R function below is consistent with the web site
output only to 5 digits after the decimal place in x and 2 digits after
the decimal place in y.  I don't know how to determine whether R or
the web site is the more accurate but at any rate either seem likely
to be well within the accuracy I need since, in truth, I don't need
high precision.

Thanks.



On 3/6/07, Roger Bivand <Roger.Bivand at nhh.no> wrote:
> On Tue, 6 Mar 2007, Gabor Grothendieck wrote:
>
> > Thanks. I looked at spTransform prior to posting but don't understand it.
> > What I want is to create a function that returns UTM coordinates:
> >
> > LL2UTM <- function(lat, long, zone = 18) {
> > ...
> > }
>
> library(rgdal)
> LL2UTM <- function(lat, long, zone=18) {
>  project(cbind(long, lat), paste("+proj=utm +zone=", zone, sep=""))
> }
>
> gets you there, for example for:
>
> lat <- seq(30,45,1)
> long <- seq(-78,-72,1)
> grd <- expand.grid(lat, long)
> LL2UTM(grd[,1], grd[,2])
>
> The project() interface to PROJ4 does not support datum transfromation,
> the spTransform() method does. The cost is having to convert the data to
> a Spatial* object:
>
> SP_grd <- SpatialPoints(cbind(grd[,2], grd[,1]),
>  proj4string=CRS("+proj=longlat +datum=NAD83"))
> SP_grd1 <- spTransform(SP_grd, CRS("+proj=utm +zone=18 +datum=NAD83"))
>
> returning a SpatialPoints object - use the coordinates() method to
> retrieve as a matrix.
>
> Assuming that your input geographical coordinates are in NAD84/WGS84, and
> the output UTM coordinates are in the same datum, project() will be OK. If
> you need datum transformation, which may lead to errors of hundreds of
> metres if not used, the PROJ4 string needs more detail. If this is at the
> block scale, datum transformation will make a difference if the input and
> output specifications vary (say placing points on a UTM map not in
> WGS84/NAD83). At the continental scale the differences are typically not
> great, but then you wouldn't use UTM anyway.
>
> Hope this helps,
>
> Roger
>
> >
> > If necessary, zone=18 can be hardcoded in the function and that
> > arg removed.
> >
> > I assume that using spTransform its just a one line body. Can
> > you give me the specific line that it should be?
> >
> > Thanks.
> >
> >
> > On 3/6/07, Edzer J. Pebesma <e.pebesma at geo.uu.nl> wrote:
> > > Gabor,
> > >
> > > package rgdal provides an interface to the PROJ.4 library for projection
> > > of geographical data. Look for the function (or rather method)
> > > spTransform. It takes any of the spatial classes provided by package sp,
> > > and afaik any of the known projection systems, including UTM and
> > > ellipsoids (LL reference "model").
> > > --
> > > Edzer
> > >
> > > Gabor Grothendieck wrote:
> > > > I am currently using this web page to convert LL to UTM (and UTM to LL).
> > > >
> > > >    http://home.hiwaay.net/~taylorc/toolbox/geography/geoutm.html
> > > >
> > > > All my points are in zone 18.  How would I accomplish the same
> > > > thing using R?
> > > >
> > > > Thanks.
> > > >
> > > > _______________________________________________
> > > > R-sig-Geo mailing list
> > > > R-sig-Geo at stat.math.ethz.ch
> > > > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> > > >
> > >
> > >
> >
> > _______________________________________________
> > R-sig-Geo mailing list
> > R-sig-Geo at stat.math.ethz.ch
> > https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> >
>
> --
> Roger Bivand
> Economic Geography Section, Department of Economics, Norwegian School of
> Economics and Business Administration, Helleveien 30, N-5045 Bergen,
> Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
> e-mail: Roger.Bivand at nhh.no
>
>