On Mon, 4 Jul 2011, Horacio Samaniego wrote:
I need to add a lat/lon axes to an equal area plot of a
SpatialLinesDataFrame object currently in units of meters. Is that
achievable using sp/maptools?
The .prj file where the object comes from is:
PROJCS["Lambert_Azimuthal_Equal_Area",GEOGCS["GCS_unnamed
ellipse",DATUM["D_unknown",SPHEROID["Unknown",6370997,0]],PRIMEM["Greenwich",0],UNIT["Degree",0.017453292519943295]],PROJECTION["Lambert_Azimuthal_Equal_Area"],PARAMETER["latitude_of_origin",45],PARAMETER["central_meridian",-100],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["Meter",1]]
I'ld hate to install arcinfo to generate 1 plot! (plus the plot is
already nice as it is!!! just need the axes!)
No, the axes are in the same coordinate system as the plot, and would be
curves here. You can overplot a long-lat grid on a projected plot:
library(sp)
data(meuse)
coordinates(meuse) <- ~x+y
proj4string(meuse) <- CRS("+init=epsg:28992")
plot(meuse)
library(rgdal)
meuse_ll <- spTransform(meuse, CRS("+proj=longlat +datum=WGS84"))
grd <- gridlines(meuse_ll)
grd_x <- spTransform(grd, CRS("+init=epsg:28992"))
plot(grd_x, add=TRUE, lty=2)
grdat_ll <- gridat(meuse_ll)
grdat_x <- spTransform(grdat_ll, CRS("+init=epsg:28992"))
text(coordinates(grdat_x), labels=parse(text=grdat_x$labels),
pos=grdat_x$pos, offset=grdat_x$offset)
Hope this helps,
Roger
thanks,
I'll summarize!
H
--
Roger Bivand
Department of Economics, NHH Norwegian School of Economics,
Helleveien 30, N-5045 Bergen, Norway.
voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: roger.biv...@nhh.no
_______________________________________________
R-sig-Geo mailing list
R-sig-Geo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo