r.external vs r.in.gdal under Linux using r.his

[Nikos Alexandris]

Get sample raster maps from QGIS' repository

wget 'https://github.com/qgis/QGIS/blob/77725c90b8a5bfe61fde5097128cdf4fb3fce610/python/plugins/processing/tests/testdata/custom/grass7/float_raster.tif?raw=true' -O float_raster.tif
wget 'https://github.com/qgis/QGIS/blob/77725c90b8a5bfe61fde5097128cdf4fb3fce610/python/plugins/processing/tests/testdata/custom/grass7/raster_4class.tif?raw=true' -O raster_4class.tif
wget 'https://github.com/qgis/QGIS/blob/77725c90b8a5bfe61fde5097128cdf4fb3fce610/python/plugins/processing/tests/testdata/custom/grass7/raster_5class.tif?raw=true' -O raster_5class.tif

# Import

r.in.gdal in=raster_4class.tif out=raster_4class
r.in.gdal in=raster_5class.tif out=raster_5class
r.in.gdal input=float_raster.tif out=float_raster

# Link to

r.external in=raster_4class.tif out=raster_4class_ext
r.external in=raster_5class.tif out=raster_5class_ext
r.external input=float_raster.tif out=float_raster_ext

# Query stats

for MAP in float_raster float_raster_ext ;do echo $MAP; r.univar -g $MAP ;echo ;done
for MAP in raster_4class raster_4class_ext ;do echo $MAP; r.univar -g $MAP ;echo ;done
for MAP in raster_5class raster_5class_ext ;do echo $MAP; r.univar -g $MAP ;echo ;done

As expected, the above commands return the same statistics for each pair of imported (r.in.gdal) and linked (r.external) map.

# Convert (supposedly) to RGB

r.his h=float_raster i=raster_4class s=raster_5class r=r g=g bl=b --o
r.his h=float_raster_ext i=raster_4class_ext s=raster_5class_ext r=r_ext g=g_ext bl=b_ext --o

# Query statistics

for MAP in r r_ext ;do echo $MAP; r.univar -g $MAP ;echo ;done
for MAP in g g_ext ;do echo $MAP; r.univar -g $MAP ;echo ;done
for MAP in b b_ext ;do echo $MAP; r.univar -g $MAP ;echo ;done

r
n=14889
null_cells=262
cells=15151
min=0
max=252
range=252
mean=71.4881456108536
mean_of_abs=71.4881456108536
stddev=51.9267840636893
variance=2696.39090319701
coeff_var=72.6369157011754
sum=1064387

r_ext
n=14889
null_cells=262
cells=15151
min=125
max=127
range=2
mean=125.927933373632
mean_of_abs=125.927933373632
stddev=0.76658997256503
variance=0.587660186037254
coeff_var=0.608752920839682
sum=1874941

g
n=14889
null_cells=262
cells=15151
min=0
max=231
range=231
mean=95.9787091141111
mean_of_abs=95.9787091141111
stddev=61.6018003150434
variance=3794.78180205448
coeff_var=64.1827764549362
sum=1429027

g_ext
n=14889
null_cells=262
cells=15151
min=125
max=127
range=2
mean=125.927933373632
mean_of_abs=125.927933373632
stddev=0.76658997256503
variance=0.587660186037254
coeff_var=0.608752920839682
sum=1874941

b
n=14889
null_cells=262
cells=15151
min=0
max=142
range=142
mean=83.9135603465646
mean_of_abs=83.9135603465646
stddev=46.9957861849984
variance=2208.60391914608
coeff_var=56.0049960827605
sum=1249389

b_ext
n=14889
null_cells=262
cells=15151
min=125
max=127
range=2
mean=125.927933373632
mean_of_abs=125.927933373632
stddev=0.76658997256503
variance=0.587660186037254
coeff_var=0.608752920839682
sum=1874941

Look at ranges only:

for MAP in r r_ext ;do echo $MAP; r.info -r $MAP ;echo ;done
r
min=0
max=252

r_ext
min=125
max=127

g
min=0
max=231

g_ext
min=125
max=127

b
min=0
max=142

b_ext
min=125
max=127

The above HIS to RGB conversion is NOT correct due to the necessity in r.his for 8-bit input maps. Nevertheless, using the "same" input maps (i.e.:

float_raster and float_raster_ext, raster_4class and raster_5class, raster_4class_ext and raster_5class_ext)

should, in theory, deliver the same outputs (i.e.:

r r_ext g g_ext b b_ext).

Looking at the numbers, there is an unexpected (?) disagreement in the statistics between the pair of same map (one map imported, one map linked to via r.external).

# With the current r.his, the right way to do this, would be:

for MAP in \
    float_raster float_raster_ext \
    raster_4class raster_4class_ext \
    raster_5class raster_5class_ext ;do \
    r.rescale \
    in=$MAP \
    to=0,255 \
    out=${MAP}_8bit
done

Convert

r.his h=float_raster_8bit i=raster_4class_8bit s=raster_5class_8bit r=r_8bit g=g_8bit bl=b_8bit
r.his h=float_raster_ext_8bit i=raster_4class_ext_8bit s=raster_5class_ext_8bit r=r_ext_8bit g=g_ext_8bit bl=b_ext_8bit

and check

r_8bit
n=14889
null_cells=262
cells=15151
min=6
max=250
range=244
mean=86.576331519914
mean_of_abs=86.576331519914
stddev=42.396468104142
variance=1797.46050770553
coeff_var=48.9700445373919
sum=1289035

r_ext_8bit
n=14889
null_cells=262
cells=15151
min=6
max=250
range=244
mean=86.576331519914
mean_of_abs=86.576331519914
stddev=42.396468104142
variance=1797.46050770553
coeff_var=48.9700445373919
sum=1289035

g_8bit
n=14889
null_cells=262
cells=15151
min=1
max=228
range=227
mean=102.257841359393
mean_of_abs=102.257841359393
stddev=47.2255504199476
variance=2230.25261246701
coeff_var=46.1828157060052
sum=1522517

g_ext_8bit
n=14889
null_cells=262
cells=15151
min=1
max=228
range=227
mean=102.257841359393
mean_of_abs=102.257841359393
stddev=47.2255504199476
variance=2230.25261246701
coeff_var=46.1828157060052
sum=1522517

b_8bit
n=14889
null_cells=262
cells=15151
min=5
max=139
range=134
mean=94.0406340251192
mean_of_abs=94.0406340251192
stddev=38.318460482594
variance=1468.30441375612
coeff_var=40.7467058041726
sum=1400171

b_ext_8bit
n=14889
null_cells=262
cells=15151
min=5
max=139
range=134
mean=94.0406340251192
mean_of_abs=94.0406340251192
stddev=38.318460482594
variance=1468.30441375612
coeff_var=40.7467058041726
sum=1400171

We suspect that something goes wrong with r.external and the range in input maps is not the same as when importing. This is not obvious at first. I.e., r.info -r reports identical ranges between imported and linked maps:

for MAP in $(glr pattern=*raster* exclude=*8bit*) ;do echo $MAP; r.info -r $MAP ;echo ;done
float_raster
min=-303
max=14847

float_raster_ext
min=-303
max=14847

raster_4class
min=1
max=4

raster_4class_ext
min=1
max=4

raster_5class
min=1
max=5

raster_5class_ext
min=1
max=5

obviously and expectedly, ranges are the same for maps rescaled to [0,255]:

float_raster_8bit
min=0
max=255

float_raster_ext_8bit
min=0
max=255

raster_4class_8bit
min=0
max=255

raster_4class_ext_8bit
min=0
max=255

raster_5class_8bit
min=0
max=255

raster_5class_ext_8bit
min=0
max=255

This came after looking at why the r.his related test in QGIS fails (see ​https://travis-ci.org/qgis/QGIS/jobs/603479075#L3608).

1. The test itself there is not correct (as far as I understand, input maps must range in [0,255] for this test to be correct) but this is a different issue.
2. The test fails after Panos' following PR ​https://github.com/qgis/QGIS/pull/32425#.

Panos and Nikos

[mmetz]

The differences are caused by different color rules. From the manual of r.his:

"First, the working color is set to the color of the corresponding cell in the map layer chosen to represent hue. Second, this color is multiplied by the red intensity of that cell in the intensity map layer. This map layer should have an appropriate gray-scale color table associated with it."

r.in.gdal reads GRASS color rules embedded in metadata, e.g. for raster_4class.tif:

  Metadata:
    COLOR_TABLE_RULES_COUNT=5
    COLOR_TABLE_RULE_RGB_0=1.000000e+00 1.600000e+00 255 255 0 0 255 0
    COLOR_TABLE_RULE_RGB_1=1.600000e+00 2.200000e+00 0 255 0 0 255 255
    COLOR_TABLE_RULE_RGB_2=2.200000e+00 2.800000e+00 0 255 255 0 0 255
    COLOR_TABLE_RULE_RGB_3=2.800000e+00 3.400000e+00 0 0 255 255 0 255
    COLOR_TABLE_RULE_RGB_4=3.400000e+00 4.000000e+00 255 0 255 255 0 0

r.external does not read these metadata and instead assigns some default color rules, or none, depending on the input datatype.

In short, you need to assign appropriate gray-scale color rules first for r.his.