[R-sig-Geo] CRS: list of datum?

Agustin Lobo alobolistas at gmail.com
Tue Oct 12 12:08:08 CEST 2010


Roger,

projInfo("datum") works fine and outputs 10 datum,
ED50 and ED87 not listed.

Anyway, looking for an operational solution, is this correct for Norway?

CRS("+proj=utm +zone=35 +ellps=intl +units=m +no_defs
+towgs84=-82.981,-99.719,-110.709,-0.5076,0.1503,0.3898")

(I'm assuming that the Kola data in StatDa follow projection standards
for Norway as in
http://www.asprs.org/resources/grids/10-99-norway.pdf).

because the fact is that the output layer produced by spTransform()
seems to be shifted.

Agus

2010/10/11 Roger Bivand <Roger.Bivand at nhh.no>:
> 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
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>
>



More information about the R-sig-Geo mailing list