r.out.gdal: add alpha channel support

[neteler]

At time the export of RGB data with nodata is rather difficult.

grass72 ~/grassdata/nc_spm_08_grass7/landsat

# test case:
# set region to be half off the NC Landsat scene:
g.region n=226831.5 s=221188.5 w=626430 e=634467 res=28.5 -p
d.mon wx0
d.rgb b=lsat7_2002_10 g=lsat7_2002_20 r=lsat7_2002_30

# yes, null is expected:
r.univar lsat7_2000_30 -g | grep null
null_cells=24750

r.out.gdal lsat7_2000 output=lsat7_2000_RGB.tif --o
Checking GDAL data type and nodata value...
 100%
Using GDAL data type <Byte>
Input raster map contains cells with NULL-value (no-data). The value 0 will
be used to represent no-data values in the input map. You can specify a
nodata value with the nodata option.
 100%
...

But in this case the 0 value is a data value and should not be used for nodata. The only solution seems to be an extra alpha channel to be stored in the resulting GDAL file (here: GeoTIFF).

[mmetz]

Replying to neteler:

At time the export of RGB data with nodata is rather difficult.

grass72 ~/grassdata/nc_spm_08_grass7/landsat

# test case:
# set region to be half off the NC Landsat scene:
g.region n=226831.5 s=221188.5 w=626430 e=634467 res=28.5 -p
d.mon wx0
d.rgb b=lsat7_2002_10 g=lsat7_2002_20 r=lsat7_2002_30

# yes, null is expected:
r.univar lsat7_2000_30 -g | grep null
null_cells=24750

r.out.gdal lsat7_2000 output=lsat7_2000_RGB.tif --o
Checking GDAL data type and nodata value...
 100%
Using GDAL data type <Byte>
Input raster map contains cells with NULL-value (no-data). The value 0 will
be used to represent no-data values in the input map. You can specify a
nodata value with the nodata option.
 100%
...

But in this case the 0 value is a data value and should not be used for nodata.

According to Landsat metadat, 0 is not a valid data value, valid range is 1 - 255. Also, r.info -r shows that 0 is not present, therefore 0 can be used for nodata.

The only solution seems to be an extra alpha channel to be stored in the resulting GDAL file (here: GeoTIFF).

You mean r.out.gdal should create this extra alpha channel internally? You could also try the PHOTOMETRIC and ALPHA creation options for GeoTIFF with a corresponding alpha channel created manually.

[neteler]

Replying to mmetz:

Replying to neteler:

...

According to Landsat metadat, 0 is not a valid data value, valid range is 1 - 255. Also, r.info -r shows that 0 is not present, therefore 0 can be used for nodata.

Well, this was just intended as simple example to reproduce the issue. Indeed I am using different data. At least there is no guarantee that 0 is never a valid dataset while having 8 bit to be exported.

The only solution seems to be an extra alpha channel to be stored in the resulting GDAL file (here: GeoTIFF).

You mean r.out.gdal should create this extra alpha channel internally?

Yes, I thought so.

You could also try the PHOTOMETRIC and ALPHA creation options for GeoTIFF with a corresponding alpha channel created manually.

Interesting - do you have an example? Would be happy to close this ticket as "worksforme" :)

[mmetz]

Replying to neteler:

Replying to mmetz:

Replying to neteler:

...

According to Landsat metadat, 0 is not a valid data value, valid range is 1 - 255. Also, r.info -r shows that 0 is not present, therefore 0 can be used for nodata.

Well, this was just intended as simple example to reproduce the issue. Indeed I am using different data. At least there is no guarantee that 0 is never a valid dataset while having 8 bit to be exported.

If there are NULL cells in the bands to be exported, a valid value according to the selected output data type needs to be selected to represent nodata, otherwise NULL cells can not be encoded, e.g. for Byte, the value representing nodata must be in the range 0-255. This is why you see "Checking GDAL data type and nodata value..." in the output of r.out.gdal. The module does not always choose 0 as nodata value, this depends on the input bands and the selected output data type.

The only solution seems to be an extra alpha channel to be stored in the resulting GDAL file (here: GeoTIFF).

You mean r.out.gdal should create this extra alpha channel internally?

Yes, I thought so.

You could also try the PHOTOMETRIC and ALPHA creation options for GeoTIFF with a corresponding alpha channel created manually.

Interesting - do you have an example? Would be happy to close this ticket as "worksforme" :)

You would need to create an alpha channel that encodes NULL cells and in the bands to be exported replace NULL cells with some value in the range 0-255. For example

RMAP="lsat7_2000_30"
GMAP="lsat7_2000_20"
BMAP="lsat7_2000_10"

OUTNAME="lsat7_2000.tif"

# extract alpha
r.mapcalc "out_a = if(isnull($RMAP) || isnull($GMAP) || isnull($BMAP), 0, 255)"

# replace NULL cells with a valid value, extract colors

# extract red
r.mapcalc "out_r = if(isnull($RMAP), 0, #$RMAP)"
# extract green
r.mapcalc "out_g = if(isnull($GMAP), 0, #$GMAP)"
# extract blue
r.mapcalc "out_b = if(isnull($BMAP), 0, #$BMAP)"

# create group
i.group group=out_rgba input=out_r,out_g,out_b,out_a

# remove any MASK because this works only if there are 
# no NULL cells in the bands to be exported 
r.mask -r

# export the group
# add PROFILE=BASELINE to createopt to produce a standard TIFF file without out GTiff extensions
r.out.gdal in=out_rgba out=$OUTNAME -cm createopt="PHOTOMETRIC=RGB,ALPHA=YES"

A nice enhancement to r.out.gdal would be to accept multiple input raster maps, this would avoid the need to create a group just for exporting multiple raster maps into a single GDAL raster dataset.

[neteler]

The example works nicely, thanks! If no objections, I'll add it to the manual.

Yet a question:

With Sentinel-2, also 0 is used as a data value (eg Computed Min/Max=0.000,13668.000). The black border mask can be retrieved from an auxiliary JP2000 file and applied as MASK in GRASS GIS. So far so nice.

Still I am not sure if it now better to use r.rescale or the # operator of r.mapcalc with some formula in order to "free up" the 0 value for paving the way of exporting with alpha as above?

[mmetz]

Replying to neteler:

With Sentinel-2, also 0 is used as a data value (eg Computed Min/Max=0.000,13668.000). The black border mask can be retrieved from an auxiliary JP2000 file and applied as MASK in GRASS GIS. So far so nice.

Still I am not sure if it now better to use r.rescale or the # operator of r.mapcalc with some formula in order to "free up" the 0 value for paving the way of exporting with alpha as above?

Using a nodata value means that any band values equal to the nodata value will be interpreted as nodata. Using an alpha channel means that all pixels with an alpha value of 0 are transparent. The alpha channel thus represents per-pixel encoding of nodata, just like the GRASS null file. That means when using an alpha channel, you do not need to "free up" any particular value, you can use any value you like to replace NULL cells, as long as the value can be represented by the Byte data type. It does not matter if that value is already present in any of the input bands.

This is about exporting RGB data, not GIS data, therefore the # operator of r.mapcalc needs to be used. The result is a regular image (like a screenshot) that can be displayed by any simple image viewer coming with your OS and/or desktop environment. When exporting, do you want to provide something like a screenshot or the actual data?

[neteler]

Replying to mmetz:

Replying to neteler:

Still I am not sure if it now better to use r.rescale or the # operator of r.mapcalc with some formula in order to "free up" the 0 value for paving the way of exporting with alpha as above?

This is about exporting RGB data, not GIS data, therefore the # operator of r.mapcalc needs to be used. The result is a regular image (like a screenshot) that can be displayed by any simple image viewer coming with your OS and/or desktop environment. When exporting, do you want to provide something like a screenshot or the actual data?

Yes exactly. Here it is for showing the data in a Web system. I have taken liberty to write up an example for the manual page from your comments, r69162. Thanks for the explanations.

(an example backport needs r69128 and r69131 - they works nicely for me in a local 7.2.svn backport).

[mmetz]

Replying to neteler:

[...] I have taken liberty to write up an example for the manual page from your comments, r69162. Thanks for the explanations.

(an example backport needs r69128 and r69131 - they works nicely for me in a local 7.2.svn backport).

r68128 has already been backported to relbr72 and relbr70 in r69129,30. r69131 is not relevant. Please feel free to backport the example.

[neteler]

Example backported in r69265. Closing.