r.out.gdal tiff output does not work

(I see the same with 6.2.1's old script version of r.out.gdal)

it works if the values are in the range of 0-255:

  #spearfish
  r.out.gdal in=fields out=fields.tif type=Byte
  qiv fields.tif   # qiv is a quick image viewer

but not if all values are out of range:

  r.info -r elevation.dem
    min=1066
    max=1840
  r.info -t elevation.dem
    datatype=CELL
  r.out.gdal in=elevation.dem out=elv.tif type=Byte

  gdalinfo -stats elv.tif | grep STATISTIC
    STATISTICS_MINIMUM=255
    STATISTICS_MAXIMUM=255
    STATISTICS_MEAN=255
    STATISTICS_STDDEV=0

  r.in.gdal in=elv.tif out=elv2
  r.info -r elv2
    min=255
    max=255

we can try with half the raster in range:

  r.mapcalc "elev_0 = elevation.dem - 1066"
  r.info -r elev_0
    min=0
    max=774
  r.out.gdal in=elev_0 out=elv0.tif type=Byte
  gdalinfo -stats elv0.tif | grep STATISTIC
    STATISTICS_MINIMUM=0
    STATISTICS_MAXIMUM=254
    STATISTICS_MEAN=138.68943338231
    STATISTICS_STDDEV=56.526296226099

If I load elv0.tif into an image viewer it comes out as all black. But in Gimp if you use the eyedropper you can see the index number changes, but the palette is black for all 256 index entries.

gdalinfo confirms palette info was not written to file:

Band 1 Block=1899x4 Type=Byte, ColorInterp=Palette
  NoData Value=255
  Metadata:
    COLOR_TABLE_RULES_COUNT=0
  Color Table (RGB with 256 entries)
    0: 0,0,0,255
    1: 0,0,0,255
    2: 0,0,0,255
    3: 0,0,0,255
    4: 0,0,0,255
...
  255: 0,0,0,255

I do not know how to apply the Gimp's Topographic.gpl color palette to the image's private "Colormap" index. It must be possible somehow, but I don't know Gimp well enough..

If you want to store integer values >255 in the Tiff, you will have to use one of Int16,UInt16,UInt32,Int32[,CInt16,CInt32?].

If you use type=Int32 etc you should preserve the CELL data for use in other GIS software, even if the saved palette is bogus.

In the above elev_0 case you can see that the 0-774 raster exported as Byte (8bit) only preserves elevation values in the range of 0-254.

Only Byte and UInt16 will attempt to write out a color palette to the GeoTIFF AFAIK.

Hamish

The data range as reported in ArcMap? 9.2 when opening a geotiff exported from r.out.gdal is (1.17549e-038,3.40282e+038).

I have noticed that when setting the files symbology, when set to a fixed number of classes (ranges) the correct data range is displayed, and the image seems to work fine. This suggests that it may be an out of range or unusual (in terms of ArcMap? logic) no-data value.

I have updated my notes here

Replying to dylan:

The data range as reported in ArcMap? 9.2 when opening a geotiff exported from r.out.gdal is (1.17549e-038,3.40282e+038).

If it were the NULL value, I'd expect to see extreme values only on one side, not both. Endian problems is another suspect, but probably unlikely if all on x86.

I have noticed that when setting the files symbology, when set to a fixed number of classes (ranges) the correct data range is displayed, and the image seems to work fine. This suggests that it may be an out of range or unusual (in terms of ArcMap? logic) no-data value.

can you do like d.historgram or r.univar in Arc to see if the values are streched over that range or are mostly correct but with weird outliers?

what if you try other output types besides Float32?

Hamish

(see also trac bug #67)

is there some other free software, either with or without GDAL that could be used as a 3rd party test point?

Any info on what Arc treats as the no-data value? Does it respect what GDAL puts in the metadata? Does Arc have the equivalent of r.null to clean it up?

Hamish

Eric mentioned several viewers that don't like the tiff produced by r.out.gdal and for me it won't work with Preview on Mac and gimp. You can try for yourself: here is an original tiff file that is readable by gimp and preview and imports to GRASS nc_spm_07 LOCATION using r.in.gdal (change tfwmf to tfw). http://skagit.meas.ncsu.edu/~helena/grasswork/grassbookdat07/ncexternal/IMG3720079200P20040216.tfwmf http://skagit.meas.ncsu.edu/~helena/grasswork/grassbookdat07/ncexternal/IMG3720079200P20040216.tif

After exporting with r.out.gdal the image is not readable. But you are right that the numbers appear to be correct but there seem to be some out of range values (nulls?). With some other rasters where I had nulls I tried to define nulls as -9999 and export to AAGrid using r.out.gdal but there was a problem too. Here is one example of exported tiff http://skagit.meas.ncsu.edu/~helena/grasswork/ortho_2001_t792_1m.tif http://skagit.meas.ncsu.edu/~helena/grasswork/ortho_2001_t792_1m.tfw

there may be some option that exports the tiffs correctly but the only one that was readable was this http://skagit.meas.ncsu.edu/~helena/grasswork/lakes.tif http://skagit.meas.ncsu.edu/~helena/grasswork/lakes.tfw

Thanks for looking into it

Helena

Helena,

I'm coming into this rather late, and just skimming things. I looked at your ortho_2001_t7921_1m.tif file (no smaller files demostrating the problem!) and I see it is 32bit integer pixel values. I'm not surprised lots of packages can't read this.

I looked at lakes.tif and it seems very similar to the ortho image in configuration - also one 32bit integer sample. So I'm surprised that it is readable and the ortho image is not in some (unnamed?) application.

PS. it would be convenient to post tiffinfo reports for some of the problem files to learn more without having to download them. For instance the ortho report is:

TIFF Directory at offset 0x23768d2 (37185746)
  Image Width: 3048 Image Length: 3048
  Bits/Sample: 32
  Sample Format: signed integer
  Compression Scheme: None
  Photometric Interpretation: min-is-black
  Samples/Pixel: 1
  Rows/Strip: 1
  Planar Configuration: single image plane
  Tag 33550: 1.000000,1.000000,0.000000
  Tag 33922: 0.000000,0.000000,0.000000,637033.000000,222504.000000,0.000000
  Tag 34735: 1,1,0,16,1024,0,1,1,1025,0,1,1,1026,34737,24,0,2048,0,1,4269,2049,34737,6,24,2054,0,1,9102,3072,0,1,32767,3074,0,1,32767,3075,0,1,8,3076,0,1,9001,3078,34736,1,2,3079,34736,1,3,3084,34736,1,1,3085,34736,1,0,3086,34736,1,4,3087,34736,1,5
  Tag 34736: 33.750000,-79.000000,36.166667,34.333333,609601.220000,0.000000
  Tag 34737: Lambert Conformal Conic|grs80|
  Tag 42112: <GDALMetadata>
  <Item name="COLOR_TABLE_RULES_COUNT" sample="0">1</Item>
  <Item name="COLOR_TABLE_RULE_RGB_0" sample="0">0.000000e+00 2.550000e+02 0 0 0 255 255 255</Item>
</GDALMetadata>

We've been having similar problems here and I've tested a bit further with current svn.

At http://geog-pc40.ulb.ac.be/grass/misc/r_out_gdal/ you can find relevant files:

- test_loc.tgz: test location with one raster file - test.tif: this same raster file exported with 'r.out.gdal type=UInt16' - test_tiffinfo.txt: tiffinfo output on test.tif (note the errors at the beginning) - test _gdalinfo.txt gdalinfo output on test.tif, which shows the correct color categories and then an all-black color table from 0 - 65535. Is it normal to have both 'COLOR_TABLE_RULE_RGB_0' entries and a color table in the tiff ?

The tif can be read and displayed by QGIS. SAGA gives a weird image (both with band and pixel interleave): test_saga.png. gvSIG crashes when trying to load the file. These are not some general picture viewer, but cartography/GIS programs, so I would hope that we should be able to solve this for reasons of interoperability.

Final test which seems to show that there is something wrong: r.in.gdal in=test.tif out=test: import takes a very long time and creates a huge colr table file:

% 0 65535 0:0 1:0 2:0 3:0 4:0 5:0 6:0 7:0 8:0 9:0 10:0 11:0 ... 65535:0 65535:0

Whereas the original colr file only contains this: % 334 460 334:0:191:191 359.2:0:255:0 359.2:0:255:0 384.4:255:255:0 384.4:255:255:0 409.6:255:127:0 409.6:255:127:0 434.8:191:127:63 434.8:191:127:63 460:20

The display of test also takes a long time, I suppose because of this long colr table.

I'm moving this up to critical as I think this is a fairly basic feature which should "just work"...

Moritz

Replying to dylan:

The data range as reported in ArcMap? 9.2 when opening a geotiff exported from r.out.gdal is (1.17549e-038,3.40282e+038).

This perfectly right insofar as this is the potential data range for float32. ArcMap? does not check the actual data range when loading a new layer the first time, even when real minimum and maximum values are given in the metadata. The display of a layer gets more reasonable when choosing e.g. histogram equalization as stretch because only now ArcMap? is looking at the real data range, gets real min and max values and adjusts the data visualization accordingly.

ArcMap? behaves like this too when opening ArcInfo? ASCII grid and ESRI .hdr labeled files with float32 data (same elevation raster exported into different file formats, keeping data type float32). Maybe this is not a GeoTIFF/GRASS/GDAL bug, but ArcMap? is simply not very user-friendly with data visualization?

Here, GRASS is much more user-friendly, because it always knows the data range of a raster and, if no color table has been specified, automatically generates a reasonable color table, usually in color and not grayscale as ArcMap? does by default...

The tiff file is handled by ArcMap? all right, data and georeferencing information are there, but data visualization needs to be modified by hand.

I have noticed that when setting the files symbology, when set to a fixed number of classes (ranges) the correct data range is displayed, and the image seems to work fine. This suggests that it may be an out of range or unusual (in terms of ArcMap? logic) no-data value. I have updated my notes here

Replying to hamish:

(I see the same with 6.2.1's old script version of r.out.gdal) it works if the values are in the range of 0-255: but not if all values are out of range:

It is equally right to say that the chosen datatype=Byte is not capable of holding the full data range and a different datatype must be chosen.

{{{ r.info -r elevation.dem min=1066 max=1840 r.info -t elevation.dem datatype=CELL r.out.gdal in=elevation.dem out=elv.tif type=Byte r.in.gdal in=elv.tif out=elv2 r.info -r elv2 min=255 max=255 }}} we can try with half the raster in range: {{{ r.mapcalc "elev_0 = elevation.dem - 1066" r.info -r elev_0 min=0 max=774 r.out.gdal in=elev_0 out=elv0.tif type=Byte

No color table was specified with r.colors, so what data visualization scheme is to be expected?

If I load elv0.tif into an image viewer it comes out as all black. But in Gimp if you use the eyedropper you can see the index number changes, but the palette is black for all 256 index entries.

Gimp and many other image viewers only support 8 bit per channel, cinepaint e.g. supports 32 bit. Most of them can not display tiff files with anything else but Byte: only tiff images, not tiff datasets. Elevation is a dataset, not an image. INT16 or UINT16 would be appropriate to hold elevation data, in case of cm accuracy with meters as units, FLOAT32 would be ok. Such tiff files hold a dataset. For elevation data, GRASS offers several color tables which can be written to a tiff file as a suggested data visualization scheme, if they have been set with r.colors. Still, datatypes other than Byte can only be read with spatial data viewers that could make use of this colour table, but not all of them are actually able to do so. They could instead choose their own color schemes.

gdalinfo confirms palette info was not written to file:

But the color palette was indeed written to the file.

{{{ Band 1 Block=1899x4 Type=Byte, ColorInterp?=Palette NoData? Value=255 Metadata: COLOR_TABLE_RULES_COUNT=0

Below here. The instructions read to display each value with R:G:B = 0:0:0 = black. A data visualization scheme should have been provided with r.colors to have control over the colortable written to the tiff file.

Color Table (RGB with 256 entries) 0: 0,0,0,255 1: 0,0,0,255 2: 0,0,0,255 3: 0,0,0,255 4: 0,0,0,255 ... 255: 0,0,0,255 }}}

If you want to store integer values >255 in the Tiff, you will have to use one of Int16,UInt16,UInt32,Int32[,CInt16,CInt32?].

It would help a lot to have a list with these data types and their ranges in the man page so that the appropriate data type for file export can be chosen.

If you use type=Int32 etc you should preserve the CELL data for use in other GIS software, even if the saved palette is bogus.

Maybe the color table is not so important, different people have different preferences e.g. for elevation data visualization. If imagery is exported, e.g. after atmospheric correction or masking, the datatype should probably be Byte as it was originally, e.g. to make postprocessing in Gimp possible. A tiff file that holds an image (e.g. GRASS display saved as tiff, a scanned image, aerial photography, LANDSAT) does not have a color table, usually. It does not need one, because the data in the file are already R:G:B or grayscale (better: single band) values.

As long as both the data in the tiff file and its georeferencing are identical to the original raster, the tiff file could be regarded as perfectly ok. Suggested data visualization with an embedded color table is nice to have, but coloring schemes are a matter of taste and not crucial for further processing (e.g. terrain analysis with a DEM).

There is still the issue of nodata: a nodata value must fall into the range of the chosen datatype, otherwise this information gets lost. r.out.gdal does not check if a supplied nodata value is valid.

Only Byte and UInt16 will attempt to write out a color palette to the GeoTIFF AFAIK.

For UINT16, the color table in the tiff file will always have 65535 entries, gdal does that. This slows down display. Have an option in r.out.gdal to choose whether a color table should be written to the tiff file?

Markus

Replying to mlennert:

We've been having similar problems here and I've tested a bit further with current svn.

I repeated the tests and it looks better now, especially if the colortable is not written to the tiff file.

Tested on Linux 64bit with gdal-1.5.2 compiled from source grass-6.4.svn_src_snapshot_2008_07_26 compiled from source QGIS-0.9.1 compiled from source (I don't get versions 0.10.0 and 0.11.0 to display ECW and JP2ECW files, only 0.9.1 does it) gvSIG-1.1.1 precompiled 32bit binary with built-in Java on Linux 64bit SAGA-2.0.3 on Windows, I don't get it to run on Linux 64bit

r.out.gdal modified: - if a nodata value is given, it is adjusted to the range of the selected datatype if it is out of range - basic statistics (min, max, mean, stddev) are calculated for the current region, replacing NULL cells with nodata value (default or specified) and then exported with GDALSetRasterStatistics() - colortable is only exported on request (new flag needs to be set) and only if there is an entry in directories colr or colr2

Additionally, a warning would be useful if the selected datatype does not cover the range of raster values to be exported (current region settings used by r.out.gdal, not full raster range as in cell_misc).

I used test_loc.tgz and exported mnt_test@PERMANENT twice as UINT16, once with colortable, once without colortable. Additionally, I used the file test.tif from the link below.

QGIS-0.9.1 loads and displays everything. test.tif appears all black, but the data are there and ok. My export with colortable is displayed exactly like in GRASS, the export without colortable in grayscale. Again, the data are ok according to stats calculated by QGIS.

SAGA-2.0.3 displays all three as grayscale, i.e. disregards any colortables. The data are ok.

gvSIG needs a proper projection specified, otherwise raster files are not loaded (user error), so I used a projection similar to the test_loc location, EPSG: 32766. I could not find the projection of test_loc in gvSIG. All three files are loaded, no crash. Only the file without colortable is displayed in grayscale, the two other files with colortable are all black. I don't know how to override the default display and assign a custom color table in gvSIG like r.colors does in GRASS. The data are ok.

Some other, commercial, software doesn't like GeoTIFFs with colortables either, so maybe it is a good idea not to export colortables, or only on request with a warning? One could also ask developers of other software if colortables in GeoTIFF are supported and read at all. Of the three packages tested here, only QGIS makes proper use of the colortable, and QGIS is too close to GRASS to be regarded as independent testing suite.

At http://geog-pc40.ulb.ac.be/grass/misc/r_out_gdal/ you can find relevant files: - test_loc.tgz: test location with one raster file - test.tif: this same raster file exported with 'r.out.gdal type=UInt16' - test_tiffinfo.txt: tiffinfo output on test.tif (note the errors at the beginning)

tiffinfo doesn't read GeoTIFF tags with projection info, error can be ignored

- test _gdalinfo.txt gdalinfo output on test.tif, which shows the correct color categories and then an all-black color table from 0 - 65535. Is it normal to have both 'COLOR_TABLE_RULE_RGB_0' entries and a color table in the tiff ?

Yes, the COLOR_TABLE_RULE entries are written by r.out.gdal in the presence of a colortable and I guess there for convenience so that GRASS users can reconstruct color rules instead of using a very long colortable. Apparently COLOR_TABLE_RULE entries are custom metadata not read by other software. Since there are colortable rules for the GRASS raster, there should be entries in the colortable in the test.tif file that are not all black. Maybe something was wrong with the GRASS svn version used by you? Can you test again?

Markus

r.out.gdal GeoTIFF compatibility testing and new r.out.gdal test version

Changes in r.out.gdal test version new flag for colortable export, default to no nodata value adjustment to range of selected data type with warning message, if user-specified nodata value is out of range data range checking against the selected data type, abort with error if range of selected data type exceeded alignment of current region extends to raster resolution (original uses current region extends and resolution), current region restored on exit raster band statistics (min, max, mean, stddev) written to metadata if supported, GDAL decides

I tested with QGIS-0.9.1 on Linux 64bit, QGIS -0.11.0-1 on XP, SAGA-2.0.3 on XP (doesn't work for me on Linux 64bit), gvSIG-1.1.2 (32bit binary) on Linux 64bit, gvSIG-1.1.1 on XP, PCI Geomatica FreeView? 9.1 on XP, ER Viewer 7.2 (ER Mapper) on XP, and ArcMap? 9.2 (ESRI) on XP.

Tests were performed with the elevation raster map layer (range within Byte, FCELL map) in the North Carolina sample dataset, creating a MASK for elevation values <= 70m, needed for testing of nodata handling with nodata=0. The elevation raster was exported as Byte, UInt16, Int32, and Float32, once with colortable (Byte and UInt16 only, colortable export for other data types not supported by gdal), once without colortable, always as GeoTIFF.

QGIS-0.9.1 on Linux 64bit

display: all ok

colortable: all ok

cell values: all ok

nodata: all ok

QGIS-0.11.1 on XP

display: all ok

colortable: all ok

cell values: all ok

nodata: all ok

gvSIG-1.1.2 (Linux 32bit binary) on Linux 64bit. Projection EPSG:32119

display: fail for UInt16 with colortable, all other ok

colortable: fail for UInt16, Byte ok

cell values: fail for Byte!!! (sometimes reports negative values), all other ok. Interestingly the display seems to use correct values. Bug in gvSIG-1.1.2

nodata: nodata not supported, actual cell value reported

gvSIG-1.1.1 on XP, Projection EPSG:32119

display: fail for UInt16 with colortable, all other ok

colortable: fail for Uint16, ok for Byte

cell values: fail for Byte!!! (sometimes reports negative values), all other ok. Interestingly the display seems to use correct values. Bug in gvSIG-1.1.1

nodata: nodata not supported, actual cell value reported

PCI Geomatica FreeView? 9.1 on XP

display: all ok

colortable: all ok

cell values: all ok

nodata: nodata not supported, actual cell value reported

ER Viewer 7.2 (ER Mapper) on XP

display: fail for UInt16 with colortable, all other ok

colortable: fail for UInt16, ok for Byte

cell values: cell value query not supported

nodata: cell value query not supported

SAGA-2.0.3 on XP (I don't get it to run on Linux 64bit)

display: all ok

colortable: all ignored

cell values: all ok

nodata: all ok

ESRI ArcMap? 9.2 on XP

display: all ok, but all without colortable and other than Byte need manual adjustment which is painfully slow, even for these small raster layers (Properties -> Symbology -> Stretch -> Minimum-Maximum)

colortable: all ok

cell values: all ok

nodata: all ok

Colortables are not allways properly rendered. GIMP and image viewers which are not using GDAL render a colortable for Byte type properly. Data types other than Byte can not be read by these 'standard' image editing and viewing applications.

If nodata specification is not supported, applications return the actual cell value, in this case 0, instead of "NoData?", "NULL", "" or similar. The files are still perfectly readable, only nodata assignment doesn't happen.

I suggest adding a flag for colortable export with default to no and nodata value checking as well as data range checking against the selected data type.

Another issue I discovered is that r.out.gdal uses the current region settings for export. This can cause implicit resampling to a different resolution which should be avoided IMHO. Instead, r.out.gdal could use the current region extends, but align the region settings used for export to the raster to be exported. This way, a subregion can still be exported, but the resolution of the GRASS raster is preserved. Also maybe check whether the current region extends are within the raster extends? Alignment of the current region to the raster to be exported by a new r.out.gdal is not crucial, it can be done manually beforehand with g.region align="myraster", but it would be very convenient if r.out.gdal does that and restores the original region settings after successful export.

The test version of gdal is available at http://markus.metz.giswork.googlepages.com/r.out.gdal.conservative.tar.gz

See enclosed README for full description and justification of all changes.

Markus

Great work! Would you mind to also test back-import into GRASS?

Markus

Replying to neteler:

Great work!

Thanks!

Would you mind to also test back-import into GRASS?

All GeoTIFFs previously exported with the test version of r.out.gdal were back-imported to GRASS-6.4.svn (snapshot_2008_08_23)

GeoTIFF files:
elevation.Byte.clrtbl.tif
elevation.Byte.noclrtbl.tif
elevation.UInt16.clrtbl.tif
elevation.UInt16.noclrtbl.tif
elevation.Int32.noclrtbl.tif
elevation.Float32.noclrtbl.tif

display: all ok, very slow for the UInt16 back-import with colortable
colortable: all ok
cell values: all ok
nodata: all ok

Back-importing elevation.UInt16.clrtbl.tif is pretty slow. Display takes quite some time too. It is a colortable with 65535 entries, one for each potential cell value, an exact copy of the GeoTIFF colortable. Removing the colortable with r.colors -r, then copying the colortable of the back-import of elevation.Byte.clrtbl.tif results in the same display, but now fast on refresh, as expected.

Back-import of all other GeoTIFF files is fast, as is display in GRASS.

I would prefer to use r.out.gmt to get the GRASS color rules into an extra file (GMT .cpt). This would also be a workaround for export file types other than Byte and UInt16 where GDAL does not support writing out a colortable. But then I don't know how widely GMT .cpt files are supported.

Markus

Markus (Metz),

what about integrating your fixes from

http://markus.metz.giswork.googlepages.com/r.out.gdal.conservative.tar.gz

?

Markus

[keep it in the bug report please!]

Comment (by neteler): what about integrating your fixes from

http://markus.metz.giswork.googlepages.com/r.out.gdal.conservative.tar.gz

?

Markus Metz wrote:

Sure, no objections from my side, I'm using this version only. But r.out.gdal is a very important module of GRASS and maybe some more testing is required. Also, the new features may not find approval by all,

yep

e.g. I've put in rather restrictive NULL cell handling: the user has to specify a nodata value that falls within the range of the selected datatype (Byte, UInt16, Int16, ...), a nodata value is no longer chosen automatically as in the original version. If a nodata value is not given, but NULL cells are present, r.out.gdal aborts with an error message.

When possible, I think it better to have a sensible default then give the user the option to override. That's what we have in the existing code. NULL=max possible is a common usage, so not bad for a default. Your change is too conservative for my vote, and generally I don't like banning users from doing things just because my imagination is limited as to why they would want to do something funny. If they really want to set nodata= out of bounds, and explicitly give a nodata= value specifying that, I'd say let them (with a warning in case they didn't realize, of course).

My version also no longer uses the current region resolution, instead the current region extends are aligned to the resolution of the raster to be exported to avoid any implicit resampling.

Nope. I see your reasons, but this is in violation of GRASS's standard raster semantics, and one of the reasons for the r.out.gdal.sh -> r.out.gdal (.c) was to "fix" this. (sort of)

Maybe it makes the GRASS<->R stats people happier?, but shouldn't be the default method, and I doubt offered at all. I would think this better handled with a note in the help page with a suggestion to use "g.region res=.. -a" if they want that.

And the colortable is only exported on request, and then with a warning message.

IMO that should be automatic with Byte maps (255 levels), perhaps needed with UInt16 maps (with 65535 levels) it should be optional...

If these changes are ok with you then integrate the changes, otherwise maybe only some, but not all changes could be integrated. More discussion needed on how to change/improve r.out.gdal?

as above. I need to reread the (long) bug history now, but fwiw I just fixed yesterday a bug in SVN where nodata= value given by the user was only used if type= was also given, otherwise the auto-type determination reset it.

I also noticed that QGIS 0.8.1 (old, but latest debianGIS package for Etch/stable) has a bug where the first palette entry always is black, even when the GeoTIFF says otherwise. (loads/views fine in other software) In my case it was a 7-color palette with cat 0 being white in the color table.

Hamish

Replying to mmetz (and others):

No color table was specified with r.colors, so what data visualization scheme is to be expected?

if none is specified (no colr/ files exists for the raster) the default "rainbow" is used. GRASS's libgis does this automatically, AFAIK the module never knows that r.colors hasn't been run.

Gimp and many other image viewers only support 8 bit per channel, cinepaint e.g. supports 32 bit.

this is a CELL map exported with type=Byte.

Most of them can not display tiff files with anything else but Byte: only tiff images, not tiff datasets.

the spearfish example I gave was Byte-only.

Elevation is a dataset, not an image.

True, but if elev data is put in the correct range, Gimp neither knows nor cares what the levels represent. I was just using it as random values, not anything meaningful.

Spearfish's elevation.dted would be a better sample dataset to use here, it is CELL 0-255. or convert spot.image into RGB bands. Other test datasets are the spearfish imagery60 data, and the North Carolina dataset lsat7 images.

Still, datatypes other than Byte can only be read with spatial data viewers that could make use of this colour table, but not all of them are actually able to do so.

Important to know when this is happening, not really our problem, and not much we can do about it as long as we are writing metadata to spec.

Hamish:

gdalinfo confirms palette info was not written to file: ... COLOR_TABLE_RULES_COUNT=0 ... 0: 0,0,0,255 ... 255: 0,0,0,255

MM:

But the color palette was indeed written to the file.

An empty one was written while within GRASS it had real colors. For the expected colors run:

seq 0 255 | r.what.color -i in=elev_0

0: 255:255:0
1: 254:255:0
2: 252:255:0
...
253: 0:255:161
254: 0:255:163
255: 0:255:165

It would help a lot to have a list with these data types and their ranges in the man page so that the appropriate data type for file export can be chosen.

Yes. Can someone in the know supply some suggested text?

Maybe the color table is not so important, different people have different preferences

The color table will be carefully constructed/chosen by some people and they may wish to keep it. e.g. r.colors.stddev addon script, r.cpt2grass addon script, and NASA color tables given for MODIS sat data (see grass's MODIS wiki page). These data are real (calibrated into sea surface temperature or mg/mL chlorophyll) and also have specially designed color tables (incl log-scale) which are nice to keep.

Suggested data visualization with an embedded color table is nice to have, but coloring schemes are a matter of taste

geotiff has a method for assigning color palettes, all I suggest is populating that correctly, apparently we are currently (or were) failing at.

There is still the issue of nodata: a nodata value must fall into the range of the chosen datatype, otherwise this information gets lost. r.out.gdal does not check if a supplied nodata value is valid.

Hopefully better in SVN now. If the user manually specifies an out of range value, I don't mind respecting that.

For UINT16, the color table in the tiff file will always have 65535 entries, gdal does that.

perhaps this should be addressed in GDAL then; only write palette entries which are given.

Have an option in r.out.gdal to choose whether a color table should be written to the tiff file?

For Uint16, how to just get the "COLOR_TABLE_RULE_RGB_0" part of this without the "(RGB with 66535 entries)" part?

r.out.gdal geotiff from 0-255 (black & white) ortho photo with grass's sepia color table:

  Metadata:
    COLOR_TABLE_RULES_COUNT=10
    COLOR_TABLE_RULE_RGB_0=0.000000e+00 2.550000e+01 0 0 0 34 23 9
    COLOR_TABLE_RULE_RGB_1=2.550000e+01 5.100000e+01 34 23 9 68 47 19
    COLOR_TABLE_RULE_RGB_2=5.100000e+01 7.650000e+01 68 47 19 101 72 33
    COLOR_TABLE_RULE_RGB_3=7.650000e+01 1.020000e+02 101 72 33 130 96 49
    COLOR_TABLE_RULE_RGB_4=1.020000e+02 1.275000e+02 130 96 49 158 122 71
    COLOR_TABLE_RULE_RGB_5=1.275000e+02 1.530000e+02 158 122 71 183 147 96
    COLOR_TABLE_RULE_RGB_6=1.530000e+02 1.785000e+02 183 147 96 206 174 126
    COLOR_TABLE_RULE_RGB_7=1.785000e+02 2.040000e+02 206 174 126 225 200 162
    COLOR_TABLE_RULE_RGB_8=2.040000e+02 2.295000e+02 225 200 162 242 227 204
    COLOR_TABLE_RULE_RGB_9=2.295000e+02 2.550000e+02 242 227 204 255 254 251
  Color Table (RGB with 256 entries)
    0: 0,0,0,255
    1: 1,0,0,255
    2: 2,1,0,255
    3: 4,2,1,255
    4: 5,3,1,255
    5: 6,4,1,255
 ...

r.out.gdal modified: - if a nodata value is given, it is adjusted to the range of the selected datatype if it is out of range

It is better to respect the user's instructions. Give them a warning if you like, but do what they say; don't quietly "fix" it for them.

basic statistics (min, max, mean, stddev) are calculated for the current region, replacing NULL cells with nodata value (default or specified) and then exported with GDALSetRasterStatistics()

how is this better than "gdalinfo -stats" ? I don't understand, you convert NULL to nodata value in stats? So for a CELL map with lots of NULL areas, if I make nodata=100 (well above the number of used palette entries) or leave it at 255, the stats will be skewed towards 100 or 255 instead of just considering the real data? huh? "gdalinfo -stats" seems to do the right thing here.

colortable is only exported on request (new flag needs to be set) and only if there is an entry in directories colr or colr2

This is wrong, libgis is in charge of handling colors. Use the API not raw file access. For one thing, otherwise you miss out on the default mapping if there isn't a colr/ file.

Additionally, a warning would be useful if the selected datatype does not cover the range of raster values to be exported

i.e. if the map was CELL 0-744 as in my elev_0 example earlier in this ticket, but the user asked for Byte (0-255). I think that's a very good idea.

Start by moving range maxes into header #defines for each type from the currently scattered "nodataval = -1E37f;" etc lines in the r.out.gdal source code.

Some other, commercial, software doesn't like GeoTIFFs with colortables either, so maybe it is a good idea not to export colortables, or only on request with a warning?

If other software ignore valid GeoTIFF metadata tags, then that's their problem and their responsibility to provide an alternative for their users. All we can do is correctly implement the spec. If other software read but incorrectly handle valid GeoTIFF metadata, then that's really their bug to deal with.

Of the three packages tested here, only QGIS makes proper use of the colortable, and QGIS is too close to GRASS to be regarded as independent testing suite.

The test suite you can't argue with is the data in the file. It either follows the geotiff spec or does not. (no idea how free-form or tight that is in this regard)

Apparently COLOR_TABLE_RULE entries are custom metadata not read by other software.

Any confirmation on this? It is custom to GRASS or custom to OSGeo family (GDAL dependents) or suggested by the geotiff spec or ...???

Even if 66535 rules are turned off by default by r.out.gdal for UInt16, these range rules would be nice to keep as they are rarely very long.

current region restored on exit

module should not have to do that, it should set WIND_OVERRIDE variable or so, so if module is interrupted by ctrl-C or whatever the temporary change does not become permanent. (but as mentioned earlier, this violates grass's raster module conventions and the feature should not be applied in SVN)

raster band statistics (min, max, mean, stddev) written to metadata

perhaps nice for casual browsing, but again what's wrong with "gdalinfo -stats" or similar scanning for that task?

if supported, GDAL decides

what does that mean? when is min/max etc not supported??

Tests were performed with the elevation raster map layer (range within Byte, FCELL map) in the North Carolina sample dataset, creating a MASK for elevation values <= 70m, needed for testing of nodata handling with nodata=0. The elevation raster was exported as Byte, UInt16, Int32, and Float32, once with colortable (Byte and UInt16 only, colortable export for other data types not supported by gdal), once without colortable, always as GeoTIFF.

Please supply map names and exact commands so we can replicate.

QGIS-0.9.1 on Linux 64bit QGIS-0.11.1 on XP ... ESRI ArcMap? 9.2 on XP

this is handy information for users to be aware of, but (again) not much we can do about 3rd party support beyond coding to geotiff spec/recomendations.

It is nice to see QGIS do well. If this is filtered through a GeoTIFF and not using the GRASS native format plugin, then it is a good sign that we are all coding to the same tune. I don't think the code itself has much cut & paste between the projects, and GDAL I/O is independently well-tested.

If nodata specification is not supported, applications return the actual cell value, in this case 0, instead of "NoData??", "NULL", "" or similar. The files are still perfectly readable, only nodata assignment doesn't happen.

this is as expected. If you don't know the tag, ignore it.

I would prefer to use r.out.gmt to get the GRASS color rules into an extra file (GMT .cpt). This would also be a workaround for export file types other than Byte and UInt16 where GDAL does not support writing out a colortable. But then I don't know how widely GMT .cpt files are supported.

It could be nice to have both internal and .cpt file saved, like current world-file support. But from GRASS or (like world files) from GDAL? FWIW, and it has been a while so I don't remember so well; I don't think I spent much time making r.out.gmt's .cpt output very robust. AFAIR it was a quick hack.

Since it has been some months since the original bug report, and the code has changed, could we please restate what this bug actually is (rather than fuzzy improved r.out.gdal niceness), and some spearfish or NC sample-dataset instructions for recreating it? After studying the ticket from beginning to end I can't say where the GRASS bug is, beyond:

• provide warning if map data overflows requested data type=
• type limits to #defines in the code, and human readable to the help page
• how to make GDAL stop writing color rules for UInt16 once Maximum= is reached (instead of writing out 66535-Max empty rules)

what else?

And we must remember that r.out.gdal is not just for GeoTIFF files..

best, Hamish

Replying to hamish:

Replying to mmetz (and others):

Several of my previous descriptions about the conservative version of r.out.gdal are out of date, e.g. the way statistics are calculated. I would prefer not to stick strictly to my version, instead use it as ideas for discussion.

Spearfish's elevation.dted would be a better sample dataset to use here, it is CELL 0-255. or convert spot.image into RGB bands. Other test datasets are the spearfish imagery60 data, and the North Carolina dataset lsat7 images.

Using sample datasets available to all is a very good idea.

Still, datatypes other than Byte can only be read with spatial data viewers that could make use of this colour table, but not all of them are actually able to do so.

Important to know when this is happening, not really our problem, and not much we can do about it as long as we are writing metadata to spec.

I agree, GRASS can not fix or accommodate other viewers/packages.

It would help a lot to have a list with these data types and their ranges in the man page so that the appropriate data type for file export can be chosen.

Yes. Can someone in the know supply some suggested text?

See the updated documentation in my version.

Maybe the color table is not so important, different people have different preferences

The color table will be carefully constructed/chosen by some people and they may wish to keep it. e.g. r.colors.stddev addon script, r.cpt2grass addon script, and NASA color tables given for MODIS sat data (see grass's MODIS wiki page). These data are real (calibrated into sea surface temperature or mg/mL chlorophyll) and also have specially designed color tables (incl log-scale) which are nice to keep.

True, but if you preprocess such data for someone else who will be running his/her own analyses, a colortable is not important. Sometimes it is desirable to keep the colortable, sometimes not, therefore I would like to have the option to keep it or not.

Suggested data visualization with an embedded color table is nice to have, but coloring schemes are a matter of taste

geotiff has a method for assigning color palettes, all I suggest is populating that correctly, apparently we are currently (or were) failing at.

I disagree, the colortable is populated correctly with the current svn version (at least in my tests), but other GIS software may not be able to make use of it. Nothing GRASS can do about that.

There is still the issue of nodata: a nodata value must fall into the range of the chosen datatype, otherwise this information gets lost. r.out.gdal does not check if a supplied nodata value is valid.

Hopefully better in SVN now. If the user manually specifies an out of range value, I don't mind respecting that.

Users may accidentally specify an out of range nodata value. That has a similar effect like replacing NULL values in the GRASS raster. I would at least like to show a warning.

For UINT16, the color table in the tiff file will always have 65535 entries, gdal does that.

perhaps this should be addressed in GDAL then; only write palette entries which are given.

File an enhancement for GDAL?

Have an option in r.out.gdal to choose whether a color table should be written to the tiff file?

For Uint16, how to just get the "COLOR_TABLE_RULE_RGB_0" part of this without the "(RGB with 66535 entries)" part?

Done in my version.

basic statistics (min, max, mean, stddev) are calculated for the current region, replacing NULL cells with nodata value (default or specified) and then exported with GDALSetRasterStatistics()

how is this better than "gdalinfo -stats" ?

I thought it could be better because other software could read the metadata. But they don't, at least the ones I tested. So there is no advantage in giving stats in the metadata.

I don't understand, you convert NULL to nodata value in stats? So for a CELL map with lots of NULL areas, if I make nodata=100 (well above the number of used palette entries) or leave it at 255, the stats will be skewed towards 100 or 255 instead of just considering the real data? huh? "gdalinfo -stats" seems to do the right thing here.

Right,I changed that and ignore nodata, stats as calculated by my version are now identical to stats calculated by gdalinfo.

colortable is only exported on request (new flag needs to be set) and only if there is an entry in directories colr or colr2

This is wrong, libgis is in charge of handling colors. Use the API not raw file access. For one thing, otherwise you miss out on the default mapping if there isn't a colr/ file.

G_read_colors() returns 0 if no entry in colr or colr2 is present, if an entry is present it returns 1. That's what I used, not raw file access.

Additionally, a warning would be useful if the selected datatype does not cover the range of raster values to be exported

i.e. if the map was CELL 0-744 as in my elev_0 example earlier in this ticket, but the user asked for Byte (0-255). I think that's a very good idea. Start by moving range maxes into header #defines for each type from the currently scattered "nodataval = -1E37f;" etc lines in the r.out.gdal source code.

Also do that for range min's?

If other software ignore valid GeoTIFF metadata tags, then that's their problem and their responsibility to provide an alternative for their users. All we can do is correctly implement the spec. If other software read but incorrectly handle valid GeoTIFF metadata, then that's really their bug to deal with.

Again, at least have the option. You can still tell people using such software that they should try GRASS...

OK, the general question is, should GRASS consider that some packages don't fully support all GeoTIFF or other file type features or not and then strictly stick to the specs?

Apparently COLOR_TABLE_RULE entries are custom metadata not read by other software.

Any confirmation on this? It is custom to GRASS or custom to OSGeo family (GDAL dependents) or suggested by the geotiff spec or ...???

My personal interpretation of the source code of r.out.gdal, needs confirmation. I can only say that the tested packages apart from QGIS don't make use of COLOR_TABLE_RULE_RGB_%d entries.

current region restored on exit

module should not have to do that, it should set WIND_OVERRIDE variable or so, so if module is interrupted by ctrl-C or whatever the temporary change does not become permanent. (but as mentioned earlier, this violates grass's raster module conventions and the feature should not be applied in SVN)

Yes, it violates GRASS raster conventions, but in this case I can't see why it would be desirable to change the resolution on the fly with g.region, i.e. nearest neighbour for raster export. At least it should be stated clearly in the documentation that r.out.gdal respects the current region extends and resolution.

if supported, GDAL decides

what does that mean? when is min/max etc not supported??

I don't no the specs of all raster types supported by GDAL. r.out.gdal can try to export all metadata available, but some other formats support more metadata, some less. It should not be of interest for GRASS, because GDAL decides whether certain metadata are supported or not.

Tests were performed with the elevation raster map layer (range within Byte, FCELL map) in the North Carolina sample dataset, creating a MASK for elevation values <= 70m, needed for testing of nodata handling with nodata=0. The elevation raster was exported as Byte, UInt16, Int32, and Float32, once with colortable (Byte and UInt16 only, colortable export for other data types not supported by gdal), once without colortable, always as GeoTIFF.

Please supply map names and exact commands so we can replicate.

I will do so in a separate reply, this one is long enough already.

QGIS-0.9.1 on Linux 64bit QGIS-0.11.1 on XP ... ESRI ArcMap? 9.2 on XP

this is handy information for users to be aware of, but (again) not much we can do about 3rd party support beyond coding to geotiff spec/recomendations.

In the end, it boils down to colortables in UInt16 files, therefore my request for an option to export colortables. That doesn't harm and can do good.

I would prefer to use r.out.gmt to get the GRASS color rules into an extra file (GMT .cpt). This would also be a workaround for export file types other than Byte and UInt16 where GDAL does not support writing out a colortable. But then I don't know how widely GMT .cpt files are supported.

Since it has been some months since the original bug report, and the code has changed, could we please restate what this bug actually is (rather than fuzzy improved r.out.gdal niceness), and some spearfish or NC sample-dataset instructions for recreating it? After studying the ticket from beginning to end I can't say where the GRASS bug is, beyond: * provide warning if map data overflows requested data type= * type limits to #defines in the code, and human readable to the help page * how to make GDAL stop writing color rules for UInt16 once Maximum= is reached (instead of writing out 66535-Max empty rules)

Ask Frank Warmerdam if this is necessary? This is GDAL code, AFAIK the GDAL API does not provide such an option.

what else?

• an option for colortable export:-)
• Nodata handling: more documentation, warning or error messages if out of range because this would be implicit raster editing.

And we must remember that r.out.gdal is not just for GeoTIFF files..

... but for many many other formats. Too many for testing? Leave as much error handling as possible to GDAL?

In the end, r.out.gdal as it is now may not be as bad as the length of this ticket suggests, and the 3 points of Hamish plus my two points would not require a lot of rewriting.

Regards, Markus

These messages get a bit long so I'll reply in full in the next message.

But quickly-

• oops, you are quite right. nodata must be set to something valid, otherwise the module doesn't know how to write out NULL values when it comes to them. It can't just skip them and nan is not an option for int maps. Thus, as you proposed, in the case where nodata= is out of range for the data type the module should exit with an error.

I was just thinking of what to do for the metadata tag; not how the module would deal with NULL values.

• if a Byte map has data 0-255 *and* NULL, it seems reasonable to me to default to nodata=255 + issue a warning. *But* if a map has 0-255 with *no* NULLs, nodata should not be set* to 255, as this would lead to a loss of fidelity. (improper removal of all cells with cat=255 for no gain)

[*] ie should not be set, unless the user requests it explicitly with the nodata= opt

sorry to be so thick, Hamish

Replying to hamish:

- if a Byte map has data 0-255 *and* NULL, it seems reasonable to me to default to nodata=255 + issue a warning. *But* if a map has 0-255 with *no* NULLs, nodata should not be set* to 255, as this would lead to a loss of fidelity. (improper removal of all cells with cat=255 for no gain)

This is the current implementation in 6.4 (unless modified recently) and my version.

[*] ie should not be set, unless the user requests it explicitly with the nodata= opt

Thicker: I disagree, this is raster editing which should be done before exporting, identical to r.null setnull=val.

As a summary, I see three issues of controversy:

1. should implicit raster editing be allowed, i.e. give and respect nodata in the absence of NULL cells, give out of range nodata value (NULL cells are replaced with min or max of data type, r.null null=nodataval gets "secretly" converted by gdal to r.null null=datatype_min or r.null null=datatype_max), respect the current region resolution

2. should the user be allowed to do strange things like giving out of range nodata values or specify a data type not covering the range of cell values

3. should r.out.gdal consider flaws in (not only) GeoTIFF support of other applications, e.g. give an option to export the colortable and don't do it by default, or should r.out.gdal (stubbornly?) stick to the specs

sorry to be so thick,

Being thick is necessary, r.out.gdal is important for data exchange with people using other applications

Markus

Replying to hamish:

No color table was specified with r.colors, so what data visualization scheme is to be expected?

if none is specified (no colr/ files exists for the raster) the default "rainbow" is used. GRASS's libgis does this automatically, AFAIK the module never knows that r.colors hasn't been run.

It wasn't always this way; prior to r30551, it would only export an explicit colour table, not the default rainbow colour table.

I have just applied in SVN (r34015,6) some updates which should help. I took a lot of your suggestions, but not all.

cleanups in pursuit of bug #73 in collaboration with Markus Metz
* rename & move export band to its own file
* flag to disable writing long (UInt16) color tables
* add type length defines (needs completion & review)
* add some TODO comments where future work is required
* help page cleanups 

Replying to mmetz:

As a summary, I see three issues of controversy: 1. should implicit raster editing be allowed, i.e. give and respect nodata in the absence of NULL cells,

AFAIK this is just the setting of a metadata tag (or not), so not really editing the raster data at all. As such I don't mind if the user decides to set it when it would otherwise not be. There's no processing involved.

give out of range nodata value (NULL cells are replaced with min or max of data type, r.null null=nodataval gets "secretly" converted by gdal to r.null null=datatype_min or r.null null=datatype_max),

It should exit with an error if the user-supplied nodata= value is out of range for the specified data type. When the module comes across a NULL cell it needs to know what to do with it.

respect the current region resolution

note added to the help page that g.region is boss.

2. should the user be allowed to do strange things like giving out of range nodata values or specify a data type not covering the range of cell values

No, G_fatal_error() for both; the module needs to know what to write out for the cell and shouldn't just make stuff up.

3. should r.out.gdal consider flaws in (not only) GeoTIFF support of other applications, e.g. give an option to export the colortable and don't do it by default, or should r.out.gdal (stubbornly?) stick to the specs

-c flag added to turn off (long) color tables for Byte and UInt16. I think it is worth checking if the GDAL code could be changed to only write out as much of the table as needed, not up to 65535-max empty rules.

AFAIAC, upon request we should offer the user the chance to export a bare-bones BASELINE GeoTIFF (we have r.out.tiff for that), but by default we should keep the GeoTIFF as self-documenting and metadata up'd as possible-- as that is what makes the data valuable and viable in the long term and support files are quickly lost.

I'd note that INTERLEAVE=PIXEL is only the default for the very latest versions of GDAL released in the last months. On my system it still defaults to BAND.

The type range limit testing stuff is still TODO. (both data max/min and nodata=)

Hamish

Replying to hamish:

Replying to mmetz:

1. should implicit raster editing be allowed, i.e. give and respect nodata in the absence of NULL cells,

AFAIK this is just the setting of a metadata tag (or not), so not really editing the raster data at all. As such I don't mind if the user decides to set it when it would otherwise not be. There's no processing involved.

This is the same like r.null setnull=nodata. Is this raster editing?

note added to the help page that g.region is boss.

Little word missing? "... if you need to realign the region settings to the original map's before export."

Also maybe remove comments on GDAL datatypes just below the ranges of datatypes, because they are not really necessary?

-c flag added to turn off (long) color tables for Byte and UInt16. I think it is worth checking if the GDAL code could be changed to only write out as much of the table as needed, not up to 65535-max empty rules.

Potentially misleading description: "Do not export long colortable"? A user might wonder whether short colortables are always exported?

The type range limit testing stuff is still TODO. (both data max/min and nodata=)

The min/max values for all integer types are identical to gdal and AFAIK universally valid. See min/max values in gdal-1.5.2 gcore/rasterio.cpp. As I understand, min/max values for float32 are indirectly determined in gdal using typecast from double to float. The min/max values for float32 and float64 as suggested by me correspond to IEE 754 to my best knowledge. Values farthest away from zero have all bits in the mantissa set and all but the last in the exponent which is equal to (1 - 1/2²⁴) * 2¹²⁸ for float32 and (1 - 1/2⁵³) * 2¹⁰²⁴ for float64, taken from http://www.cs.berkeley.edu/~wkahan/ieee754status/IEEE754.PDF.

There is a cplusplus style comment in local_proto.h. OK, it's in a c style comment, but that could be uncommented soon.

In the TODO line 386 rather use GDALGetDataTypeName(datatype) than type->answer because type is not required and datatype is known by then.

Even looking very hard I can't find more to criticize;-)

I guess once the type ranges are confirmed the TODOs can be done.

Markus

Replying to mmetz:

Hamish:

AFAIK this is just the setting of a metadata tag (or not), so not really editing the raster data at all. As such I don't mind if the user decides to set it when it would otherwise not be. There's no processing involved.

Replying to mmetz:

This is the same like r.null setnull=nodata. Is this raster editing?

It is not the same because the raster cells are intact. You can "cleanse" away the metadata with tifftopnm file.gtiff | pnmtotiff > file.tif and the metadata "nodata" tag is gone and your data is as it once was.

Little word missing? "... if you need to realign the region settings to the original map's before export." Also maybe remove comments on GDAL datatypes just below the ranges of datatypes, because they are not really necessary?

thanks, done.

-c flag added

Potentially misleading description: "Do not export long colortable"? A user might wonder whether short colortables are always exported?

they'd be right. Short tables are still written.

• Long: Byte and UInt16 maps write out a RGB value for each 256/66536 possible int.

• Short: These are the COLOR_TABLE_RULE_RGB_n metadata created by

GRASS by copying over the contents of the colr/ file. They are typically 4-5 lines long and work by ranges not per-possibility.

I agree it is unpleasant if the user specifies createopt="PROFILE=BASELINE" for a bare-bones GeoTIFF and GRASS still writes out COLOR_TABLE_RULE metadata anyway.

But typically the user will use that flag to avoid the list of 66536 RGB values, not the list of 5 rules.

The type range limit testing stuff is still TODO. (both data max/min and nodata=)

The min/max values for all integer types are identical to gdal and AFAIK universally valid. See min/max values in gdal-1.5.2 gcore/rasterio.cpp. As I understand, min/max values for float32 are indirectly determined in gdal using typecast from double to float. The min/max values for float32 and float64 as suggested by me correspond to IEE 754 to my best knowledge. Values farthest away from zero have all bits in the mantissa set and all but the last in the exponent which is equal to (1 - 1/2²⁴) * 2¹²⁸ for float32 and (1 - 1/2⁵³) * 2¹⁰²⁴ for float64, taken from http://www.cs.berkeley.edu/~wkahan/ieee754status/IEEE754.PDF.

thanks for the refs., most useful

I notice gdal's gcore/rasterio.cpp has its Int32 min at -2147483647 not -2147483648. That might make a diff if we wanted to use the min as the default nodata value. ?? (or is that < not <= so ok?)

There is a cplusplus style comment in local_proto.h. OK, it's in a c style comment, but that could be uncommented soon.

just a temporary thing, while we sort this out, and less dangerous than a c-comment nested in another c-comment.

In the TODO line 386 rather use GDALGetDataTypeName(datatype) than type->answer because type is not required and datatype is known by then.

I've copied your comment into my local copy of main.c.

I guess once the type ranges are confirmed the TODOs can be done.

right.

I'm no expert on when it is best to terminate a big number with "UL", "L", ".0", "f" or use the likes of "(double)(unsigned int)0xFFFFFFFFu"

I notice with the current code that if I make a UInt16 map and don't specify a nodata value the G_message reports it defaults to -32768, but this may just be %d speaking, not the variable.

Hamish

Replying to hamish:

Replying to mmetz: - Long: Byte and UInt16 maps write out a RGB value for each 256/66536 possible int. - Short: These are the COLOR_TABLE_RULE_RGB_n metadata created by GRASS by copying over the contents of the colr/ file. They are typically 4-5 lines long and work by ranges not per-possibility. I agree it is unpleasant if the user specifies createopt="PROFILE=BASELINE" for a bare-bones GeoTIFF and GRASS still writes out COLOR_TABLE_RULE metadata anyway. But typically the user will use that flag to avoid the list of 66536 RGB values, not the list of 5 rules.

Yes, but for GeoTIFF these short colortables are not colortables in the GeoTIFF sense, but custom metadata. The equivalent in gdal_translate would be gdal_translate -mo "COLOR_TABLE_RULE_0=..." Qgis reads these rules instead of the colortable if rules are present, but qgis is the only application doing so. Grass doesn't do that itself, i.e. it exports COLOR_TABLE_RULE, but ignores these entries on import (enhancement for r.in.gdal?).

I notice gdal's gcore/rasterio.cpp has its Int32 min at -2147483647 not -2147483648. That might make a diff if we wanted to use the min as the default nodata value. ?? (or is that < not <= so ok?)

Hmm, interesting. limits.h defines INT_MAX as 2147483647 and INT_MIN as -INTMAX - 1 Also true: USHRT_MAX (UInt16 max) = SHRT_MAX - SHRT_MIN (Int16 max and min)

UINT_MAX (UInt32 max) = INT_MAX - INT_MIN (Int32 max and min) or INT_MAX + INT_MAX + 1

But if gdal decides Int32 min is -2147483647, then it be.

I'm no expert on when it is best to terminate a big number with "UL", "L", ".0", "f" or use the likes of "(double)(unsigned int)0xFFFFFFFFu"

Not that I am, I would suggest using the same like gdal or limits.h (UInt32 max = 4294967295U)

Float23 and Float 64 limits: on my system (Linux 64 with gcc 4.1.2) there are FLT_MAX and DBL_MAX defined and can be used for TYPE_FLOAT32_MIN/MAX and TYPE_FLOAT64_MIN/MAX. Are FLT_MAX and DBL_MAX defined on other systems (in float.h) too? This would be more accurate than the current approximation in local_proto.h.

I notice with the current code that if I make a UInt16 map and don't specify a nodata value the G_message reports it defaults to -32768, but this may just be %d speaking, not the variable.

If a nodata value is chosen automatically in the presence of NULL cells and absence of user-specified nodata option, this must be correct and consistent, i.e. always max or always min of the datatype which is not the case right now. Should be corrected when the range testing is implemented.

Markus

The attached patch r.out.gdal.rangecheck.patch would add a check whether the user specified nodata value is within the range of the selected GDAL datatype and whether the selected GDAL datatype covers the range of each input raster band. The default nodata value is now set to the minimum of the selected GDAL datatype for all datatypes. Writing out an all black colortable in case that no colortable could be read is disabled. BTW, it seems that some changes in the patch are there because the original code was not indented with grass_indent.sh.

Markus

Replying to mmetz:

BTW, it seems that some changes in the patch are there because the original code was not indented with grass_indent.sh.

It's best if you can avoid mixing substantial changes and formatting changes in the same patch. If the code wasn't formatted either preserve the original formatting (or use e.g. "svn diff -x -b" to ignore whitespace changes) or re-format the code first.

FWIW, I've reformatted r.out.gdal in r34914.

Replying to glynn:

FWIW, I've reformatted r.out.gdal in r34914.

FWIW, it was just minor whitespace on roughly 5 lines. I don't think it is necessary to run intent compulsively as long as the whitespace doesn't get too abused.

Hamish

Replying to hamish:

FWIW, it was just minor whitespace on roughly 5 lines. I don't think it is necessary to run intent compulsively as long as the whitespace doesn't get too abused.

The question is, keep the GRASS source code in standard format yes or no? According to the discussion on it and to SUBMITTING, the answer is something like "YES, we had a long discussion on it, reformatting the whole source code was a major change, now please be so kind as to keep the code in the standard format" and not "Yes, but a little deviation here and there doesn't harm". No offense! My patch would introduce some new code, that's why I ran indent on it, otherwise I would have been told to do so first before proposing a patch.

Justification for the patch: the TODO's in the code have been attended to; that would avoid quite a few of the problems reported in this very long ticket and IMHO the module could be regarded as working all right if these TODO's are solved.

I replaced the attached patch with a new patch against the newly formatted r.out.gdal

Markus M

TODO's attended, patch against r34914

Is this still an issue, or have Markus' patches solved the problem ?

Replying to mlennert:

Is this still an issue, or have Markus' patches solved the problem ?

My patch only adds safety checks to make sure that the selected GDAL data type covers the range of the raster band to be exported and that the given null value is within the range of the selected GDAL data type. That would be an enhancement, not a bugfix.

AFAICT, r.out.gdal as it is now (without the patch) can produce perfectly readable output, but it is possible to produce bogus output without these two safety checks. This is now a trade-off between code complexity and additional checks of user-specified options. I would opt to keep the code simple, also because the manual has been updated.

In any case, it would be great if this ticket could be solved before the final release of 6.4.

Replying to mmetz:

Replying to mlennert:

Is this still an issue, or have Markus' patches solved the problem ?

My patch only adds safety checks to make sure that the selected GDAL data type covers the range of the raster band to be exported and that the given null value is within the range of the selected GDAL data type. That would be an enhancement, not a bugfix. AFAICT, r.out.gdal as it is now (without the patch) can produce perfectly readable output, but it is possible to produce bogus output without these two safety checks. This is now a trade-off between code complexity and additional checks of user-specified options. I would opt to keep the code simple, also because the manual has been updated. In any case, it would be great if this ticket could be solved before the final release of 6.4.

Yes, but what I'm trying to understand is what exactly needs to be solved. From the discussion, I gather that most of the problem is not a bug in r.out.gdal, but rather the incapacity of certain viewers to display certain types of data. r.out.gdal does export an image to be viewed, but raster data for analysis and a viewer cannot handle this.

Or am I misunderstanding something ?

Moritz

Replying to mlennert:

From the discussion, I gather that most of the problem is not a bug in r.out.gdal, but rather the incapacity of certain viewers to display certain types of data.

Yes, that's one problem, but a problem of the viewers. Standard image viewers or editors like GIMP can only read type BYTE, and probably only 1 band or 3 bands in e.g. a tiff. Other types (e.g. INT16, FLOAT32) are not supported by these image viewers. Not sure about the number of bands.

GIS software can usually both read, analyse, and visualize georeferenced raster data as produced by r.out.gdal. The tiff colortable problem is already solved. Some tests are in comment 13:

http://trac.osgeo.org/grass/ticket/73#comment:13

Replying to mmetz:

Replying to mlennert:

From the discussion, I gather that most of the problem is not a bug in r.out.gdal, but rather the incapacity of certain viewers to display certain types of data.

Yes, that's one problem, but a problem of the viewers. Standard image viewers or editors like GIMP can only read type BYTE, and probably only 1 band or 3 bands in e.g. a tiff. Other types (e.g. INT16, FLOAT32) are not supported by these image viewers. Not sure about the number of bands. GIS software can usually both read, analyse, and visualize georeferenced raster data as produced by r.out.gdal. The tiff colortable problem is already solved. Some tests are in comment 13: http://trac.osgeo.org/grass/ticket/73#comment:13

Well, is this still an open bug then, or can we close it ? Or at least downgrade it so that it doesn't block 6.4 release ?

Moritz

Replying to mlennert:

Well, is this still an open bug then, or can we close it ? Or at least downgrade it so that it doesn't block 6.4 release ?

Wait for feedback from Helena, Eric, Dylan, and Hamish before downgrading or closing the ticket ?

This was and is my main concern. MarkusM's patch could address this, I haven't run exaustive tests so don't know.

Replying to hamish:

- oops, you are quite right. nodata must be set to something valid, otherwise the module doesn't know how to write out NULL values when it comes to them. It can't just skip them and nan is not an option for int maps. Thus, as you proposed, in the case where nodata= is out of range for the data type the module should exit with an error. I was just thinking of what to do for the metadata tag; not how the module would deal with NULL values. - if a Byte map has data 0-255 *and* NULL, it seems reasonable to me to default to nodata=255 + issue a warning. *But* if a map has 0-255 with *no* NULLs, nodata should not be set* to 255, as this would lead to a loss of fidelity. (improper removal of all cells with cat=255 for no gain) [*] ie should not be set, unless the user requests it explicitly with the nodata= opt sorry to be so thick, Hamish

As for the other matter, I'm not worried if other software do not respect the GeoTiff? spec, as long as we do. (it's all we have control over anyway)

more in next post..

Replying to mmetz:

My patch only adds safety checks to make sure that the selected GDAL data type covers the range of the raster band to be exported and that the given null value is within the range of the selected GDAL data type. That would be an enhancement, not a bugfix.

Pretend you are r.out.gdal iterating through the raster cells to output. You come across a NULL cell. You have to write down some value, and are completely constrained to options available in the output type. What do you write out? Overflow/wrap the type? type min? type max? well (as things were) whatever you do is going to be arbirtary, and thus could lead to data loss, which is grounds for a most serious blocker bug.

also, because the manual has been updated.

that is easily modified, so shouldn't stop us.

In any case, it would be great if this ticket could be solved before the final release of 6.4.

yup

Hamish

Replying to hamish:

Pretend you are r.out.gdal iterating through the raster cells to output. You come across a NULL cell. You have to write down some value, and are completely constrained to options available in the output type. What do you write out? Overflow/wrap the type? type min? type max? well (as things were) whatever you do is going to be arbirtary, and thus could lead to data loss, which is grounds for a most serious blocker bug.

That would require an additional check that the nodata value is not present in the raster band to be exported. This is not too difficult.

What would be more difficult is to select a nodata value that is not present in the raster band to be exported if type min and type max are present. Offhand, I can only think of a helper data structure (sorted heap/binary tree/hash table) where all raster values without duplicates are kept which can consume a lot of memory for large floating point rasters. This would also require two passes over each raster band, first investigating the values and selecting a proper nodata value, second replacing NULLs with this nodata value. Then, what to do if a user-specified nodata value is within the range of the selected GDAL data type but this value is present in the raster? Give a warning but let the user override the check, causing data loss?

Sorry, I haven't thought about this problem of nodata value causing data loss, somehow left this to the user.

Markus M

Replying to mmetz:

Sorry, I haven't thought about this problem of nodata value causing data loss, somehow left this to the user.

And I think it should be. Just give a big fat warning in the man page, and possibly as output of the module saying which value is used for null data and that if there are non-null pixels which have the value given as nodata value, then you will lose the information contained in them.

Moritz

Replying to mmetz:

Replying to mlennert:

Well, is this still an open bug then, or can we close it ? Or at least downgrade it so that it doesn't block 6.4 release ?

Wait for feedback from Helena, Eric, Dylan, and Hamish before downgrading or closing the ticket ?

I suppose there's not much we can do if other image viewers are not configured to read proper geotiffs. For the record, the only tiff output I can create that actually is readable (i.e., color table remains intact) is r.out.tiff. Maybe allow a flag in r.out.gdal that outputs r.out.tiff-style tiffs and retire r.out.tiff in Grass 7?

~ Eric.

as for me - I did not file it as critical, because I thought it was partially problem on my side (not using the correct options) - it was upgraded to critical later.

I think we just need to keep around r.out.arc, r.out.tiff in case the software that you are trying to import the data to does not like the r.out.gdal geotiff as we cannot really solve other software issues with geotiff. I also agree that the warning about possible loss of data with badly chosen NULL as suggested by Moritz should take care of the null problem.

All I really tried to do was convert the nc_spm raster data into something readable by other software tools and r.out.gdal worked but some rasters were not readable or came out black so it required change in color table etc. (e.g., in ArcGIS, gimp). We ended up using r.out.arc - and lost the color tables and projection info but otherwise it worked fine. I haven't tried the very latest version of r.out.gdal though.

Helena

Replying to helena:

All I really tried to do was convert the nc_spm raster data into something readable by other software tools and r.out.gdal worked but some rasters were not readable or came out black so it required change in color table etc. (e.g., in ArcGIS, gimp). We ended up using r.out.arc - and lost the color tables and projection info but otherwise it worked fine.

Just out of curiosity: how is losing color table and projection info (r.out.arc) better than just a lost color table (r.out.gdal) ?

Moritz

An attempt to solve this ticket: updated r.out.gdal in trunk r36750 and develbranch_6 r36751

The manual has been updated considerably, trying to find a compromise between useful information and a lengthy manual. Amongst other additions, there is a new paragraph on GeoTIFF caveats in NOTES, I hope that helps a bit.

Markus M

Replying to epatton:

I suppose there's not much we can do if other image viewers are not configured to read proper geotiffs. For the record, the only tiff output I can create that actually is readable (i.e., color table remains intact) is r.out.tiff. Maybe allow a flag in r.out.gdal that outputs r.out.tiff-style tiffs and retire r.out.tiff in Grass 7?

r.out.gdal and r.out.tiff do two fundamentally different things and IMHO should be neither compared nor merged.

r.out.tiff exports an image that can be viewed with any image viewer and looks like what you see on the grass display. The pixel values are most likely very different from the raster map that is exported, because with the -p flag everything is converted to Byte. Without the -p flag, the pixel values are the colour values of the color table, not the raster values. r.out.tiff output is thus not meant to be suitable for spatial analysis (like a screenshot of the display).

r.out.gdal exports data, by default with a color table for data visualization if supported by the selected output format. The color table is just for convenience and not needed for analysis. Data export was successful if r.univar or some equivalent in another GIS software package gives results identical to r.univar on the grass raster. There is a difference between "can read the GeoTIFF" and "displays the GeoTIFF like grass displays the original raster", therefore the display is not suitable for testing because there are too many reasons why the display can be different (different monitor with different settings, different OS with different gamma, file format does not support color tables, conversion from color rules to color tables can introduce slight modifications, other GIS software uses different graphics engine, other GIS software can not read embedded color table).

Sorry for being so thick, but I want to make the point of distinguishing images and data. Unfortunately, a .tif file can be both.

Markus M

Replying to mmetz:

r.out.gdal and r.out.tiff do two fundamentally different things and IMHO should be neither compared nor merged. r.out.tiff exports an image that can be viewed with any image viewer and looks like what you see on the grass display. The pixel values are most likely very different from the raster map that is exported, because with the -p flag everything is converted to Byte. Without the -p flag, the pixel values are the colour values of the color table, not the raster values. r.out.tiff output is thus not meant to be suitable for spatial analysis (like a screenshot of the display). r.out.gdal exports data, by default with a color table for data visualization if supported by the selected output format. The color table is just for convenience and not needed for analysis. Data export was successful if r.univar or some equivalent in another GIS software package gives results identical to r.univar on the grass raster. There is a difference between "can read the GeoTIFF" and "displays the GeoTIFF like grass displays the original raster", therefore the display is not suitable for testing because there are too many reasons why the display can be different (different monitor with different settings, different OS with different gamma, file format does not support color tables, conversion from color rules to color tables can introduce slight modifications, other GIS software uses different graphics engine, other GIS software can not read embedded color table). Sorry for being so thick, but I want to make the point of distinguishing images and data. Unfortunately, a .tif file can be both.

thanks for the clarification. I didn't know that. Probably it should be explained more deeply in the manual page.

Martin

my 2c o'today,

• if type is UInt8 and input map has values outside of that, just report an error and exit. you don't have to search every cell first, just use G_read_range()+G_get_range_min_max() or G_read_fp_range()+G_get_fp_range_min_max(). tricks like "if(val<type_min) val=(type_min||NULL);" are evil. the user can use r.mapcalc or r.reclass if they really want that. better we refuse to honour broken requests.

• exit with an error if user supplied nodata= is outside of range for the given type. if no value is given use the following rules:
  • signed int types: minimum in range
  • unsigned int types: maximum in range
  • all floating point: IEEE's NaN

• ?? can we wait to tell gdal the nodata value (for the metatag) until after we have processed all rows?? if so we can just not set one with GDAL at all if we didn't come across any NULL cells in the write loop. We must decide what to use if we do come across one before we start the loop (see previous point) that would allow for a 0-255 GRASS CELL map with no NULLs to preserve all data, while allowing a CELL map with 0-255+NULL, type=UInt8, and the user specified a nodata= value in the range of 0-255 to correctly follow the user's wishes. If we do that then we don't need a -f flag (which I'd rather we didn't have).

• drop the complex types in grass7. we can easily re-add them if the next-gen raster format can store multiband i.fft output maps or whatever. has anyone ever once come across in the wild a map using these types? or are they on the list just because they can be?

• if writing a message saying "NULL cells were found in the map but the nodata option was not specified. Using a default value of xxxxx." beware that what is printed by the G_warning() or G_message() text may not be the internally stored value, as it has been passed through printf %{}. so don't trust your own eyes looking at that output + we should carefully test that. I guess it means supplying a different %ld %u or whatever for each individual type, but so be it.

[crossing-threads]

• 999 as NULL is dangerous, it can easily be an elevation value. -9999 is bad too, you just hit it less often (& data isn't always elevation). But even though you hit it less often doesn't mean it's still not going to bite you one day, so it's a bug. bugs which don't happen all the time could be considered more evil, because you have come to trust the program and are not expecting them. but alas, probably for r.out.arc we are stuck with whatever values Arc will take for input.

Hamish

Replying to hamish:

my 2c o'today, * if type is UInt8 and input map has values outside of that, just report an error and exit. you don't have to search every cell first, just use G_read_range()+G_get_range_min_max() or G_read_fp_range()+G_get_fp_range_min_max().

Note that the current region is respected and the actual data range to be exported can be smaller than reported by G_get_fp_range_min_max(). OTOH, taking that into consideration may be overdoing it a bit.

* exit with an error if user supplied nodata= is outside of range for the given type.

Glynn suggested interpretation: if (nodata != (double) (GDAL datatype) nodata) -> warning, nodata = (double) (GDAL datatype) nodata, proceed. (I hope I got that right for a change...) Result would be a valid raster, granted that subsequent checks are passed.

* ?? can we wait to tell gdal the nodata value (for the metatag) until after we have processed all rows?? if so we can just not set one with GDAL at all if we didn't come across any NULL cells in the write loop.

We have waited all the time :-) That was always done in the C version and AFAICT nobody wants to change that.

We must decide what to use if we do come across one before we start the loop (see previous point) that would allow for a 0-255 GRASS CELL map with no NULLs to preserve all data, while allowing a CELL map with 0-255+NULL, type=UInt8, and the user specified a nodata= value in the range of 0-255 to correctly follow the user's wishes.

0-255+NULL, type=UInt8 causes data loss unless the user found a free unused slot within 0-255. IMHO, big warning if not error on data loss i.e. when nodata value encountered *and* NULLs present.

I'll try to summarize the suggestions so far into a description of the new design for r.out.gdal, but I need a bit more time to read again through all contributions.

Markus M

Trying to summarize the suggestions and discussions for a new design:

• before actual export new precision test for each band
  • FCELL/DCELL to be exported as some GDAL integer type -> warning
  • CELL exported as GDT_FLOAT32 and raster_min < -2²⁴ or raster_max > 2²⁴ -> warning
  • DCELL exported as GDT_Float32 -> warning

• before actual export new range test for each band
  • if raster_min >= type_min and raster_max =< type_max -> fine, proceed
  • else if raster_min > type_max or raster_max < type_min -> complete data loss, error
  • else there is a partial data range overlap -> check data range during actual export

• before actual export, get a reasonable default nodata value if none given
  • NaN for all GDAL floating point datatypes. Is NaN constructed with 0.0/0.0 ok?
  • GDAL signed int types: first try minimum in potential range, if raster_min <= type_min, try maximum in potential range, if raster_max >= type_max, use minimum (would be (GInt32) 0x80000000 for GDT_Int32). This can re-use results of range check above.
  • GDAL unsigned int types: first try maximum in potential range, if raster_max >= type_max, try minimum in potential range, if raster_min <= type_min, use maximum. This can re-use results of range check above.
  • no warnings or errors at this stage because it is yet unknown if there are any NULL cells in the (with g.region selected part of the) raster maps to be exported

• before actual export, in case of custom nodata make sure the metadata nodata value and the raster nodata value are identical
  • if (nodata != (double) (GDAL datatype) nodata) -> warning and nodata = (double) (GDAL datatype) nodata
  • This would give feedback to the user about what GDAL has to do later for the selected export datatype and nodata value. And the metadata nodata value would always be identical to the value replacing NULL cells.

• during actual export, check for presence of cells with nodata value
  • NULL cells present: message about nodata value used to replace NULLs
  • if there are cells == nodata value and NULL cells were assigned that nodata value -> -f flag and user nodata value: warning, else error

• during actual export, check actual raster range against range of GDAL datatype if needed (depends on result of range check above)
  • raster_min < type_min or raster_max > type_max -> data loss, error

as in previous C versions of r.out.gdal, only write GDAL nodata to metadata if NULLs were encountered. It seems that GDAL metadata nodata is raster band specific, not globally set for the whole file; r.out.gdal behaves accordingly.

All warnings and errors should clearly explain what is the problem and give tips on how it can be solved.

Hamish, what exactly should this compatibility flag do? There is all sorts of software with all sorts of different deficiencies out there...

I hope I didn't forget something

Markus M

All changes described in comment 49 are submitted in trunk r37764. The -f flag overrides the precision test before actual export and the nodata value test during actual export.

Markus M

Replying to mmetz:

- NaN for all GDAL floating point datatypes. Is NaN constructed with 0.0/0.0 ok?

yes (AFAIK)

- GDAL signed int types: first try minimum in potential range, if raster_min <= type_min, try maximum in potential range, if raster_max >= type_max, use minimum (would be (GInt32) 0x80000000 for GDT_Int32). This can re-use results of range check above. - GDAL unsigned int types: first try maximum in potential range, if raster_max >= type_max, try minimum in potential range, if raster_min <= type_min, use maximum. This can re-use results of range check above.

it all seems a bit complicated, but ok.

* before actual export, in case of custom nodata make sure the metadata nodata value and the raster nodata value are identical

why? if custom nodata then export NULLs in the map to be the custom value and clobber any real data which had that value. (perhaps I don't understand something here..)

- if (nodata != (double) (GDAL datatype) nodata) -> warning and nodata = (double) (GDAL datatype) nodata

if the user asks to use a certain nodata value and it is illegal for the data type then exit with an error, probably giving the available range in the error message. don't automagically correct it for them and continue (ie override their expressed wishes). It is a recipe for pain.

- if there are cells == nodata value and NULL cells were assigned that nodata value -> -f flag and user nodata value: warning, else error

ok

Hamish, what exactly should this compatibility flag do? There is all sorts of software with all sorts of different deficiencies out there...

umm, I forget what that was in reference to. ??minimalistic metadata output??

is everyone happy with the colortable export stuff now? (my only issue with it is that if you pass the no-metadata create option to GDAL, GRASS adds its stuff anyway)

I shifted the bug target to 6.5 as this needs testing before going into 6.4.0 (this is a substantial last-minute change to a core module & we can't keep on resetting the RC cycle + we seem to have survived with it as-was for this long...). Strong candidate for 6.4.1.

thanks, Hamish

Replying to hamish:

Replying to mmetz:

- GDAL signed int types: (...)

- GDAL unsigned int types: (...)

it all seems a bit complicated, but ok.

that's my understanding of the discussion with Glynn: http://lists.osgeo.org/pipermail/grass-dev/2009-April/043464.html and following.

* before actual export, in case of custom nodata make sure the metadata nodata value and the raster nodata value are identical

why? if custom nodata then export NULLs in the map to be the custom value and clobber any real data which had that value. (perhaps I don't understand something here..)

What I meant here was to make sure that the value that replaces NULL cells is identical to the metadata nodata value, otherwise nodata info gets lost. Clobbering real data with the nodata value is a different issue (that is taken care of).

- if (nodata != (double) (GDAL datatype) nodata) -> warning and nodata = (double) (GDAL datatype) nodata

if the user asks to use a certain nodata value and it is illegal for the data type then exit with an error, probably giving the available range in the error message. don't automagically correct it for them and continue (ie override their expressed wishes). It is a recipe for pain.

ok, will change it. The first warning will tell what would happen to the nodata value, the second warning will report the selected GDAL datatype and its range, then raster export aborts. E.g.

WARNING: Mismatch between metadata nodata value and actual nodata value in
         exported raster: specified nodata value -9999.000000 gets
         converted to 241 by selected GDAL datatype.
WARNING: GDAL datatype: Byte, valid range: 0 - 255
ERROR: Raster export aborted.

BTW, it would be the selected GDAL datatype that overrides the user's wishes by type casting, not r.out.gdal.

Hamish, what exactly should this compatibility flag do? There is all sorts of software with all sorts of different deficiencies out there...

umm, I forget what that was in reference to. ??minimalistic metadata output??

There are already GeoTIFF compatibility hints in the manual, I would like to leave it like that.

is everyone happy with the colortable export stuff now? (my only issue with it is that if you pass the no-metadata create option to GDAL, GRASS adds its stuff anyway)

Please note that nothing much changed there, colortable export can be disabled from 6.4 onwards. AFAIKT, that GRASS stuff is currently only used by QGIS, not even GRASS i.e. r.in.gdal or r.external use that information.

this is a substantial last-minute change to a core module

Unfortunately, but I would really like to get this ticket solved, please test!

Markus M

see also #788.