[R-sig-Geo] confusing format of shapefile, change coordinate system to long/lat

[Roger Bivand]

On Thu, 11 Jun 2009, Demmler J. wrote:

> Dear R fellow-users,
>
> I'm new to using the "maptools" and "rgdal" libraries in R, so I'm very 
> sorry if this is something really obvious. I'm trying to get an ESRI 
> shape file (Wales lower super output areas) into the right format. I 
> think a "proper" shapefile should just have two columns with the 
> coordinates, so this one is a bit odd. I'm not quite sure what the 
> problem is, if it is the dimensions or the class of the dataset (or 
> both).

Your shapefile contains 1935 observations, and the geometric form 
(support) of the observations is as polygons. As read, they are in the 
OSGB National Grid coordinate reference system. It seems that you believe 
that the data is in geographical coordinates, but it is not. If the two 
attributes GEOEAST and GEONORTH are geographical coordinates, you need to 
extract them, cbind() them to make a two column matrix, then make a 
SpatialPoints object - watch the datum of the CRS(), it will bite - and 
spTransform() that.

But it isn't obvious that you want to do that, or that GEOEAST, GEONORTH 
are in geographical coordinates (or POPEAST, POPNORTH for that matter) - 
do summary() of the object you have read and look at those variables. Find 
out from the data provider what they are, I suspect that the GEOEAST, 
GEONORTH pair is the centroid of the polygon, and that POPEAST, POPNORTH 
is the population-weighted centroid of the polygon for the population 
distribution at a known date.

Hope this helps,

Roger

>
> library(maptools)
> library(rgdal)
>
>> ogrInfo(".","Wales_low_soa_2001")
> Driver: ESRI Shapefile number of rows 1935
> Feature type: wkbPolygon with 2 dimensions
> Number of fields: 7
>      name type length typeName
> 1 GEONORTH    2     31     Real
> 2 POPNORTH    2     31     Real
> 3  GEOEAST    2     31     Real
> 4     NAME    4    254   String
> 5    LABEL    4    254   String
> 6  POPEAST    2     31     Real
> 7 ZONECODE    4    254   String
>
>> wales_ll <- readOGR(".","Wales_low_soa_2001")
> OGR data source with driver: ESRI Shapefile
> Source: ".", layer: "Wales_low_soa_2001"
> with  1935  rows and  7  columns
> Feature type: wkbPolygon with 2 dimensions
>
> I also tried ReadShapePoly, but I can't see any noticable difference.
>
>> class(wales_ll)
> [1] "SpatialPolygonsDataFrame"
> attr(,"package")
> [1] "sp"
>
> I would like to do something like this (I also tried to just use the GEONORTH and GEOEAST columns, but wasn't very successful either):
>
>> proj4string(wales_ll) <- CRS("+proj=longlat ellps=WGS84")
> Error in `proj4string<-`(`*tmp*`, value = <S4 object of class "CRS">) :
>  Geographical CRS given to non-conformant data
>
>> summary(wales_ll)
> Object of class SpatialPolygonsDataFrame
> Coordinates:
>        min      max
> r1 146597.1 355312.8
> r2 164538.4 395994.1
> Is projected: TRUE
> proj4string :
> [+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000
> +y_0=-100000 +ellps=airy +units=m +no_defs]
> etc.
>
>
> and then later I would like to do:
>
>> EPSG <- make_EPSG()
>> EPSG[grep("British National Grid", EPSG$note), 1:2]
>      code                                note
> 2791 27700 # OSGB 1936 / British National Grid
>
>> wales.gb <- spTransform(wales_ll, CRS("+init=epsg:27700"))
>
> and possibly at some point writeOGR to kml...
>
> Any help with this is highly appreciated!
>
> Thanks  Joanne
>
> -------------------------------------------------------------------
> Linux (Ubuntu 9.04), installed from pre-compiled binary
>
>> sessionInfo()
> R version 2.8.1 (2008-12-22)
> i486-pc-linux-gnu
>
> locale:
> LC_CTYPE=en_GB.UTF-8;LC_NUMERIC=C;LC_TIME=en_GB.UTF-8;LC_COLLATE=en_GB.UTF-8;LC_MONETARY=C;LC_MESSAGES=en_GB.UTF-8;LC_PAPER=en_GB.UTF-8;LC_NAME=C;LC_ADDRESS=C;LC_TELEPHONE=C;LC_MEASUREMENT=en_GB.UTF-8;LC_IDENTIFICATION=C
>
> attached base packages:
> [1] stats     graphics  grDevices utils     datasets  methods   base
>
> other attached packages:
> [1] rgdal_0.6-8     maptools_0.7-23 sp_0.9-37       foreign_0.8-30
>
> loaded via a namespace (and not attached):
> [1] grid_2.8.1      lattice_0.17-25
>
> 	[[alternative HTML version deleted]]
>
> _______________________________________________
> 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