1

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

6
  • What does summary(r) say about it in terms of its extent?
    – Spacedman
    May 26 at 9:03
  • Does a simple plot(r) show anything meaningful? Are you sure its not all -9999 missing data?
    – Spacedman
    May 26 at 9:10
  • please see the eddits
    – Yves
    May 26 at 9:21
  • I think I meant you to show print(r) which gives the full details of the geometry, but I can see from the plot that the extents etc are correct. Can you edit and add the version number of tmap? What do you get with your tmap code? Just a blank web page? Tried tm_shape(r) + tm_raster() (ie without the "red") to map the data?
    – Spacedman
    May 26 at 9:32
  • The tmap version is 3.3-1. without the red, a legend occurs, but the data doesn't seem to be projected at the correct scale. The issue also occurs when I try to import it in ArcGIS Pro.
    – Yves
    May 26 at 9:37
2

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).

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.