2

Here is an updated sense of where I'm at. I decided to re-read in my data file using ReadOGR.

MyDataFile <- read.csv("reduced1.csv")
class(MyDataFile)
coordinates(MyDataFile)<-~X+Y
class(MyDataFile)

writeOGR(MyDataFile, "Folder", "File", driver = "ESRI Shapefile")
proj4string(MyDataFile)

Checking the CRS at this point gives me 'NA' so I try and assign it

Redefined <- readOGR(dsn="Folder", layer = "File")
proj4string(Redefined) <- CRS("+init=epsg:4326") # WGS 84

[1] "+init=epsg:4326 +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"

Now for reading in my boundaries file:

Map <-readOGR(".","infuse_dist_lyr_2011")
proj4string(Map)

[1] "+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.9996012717 +x_0=400000 +y_0=-100000 +datum=OSGB36 +units=m +no_defs +ellps=airy +towgs84=446.448,-125.157,542.060,0.1502,0.2470,0.8421,-20.4894"

So now to transform.

NewMap <- spTransform(Map, proj4string(Redefined))
  • The CRS function doesn't apply to sp-class objects, its applied to strings to create CRS objects. Use proj4string(map) or raster::projection(map) to get the projection string of an sp-class object. That's where your warning comes from. – Spacedman Feb 15 '17 at 13:48
  • And does CRS(map) really give an error about CRS(h)? We care about precision here... – Spacedman Feb 15 '17 at 13:50
  • Apolgies you're correct - I changed the name to map for the purposes of this post but omitted to do so consistently - sorry! – Henry Cann Feb 15 '17 at 13:56
  • Also proj4string returns strings which seem to suggest my files are in the right CRS. However the exact string is slightly different each time - wondering if this matters? I.e for my datafile "+proj=utm +zone=10 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0" whereas for my shapefile +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0" - the latter being shorted and without the zone and unit bits. – Henry Cann Feb 15 '17 at 13:58
  • 1
    That sounds like you've got the wrong proj4string assigned to some data. If you read with readOGR or raster::shapefile it should pick up the projection from the .prj file. In your rewrite, show the output of summary on your spatial objects. – Spacedman Feb 15 '17 at 17:58
3

Try reading the layers using rgdal::readOGR , this reads the projection automatically and is (IMO) more robust that readShapePoly.

Then if needed, run spTransform with one of your layers as the "to" option.

  • Tom, I will try readOGR for my shapefile but what about my datafile initially read in from csv? Can readOGR be used? – Henry Cann Feb 15 '17 at 13:44
  • Just to add that I followed the advice of using ReadOGR with both files. I then was able to use sp Transform into epsg:4238 for both. However for the datafile, which already had the right coordinates, the coordinate column has now gone awry. E.g the first entry reads -127.4887 0.00042064 where it used to be in degrees decimal – Henry Cann Feb 15 '17 at 14:33
  • I'm slightly confused over what your problem is now. The co-ordinates will change to match to new projection system. Are the data appearing where you would expect when plotted? – Tom Higginbottom Feb 16 '17 at 15:16

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.