Geographical CRS given to non-conformant data in R

Question

sorry for bothering you for this, but it is maybe 5 hours I am getting mad on this issue and I don't manage to sort it out.

I have a dataset of about 37,000 entries. Each of them has its own Lon and Lat coordinate values. Checking the overall values, they range respectively as follow: Latitude (-54.4871,70.66344) and Longitude (-177.375, 178.4419). This is absolutely reasonable.

I created a shapefile with these 37 thousands points using ArcGIS: everything works fine.

I need then to process these data using R, the command I used for my code is (maptools package):

cells <- readShapeSpatial('RES',IDvar="id_obj", 
                          proj4string=CRS("+proj=longlat +datum=WGS84"))

but R gives an error:

Error in validityMethod(as(object, superClass)) : Geographical CRS given to non-conformant data: 2.76663393422e+145

(I have no idea on where this number comes from, it is not part of my dataset...)

Reading other posts on this blog seems that the cause should be an invalid data for lon or lat, but as I mentioned above, this is not the case for my dataset.

I tried to create different shapefiles, the first one was not projected, using several projections (WGS84 Mercator, web mercator...), but the error is always the same...

Thanks for your help.

Answer 1

The bottom line is that your shapefile seems to be corrupt.

The points shapefile has two main sections, a coords section that has the coordinates of the points, and a data section that has "attribute" data (information about the points, like region and country in your case). Your shapefile has Lon and Lat in the data section as well, but they don't match:

library(rgdal)
setwd("<directory with shapefile...>")
map <- readOGR(dsn=".", layer="test")

range(map@data$Lat)
# [1] -54.48708  70.66344

range(map@coords[,2])
# [1]  -5.448708e+01  2.766634e+145

Reprojection involves transforming the information in the coords section, which is why it failed.

Here is a workaround, but hacking the SpatialPointsDataFrame is not a good idea:

map@coords <- as.matrix(map@data[c("Lon","Lat")])
map@bbox   <- rbind(range(map@coords[,1]),range(map@coords[,2]))
wgs.84 <- "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(map) <- CRS(wgs.84)

library(ggplot2)
gg <- data.frame(map@coords)
ggplot(gg) + 
  geom_point(aes(x=Lon,y=Lat), size=1, alpha=0.5, colour="blue") + 
  coord_fixed()

mercator <- "+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs"
map.mercator <- spTransform(map,CRS=CRS(mercator))
gg <- data.frame(map.mercator@coords)
ggplot(gg) + 
  geom_point(aes(x=Lon,y=Lat), size=1, alpha=0.5, colour="green") + 
  coord_fixed()

I'd recommend you re-create the shapefile and try again.