[R-sig-Geo] CRS: list of datum?
Roger Bivand
Roger.Bivand at nhh.no
Mon Oct 11 23:01:03 CEST 2010
On Mon, 11 Oct 2010, Agustin Lobo wrote:
> Casper,
>
> Thanks, this is very useful, although does not solve my problem.
The list of datums known to proj4 is given by projInfo("datum") - untried,
please try yourself. If they are not given, it is because they are
ambiguous, so proj4 sensibly errs on the side of caution. You have to give
the towgs84= values yourself for your area.
Roger
> What I would need is to define:
>
> proj4string(lakes) <- CRS("+proj=utm +zone=35 +datum=ED87
> +towgs84=-82.981,-99.719,-110.709,-0.5076,0.1503,0.3898")
>
> according to
> http://www.asprs.org/resources/grids/10-99-norway.pdf
>
> My question is why +datum=WGS84 works but +datum=ED87 or +datum=ED50
> do not work?
>
> Using A <- make_EPSG() I've been able to find:
>
>> A$note[grep("ED",A$note)]
>> A[160,]
> code note prj4
> 160 4231 # ED87 +proj=longlat +ellps=intl +no_defs
>
> but this is a lonlat projection, what I need is a UTM projection.
>
> I've also found:
>
>> A[2540,]
> code note
> prj4
> 2540 23035 # ED50 / UTM zone 35N +proj=utm +zone=35 +ellps=intl
> +units=m +no_defs
>
> and have used:
>
>> proj4string(lakes) <-CRS(paste("+init=epsg:23035",
> "+towgs84=-82.981,-99.719,-110.709,-0.5076,0.1503,0.3898"))
>> lakesll <- spTransform(lakes, CRS("+proj=longlat +datum=WGS84"))
>> writeOGR(lakesll,dsn="Kolalakesll",layer="Kolalakesll",driver="ESRI Shapefile")
>
> but lakesll is shifted if overlaid on an image with lonlat and WGS84
> "projection" and
> I suspect the reason is using ED50 instead of ED87.
>
> Note1. lakes comes from:
>
> library(StatDA)
> library(maptools)
> xy <- kola.background$lakes
> names(xy) <- c("x", "y")
> lakes <- map2SpatialLines(xy)
> lakesdata <- data.frame(ID=1:131,lakename=rep("",131),stringsAsFactors=F)
> row.names(lakesdata) <- 1:131
> lakes <- SpatialLinesDataFrame(sl=lakes, data=lakesdata, match.ID = TRUE)
>
> Note2: sessionInfo()
> R version 2.11.1 (2010-05-31)
> i486-pc-linux-gnu
>
> locale:
> [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
> LC_TIME=en_US.UTF-8
> [4] LC_COLLATE=en_US.UTF-8 LC_MONETARY=en_US.UTF-8
> LC_MESSAGES=en_US.UTF-8
> [7] LC_PAPER=en_US.UTF-8 LC_NAME=en_US.UTF-8
> LC_ADDRESS=en_US.UTF-8
> [10] LC_TELEPHONE=en_US.UTF-8 LC_MEASUREMENT=en_US.UTF-8
> LC_IDENTIFICATION=en_US.UTF-8
>
> attached base packages:
> [1] stats graphics grDevices utils datasets methods base
>
> other attached packages:
> [1] rgdal_0.6-27 raster_1.1.2 sp_0.9-62 rkward_0.5.4
>
> loaded via a namespace (and not attached):
> [1] grid_2.11.1 lattice_0.19-11 tools_2.11.1
>
> Thanks!
> Agus
>
> 2010/10/11 caspar hallmann <caspar.hallmann at gmail.com>:
>> HI Agus,
>>
>> see
>> A=make_EPSG()
>> A
>>
>> in rgdal
>>
>> Caspar
>>
>>
>>
>>
>>
>> On Sun, Oct 10, 2010 at 10:36 PM, Agustin Lobo <alobolistas at gmail.com> wrote:
>>> HI!
>>>
>>> While
>>>> CRS("+proj=utm +zone=32 +datum=WGS84")
>>> CRS arguments:
>>> +proj=utm +zone=32 +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
>>>
>>> is correct, the same with ED50 is not:
>>>
>>>> CRS("+proj=utm +zone=32 +datum=ED50")
>>> Error in CRS("+proj=utm +zone=32 +datum=ED50") :
>>> unknown elliptical parameter name
>>>
>>> Is there any way to get the list of available datum acronyms?
>>>
>>> Thanks
>>>
>>> Agus
>>>
>>> _______________________________________________
>>> 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
>>
>
> _______________________________________________
> 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
More information about the R-sig-Geo
mailing list