r-tmap / tmap

R package for thematic maps
https://r-tmap.github.io/tmap
GNU General Public License v3.0
865 stars 121 forks source link

Feature request: have tm_rgb able to work from 3 color bands whitout color table #321

Closed prosoitos closed 5 years ago

prosoitos commented 5 years ago

I am new to GIS and this issue may simply be a consequence of my lack of knowledge. My apologies if this is the case.

I have .tif images with 3 bands (for red, green, and blue). RasterBrick objects created from these files with raster::brick() are used successfully by raster::plotRGB(). However, when trying to use them with tmap::tm_rgb() (tm_shape(<my brick object>) + tm_rgb()), I get the error message:

Error: Raster object does not have a color table, nor numeric data that can be converted to colors. Use tm_raster to visualize the data.

And indeed, there is no color table (raster::colortable() on such object returns logical(0)).

So, why is raster::plotRGB() able to create an image without a color table, simply from the color bands in the RasterBrick object while tmap::tm_rgb() is not?

Is there a way to create the missing color table from the bands? I looked at the tmap manual, but could not find any information on this. Did I totally miss it somewhere? It may be something utterly obvious to anyone with a little GIS background and this may be a very silly question.

Is it possible to change tmap::tm_rgb() to make it able to work from the kind of RasterBricks raster::plotRGB() seems to be happy with?

I would like to have all my maps produced with tmap, but for maps that require an rgb image basemap, I am still unable to do so.

mtennekes commented 5 years ago

Could you share the tif file? (Can also in private)

prosoitos commented 5 years ago

It is very heavy...

In case it might be enough, here is the output of running str() on the BrickLayer object created by raster::brick():

Formal class 'RasterBrick' [package "raster"] with 12 slots
  ..@ file    :Formal class '.RasterFile' [package "raster"] with 13 slots
  .. .. ..@ name        : chr "po_118030_rgb_0010002.tif"
  .. .. ..@ datanotation: chr "INT2U"
  .. .. ..@ byteorder   : chr "little"
  .. .. ..@ nodatavalue : num -Inf
  .. .. ..@ NAchanged   : logi FALSE
  .. .. ..@ nbands      : int 3
  .. .. ..@ bandorder   : chr "BIL"
  .. .. ..@ offset      : int 0
  .. .. ..@ toptobottom : logi TRUE
  .. .. ..@ blockrows   : int [1:3] 1 1 1
  .. .. ..@ blockcols   : int [1:3] 14201 14201 14201
  .. .. ..@ driver      : chr "gdal"
  .. .. ..@ open        : logi FALSE
  ..@ data    :Formal class '.MultipleRasterData' [package "raster"] with 14 slots
  .. .. ..@ values    : logi[0 , 0 ] 
  .. .. ..@ offset    : num 0
  .. .. ..@ gain      : num 1
  .. .. ..@ inmemory  : logi FALSE
  .. .. ..@ fromdisk  : logi TRUE
  .. .. ..@ nlayers   : int 3
  .. .. ..@ dropped   : NULL
  .. .. ..@ isfactor  : logi [1:3] FALSE FALSE FALSE
  .. .. ..@ attributes:List of 3
  .. .. .. ..$ : NULL
  .. .. .. ..$ : NULL
  .. .. .. ..$ : NULL
  .. .. ..@ haveminmax: logi TRUE
  .. .. ..@ min       : num [1:3] 0 0 0
  .. .. ..@ max       : num [1:3] 65535 65535 65535
  .. .. ..@ unit      : chr ""
  .. .. ..@ names     : chr [1:3] "po_118030_rgb_0010002.1" "po_118030_rgb_0010002.2" "po_118030_rgb_0010002.3"
  ..@ legend  :Formal class '.RasterLegend' [package "raster"] with 5 slots
  .. .. ..@ type      : chr(0) 
  .. .. ..@ values    : logi(0) 
  .. .. ..@ color     : logi(0) 
  .. .. ..@ names     : logi(0) 
  .. .. ..@ colortable: chr(0) 
  ..@ title   : chr(0) 
  ..@ extent  :Formal class 'Extent' [package "raster"] with 4 slots
  .. .. ..@ xmin: num 714019
  .. .. ..@ xmax: num 728220
  .. .. ..@ ymin: num 8127104
  .. .. ..@ ymax: num 8133054
  ..@ rotated : logi FALSE
  ..@ rotation:Formal class '.Rotation' [package "raster"] with 2 slots
  .. .. ..@ geotrans: num(0) 
  .. .. ..@ transfun:function ()  
  ..@ ncols   : int 14201
  ..@ nrows   : int 5950
  ..@ crs     :Formal class 'CRS' [package "sp"] with 1 slot
  .. .. ..@ projargs: chr NA
  ..@ history : list()
  ..@ z       : list()

and this is the output of calling the BrickLayer object itself:

class      : RasterBrick 
dimensions : 5950, 14201, 84495950, 3  (nrow, ncol, ncell, nlayers)
resolution : 1, 1  (x, y)
extent     : 714018.9, 728219.9, 8127104, 8133054  (xmin, xmax, ymin, ymax)
crs        : NA 
source     : po_118030_rgb_0010002.tif 
names      : po_118030_rgb_0010002.1, po_118030_rgb_0010002.2, po_118030_rgb_0010002.3 
min values :                       0,                       0,                       0 
max values :                   65535,                   65535,                   65535 
mtennekes commented 5 years ago

Quite difficult to see what's going on without able to reproduce it.

You could try to upload the tif file to a site like https://www.filedropper.com/ What also could be an option is to crop the raster object (with the crop function) to make it smaller.

prosoitos commented 5 years ago

OK. I wasn't sure whether the structure of the file would be sufficient or not. I will raster::crop() a small section of the file and see if I can upload it here.

prosoitos commented 5 years ago

sample_brick

prosoitos commented 5 years ago

This is a .rds file. But GitHub doesn't allow for this type of file to be uploaded, so, as a workaround, I renamed it as a .png file.

If you download it and replace .png by .rds, I think you should have a non-corrupted rds file.

Let me know if this doesn't work and I will try something else.

prosoitos commented 5 years ago

To get the BrickLayer object contained in this rds file, I:

  1. read my .tif file with raster::brick()
  2. gave it a projection as it was missing with:
    projection(sample_bri) <- "+proj=utm +zone=6 +south +datum=WGS84 +units=m +no_defs"
  3. cropped it with raster::crop()

I don't think that any of these steps is affecting the problem at hand, but at least you know where this is coming from.

mtennekes commented 5 years ago

This

tm_shape(x) + 
    tm_rgb(max.value = 65535)

gives

Rplot

Is this what you expected?

prosoitos commented 5 years ago

Oh! Yes! Great! Thank you!!

Would you mind explaining to me how you figured that out and how this works? I can't find anything in the tmap man pages that helps here.

prosoitos commented 5 years ago

I can see the min (0 3 times) and max (65535 3 times) values both in the object and in its structure.

Are those .tif files weird (I didn't make them and I don't know how they were made) and is this a bit of an unusual thing?

[max.value is not even listed in the tm_rgb() man page and must be part of the ... arguments. I guess if this were a classic scenario, this would be documented???] (edit: I should have installed the devel version before writing this!)

In any event, thank you! raster::plotRGB() was not a great workaround as I built all the layers I want to add to this with tmap and it is a lot easier to have it all plotted with your great package :slightly_smiling_face:

prosoitos commented 5 years ago

I inspected all the RasterBrick objects created from all my .tif files and they all have the same min and max values, except for one.

That one has 3 different max values... What would be the proper way to handle this? I could do trial and error (like passing a vector with all 3 values to the max.value argument for instance, or trying one of those values as they are not hugely far apart, etc.) until I find something that works, but since this is all very unclear to me, I thought that I would ask. I'd like to understand tm_rgb() better.

Thank you!

And here is that different RasterBrick object:

class      : RasterBrick 
dimensions : 7520, 12759, 95947680, 3  (nrow, ncol, ncell, nlayers)
resolution : 1, 1  (x, y)
extent     : 735667.9, 748426.9, 8116643, 8124163  (xmin, xmax, ymin, ymax)
crs        : +proj=utm +zone=6 +south +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
source     : po_118030_rgb_0050003.tif 
names      : po_118030_rgb_0050003.1, po_118030_rgb_0050003.2, po_118030_rgb_0050003.3 
min values :                       0,                       0,                       0 
max values :                   64186,                   65535,                   65328
prosoitos commented 5 years ago

(edit: ignore this comment: I should have installed the devel version before writing this!)

Euh... I actually tried running

tm_shape(x) + 
  tm_rgb(max.value = 65535)

and I get:

Error in tm_raster(alpha = NA, saturation = 1, interpolate = TRUE, max.value = 65535) : 
  unused argument (max.value = 65535)

Not sure why the code that worked for you is not working on my end...

I am running tmap 2.2 which should be the latest on CRAN. I did not try the devel version. I suspect that is the explanation of the difference?

prosoitos commented 5 years ago

OK. Yes, I tried the devel version and I now have it working. Sorry for the previous comment. I should have checked first.

Oh, and the devel version also has a max.value argument to tm_rgb()... I definitely should have installed it before writing here... Sorry about that.

prosoitos commented 5 years ago

When omitting the max.value in tm_rgb(), the new tmap version also has a much more informative warning message instead of the error message I used to get :slightly_smiling_face:

And the function seems to work anyway now (while it wasn't previously). Except that in my case of course almost all my values are out of range, so I am getting an empty map.

Warning message:
Raster values found that are outside the range [0, 255] 

Thanks!

prosoitos commented 5 years ago

And that answers my previous question about the RasterBrick which has 3 different max values: I have to put the largest of the 3 or I am getting a warning.

With the devel version, things make much more sense!

mtennekes commented 5 years ago

Nice to hear. Are all issues solved now? BTW: it is only possible to specify one value for max.value. Otherwise you'll need to rescale certain raster layers manually.

prosoitos commented 5 years ago

Are all issues solved now?

Yes! Thank you very much @mtennekes. I will close this issue.

BTW: it is only possible to specify one value for max.value. Otherwise you'll need to rescale certain raster layers manually.

OK, thank you for confirming. Is there any side-effect when the brick contains several different max values? In the case I have, I can't see any, but the values are very close. I was wondering if there may be a problem if the values are much further apart.