When PIXELTYPE=SIGNEDBYTE, band.GetStatistics returns wrong values

[jjr8]

Evan R. requested that I open a ticket on this issue that I raised on gdal-dev.

With GDAL 1.6.2, I can successfully create a signed 8-bit raster with the HFA or GTiff drivers using the PIXELTYPE=SIGNEDBYTE option:

>>> a = numpy.array([[-128, -1, 0, 1, 127]], dtype='int8')
>>> a
array([[-128,   -1,    0,    1,  127]], dtype=int8)
>>> b = numpy.cast['uint8'](a)
>>> b
array([[128,   255,    0,    1,  127]], dtype=uint8)
>>>
>>> gdal.UseExceptions()
>>> driver = gdal.GetDriverByName('HFA')
>>> dataset = driver.Create(r'C:\Temp\test_int8.img', 5, 1, 1, gdal.GDT_Byte, ['PIXELTYPE=SIGNEDBYTE'])
>>> band = dataset.GetRasterBand(1)
>>> band.WriteArray(b)
0
>>> del band, dataset, driver

When I look at that raster with ArcGIS, it shows it is signed 8-bit integer and the five values show up properly. But when I open it with GDAL and call band.GetStatistics(), it looks like the calculation is performed using unsigned 8-bit integers:

>>> dataset = gdal.Open(r'C:\Temp\test_int8.img', gdalconst.GA_ReadOnly)
>>> band = dataset.GetRasterBand(1)
>>> band.GetMetadata('IMAGE_STRUCTURE')
{'PIXELTYPE': 'SIGNEDBYTE'}
>>> band.GetStatistics(False, True)
[0.0, 255.0, 102.2, 95.199579831005551]

Evan commented that "this should be fixable within GDAL itself by testing the presence of PIXELTYPE=SIGNEDBYTE and by casting the value to a signed byte in that case."

Thank you very much for looking at this. Please adjust Trac fields appropriately. It would be nice if this was fixed very soon, so I set the milestone to 1.6.3. But if you must fix it later I can implement a workaround of having my application calculate the statistics itself for SIGNEDBYTE rasters.

One final note: band.SetNoDataValue() and band.GetNoDataValue() appear to work correctly:

>>> driver = gdal.GetDriverByName('HFA')
>>> dataset = driver.Create(r'C:\Temp\test_int8d.img', 5, 1, 1, gdal.GDT_Byte, ['PIXELTYPE=SIGNEDBYTE'])
>>> band = dataset.GetRasterBand(1)
>>> band.SetNoDataValue(-20)
0
>>> band.WriteArray(b)
0
>>> del band, dataset, driver
>>> dataset = gdal.Open(r'C:\Temp\test_int8d.img', gdalconst.GA_ReadOnly)
>>> band = dataset.GetRasterBand(1)
>>> band.GetMetadata('IMAGE_STRUCTURE')
{'PIXELTYPE': 'SIGNEDBYTE'}
>>> band.GetNoDataValue()
-20.0

I tried several different values (not just -20) and also verified that ArcGIS reported the same NoData value and that it showed NoData pixels as NoData.

[jjr8]

Even - sorry for misspelling your name in this report!

[jjr8]

Here is another problem which could be related. Please test this when you fix the original problem. Please let me know if you'd like me to open a new ticket on this issue.

The problem is that band.GetStatistics() throws an exception if all of the following are true:

• HFA driver
• PIXELTYPE=SIGNEDBYTE
• COMPRESSED=YES
• band.FlushCache() is called prior to band.GetStatistics()

Here is an example:

>>> a = numpy.array([[-128, 5, 127]]*3, dtype='int8')
>>> a
array([[-128,    5,  127],
       [-128,    5,  127],
       [-128,    5,  127]], dtype=int8)
>>> b = numpy.cast['uint8'](a)
>>> b
array([[128,   5, 127],
       [128,   5, 127],
       [128,   5, 127]], dtype=uint8)
>>>
>>> gdal.ErrorReset()
>>> gdal.UseExceptions()
>>> driver = gdal.GetDriverByName('HFA')
>>> dataset = driver.Create(r'C:\Temp\test_int8d.img', 5, 3, 3, gdal.GDT_Byte, ['COMPRESSED=YES', 'PIXELTYPE=SIGNEDBYTE'])
>>> band = dataset.GetRasterBand(1)
>>> band.WriteArray(b)
0
>>> band.FlushCache()
>>> band.GetStatistics(False, True)
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "C:\Python25\Lib\site-packages\GeoEco\AssimilatedModules\osgeo\gdal.py", line 725, in GetStatistics
    return _gdal.Band_GetStatistics(*args)
RuntimeError: Failed to compute statistics, no valid pixels found in sampling.

If any one of the four conditions noted above is not true, !band.GetStatistics?() does not throw the exception. (Of course, for PIXELTYPE=SIGNEDBYTE, it still returns the wrong values.)

[jjr8]

Replying to jjr8:

Ack! I made a copy and paste mistake. The problem does exist, but I miscopied one line of code in the example. This:

>>> dataset = driver.Create(r'C:\Temp\test_int8d.img', 5, 3, 3, gdal.GDT_Byte, ['COMPRESSED=YES', 'PIXELTYPE=SIGNEDBYTE'])

should be this:

dataset = driver.Create(r'C:\Temp\test_int8d.img', 3, 3, 1, gdal.GDT_Byte, ['COMPRESSED=YES', 'PIXELTYPE=SIGNEDBYTE'])

The 5, 3, 3 should have been 3, 3, 1. Sorry for any confusion.

[rouault]

First issue fixed in trunk (r17676) and branches/1.6 (r17677). Test added in r17678

The second issue you raise has nothing to do with GetStatistics?() but is in the HFA driver itself. I've filed it as #3152