TransWikia.com

Wrong projection of raster using R

Geographic Information Systems Asked on June 1, 2021

I have a raster layer with the following structure:

ncols 4320
nrows 2160
xllcorner -180.0
yllcorner -90.0
cellsize 0.0833333
NODATA_value -9999.0
-9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0 -9999.0

If I then try to read it and project, it does not appear:

library(raster)
library(tmap)

r <- raster("baseline/baseline/zip/1870AD_pop/popc_1870AD.asc")
crs(r) <- "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0" 

tmap_mode("view")
tm_shape(r)+tm_raster("red")

For some reason, I’m unable to display it properly, even with the added crs.

plot(r) does give me an expected output:

plot

However, when trying to plot it on another map, it fails to do it correctly.

summary(r) gives the following output:

> summary(r)
         popc_1870AD
Min.    0.000000e+00
1st Qu. 0.000000e+00
Median  7.089942e+00
3rd Qu. 1.621315e+02
Max.    1.331152e+05
NA's    7.092085e+06
Warnmeldung:
In .local(object, ...) :
  summary is an estimate based on a sample of 1e+05 cells (1.07% of all cells)

The file can be downloaded here:
http://www.filedropper.com/popc1870ad

One Answer

Your data goes up/down to 90N/S. On a Mercator projection, like on most web mapping software, this goes to +Infinity/-Infinity.

tmap seems to have even tighter limits on what it can display. If I crop to +/-50 I can see your data within that latitude band:

> rc50 = crop(r,extent(c(-180,180,-50,50)))
> tm_shape(rc50) + tm_raster()

enter image description here

but if I change it to -51 and 51, I see nothing!

If you reproject the data to the web mercator projection, then you can display more, but you still have to crop it first because otherwise you go off to infinity:

> rc = crop(r,extent(c(-180,180,-80,85)))
> rp = projectRaster(rc, crs="+init=epsg:3857")

enter image description here

It doesn't matter if you set the CRS like you do or not, since you are using the default CRS anyway (but its a Good Thing to set it anyway).

Correct answer by Spacedman on June 1, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP