[GRASS-SVN] r40120 - grass/branches/develbranch_6/raster/r.out.gdal

svn_grass at osgeo.org svn_grass at osgeo.org
Wed Dec 23 07:51:15 EST 2009


Author: mmetz
Date: 2009-12-23 07:51:14 -0500 (Wed, 23 Dec 2009)
New Revision: 40120

Modified:
   grass/branches/develbranch_6/raster/r.out.gdal/description.html
   grass/branches/develbranch_6/raster/r.out.gdal/export_band.c
   grass/branches/develbranch_6/raster/r.out.gdal/local_proto.h
   grass/branches/develbranch_6/raster/r.out.gdal/main.c
Log:
check options before actual export, backport from trunk r40119

Modified: grass/branches/develbranch_6/raster/r.out.gdal/description.html
===================================================================
--- grass/branches/develbranch_6/raster/r.out.gdal/description.html	2009-12-23 12:49:43 UTC (rev 40119)
+++ grass/branches/develbranch_6/raster/r.out.gdal/description.html	2009-12-23 12:51:14 UTC (rev 40120)
@@ -108,7 +108,7 @@
 <em>r.out.gdal</em> exports may appear all black or gray on initial
 display in other GIS software. This is not a bug of <em>r.out.gdal</em>,
 but often caused by the default color table assigned by that software.
-The dafault color table may be grayscale covering the whole range of
+The default color table may be grayscale covering the whole range of
 possible values which is very large for e.g. Int32 or Float32. E.g.
 stretching the color table to actual min/max would help (sometimes under
 symbology).
@@ -132,7 +132,7 @@
 
 <h4>Improving GeoTIFF compatibility</h4>
 
-To create a GeoTIFF that is highly compatibility with various other GIS
+To create a GeoTIFF that is highly compatible with various other GIS
 software packages, it is recommended to keep the GeoTIFF file as simple
 as possible. You will have to experiment with which options your
 software is compatible with, as this varies widely between vendors and

Modified: grass/branches/develbranch_6/raster/r.out.gdal/export_band.c
===================================================================
--- grass/branches/develbranch_6/raster/r.out.gdal/export_band.c	2009-12-23 12:49:43 UTC (rev 40119)
+++ grass/branches/develbranch_6/raster/r.out.gdal/export_band.c	2009-12-23 12:51:14 UTC (rev 40120)
@@ -6,7 +6,7 @@
 * PURPOSE:      Exports GRASS raster to GDAL suported formats;
 *               based on GDAL library.
 *
-* COPYRIGHT:    (C) 2006-2008 by the GRASS Development Team
+* COPYRIGHT:    (C) 2006-2009 by the GRASS Development Team
 *
 *               This program is free software under the GNU General Public
 *   	    	License (>=v2). Read the file COPYING that comes with GRASS
@@ -15,20 +15,225 @@
 *****************************************************************************/
 
 #include <grass/gis.h>
+#include <grass/raster.h>
 #include <grass/glocale.h>
 
 #include "cpl_string.h"
 #include "gdal.h"
 #include "local_proto.h"
 
+int exact_range_check(double, double, GDALDataType, const char *);
 
-int export_band(GDALDatasetH hMEMDS, int band, const char *name,
-		const char *mapset, struct Cell_head *cellhead,
-		RASTER_MAP_TYPE maptype, double nodataval,
-		const char *nodatakey, int suppress_main_colortable,
+/* exact check for each band
+ * returns 0 on success
+ * -1 if given nodata value was present in data
+ * -2 if selected GDAL datatype could not hold all values
+ * */
+int exact_checks(GDALDataType export_datatype,
+		const char *name, const char *mapset,
+		struct Cell_head *cellhead, RASTER_MAP_TYPE maptype,
+		double nodataval, const char *nodatakey,
 		int default_nodataval)
 {
+    int bHaveMinMax;
+    double dfCellMin;
+    double dfCellMax;
+    struct FPRange sRange;
+    int fd;
+    int cols = cellhead->cols;
+    int rows = cellhead->rows;
+    int ret = 0;
 
+    /* Open GRASS raster */
+    fd = G_open_cell_old(name, mapset);
+    if (fd < 0) {
+	G_warning(_("Unable to open raster map <%s>"), name);
+	return -1;
+    }
+
+    /* Get min/max values. */
+    if (G_read_fp_range(name, mapset, &sRange) == -1) {
+	bHaveMinMax = FALSE;
+    }
+    else {
+	bHaveMinMax = TRUE;
+	G_get_fp_range_min_max(&sRange, &dfCellMin, &dfCellMax);
+    }
+
+    /* Create GRASS raster buffer */
+    void *bufer = G_allocate_raster_buf(maptype);
+
+    if (bufer == NULL) {
+	G_warning(_("Unable to allocate buffer for reading raster map"));
+	return -1;
+    }
+    char *nulls = (char *)G_malloc(cols);
+
+    if (nulls == NULL) {
+	G_warning(_("Unable to allocate buffer for reading raster map"));
+	return -1;
+    }
+
+    /* Copy data form GRASS raster to GDAL raster */
+    int row, col;
+    int n_nulls = 0, nodatavalmatch = 0;
+
+    dfCellMin = TYPE_FLOAT64_MAX;
+    dfCellMax = TYPE_FLOAT64_MIN;
+
+    /* Better use selected GDAL datatype instead of 
+     * the best match with GRASS raster map types ? */
+
+    if (maptype == FCELL_TYPE) {
+
+	FCELL fnullval = (FCELL) nodataval;
+
+	G_debug(1, "FCELL nodata val: %f", fnullval);
+
+	for (row = 0; row < rows; row++) {
+
+	    if (G_get_raster_row(fd, bufer, row, maptype) < 0) {
+		G_warning(_("Unable to read raster map <%s> row %d"),
+			  name, row);
+		return -1;
+	    }
+	    G_get_null_value_row(fd, nulls, row);
+	    for (col = 0; col < cols; col++) {
+		if (nulls[col]) {
+		    ((FCELL *) bufer)[col] = fnullval;
+		    n_nulls++;
+		}
+		else {
+		    if (((FCELL *) bufer)[col] == fnullval) {
+			nodatavalmatch = 1;
+		    }
+		    if (dfCellMin > ((FCELL *) bufer)[col])
+			dfCellMin = ((FCELL *) bufer)[col];
+		    if (dfCellMax < ((FCELL *) bufer)[col])
+			dfCellMax = ((FCELL *) bufer)[col];
+		}
+	    }
+	    G_percent(row + 1, rows, 2);
+	}
+    }
+    else if (maptype == DCELL_TYPE) {
+
+	DCELL dnullval = (DCELL) nodataval;
+
+	G_debug(1, "DCELL nodata val: %f", dnullval);
+
+	for (row = 0; row < rows; row++) {
+
+	    if (G_get_raster_row(fd, bufer, row, maptype) < 0) {
+		G_warning(_("Unable to read raster map <%s> row %d"),
+			  name, row);
+		return -1;
+	    }
+	    G_get_null_value_row(fd, nulls, row);
+	    for (col = 0; col < cols; col++) {
+		if (nulls[col]) {
+		    ((DCELL *) bufer)[col] = dnullval;
+		    n_nulls++;
+		}
+		else {
+		    if (((DCELL *) bufer)[col] == dnullval) {
+			nodatavalmatch = 1;
+		    }
+		    if (dfCellMin > ((DCELL *) bufer)[col])
+			dfCellMin = ((DCELL *) bufer)[col];
+		    if (dfCellMax < ((DCELL *) bufer)[col])
+			dfCellMax = ((DCELL *) bufer)[col];
+		}
+	    }
+	    G_percent(row + 1, rows, 2);
+	}
+    }
+    else {
+
+	CELL inullval = (CELL) nodataval;
+
+	G_debug(1, "CELL nodata val: %d", inullval);
+
+	for (row = 0; row < rows; row++) {
+
+	    if (G_get_raster_row(fd, bufer, row, maptype) < 0) {
+		G_warning(_("Unable to read raster map <%s> row %d"),
+			  name, row);
+		return -1;
+	    }
+	    G_get_null_value_row(fd, nulls, row);
+	    for (col = 0; col < cols; col++) {
+		if (nulls[col]) {
+		    ((CELL *) bufer)[col] = inullval;
+		    n_nulls++;
+		}
+		else {
+		    if (((CELL *) bufer)[col] == inullval) {
+			nodatavalmatch = 1;
+		    }
+		    if (dfCellMin > ((CELL *) bufer)[col])
+			dfCellMin = ((CELL *) bufer)[col];
+		    if (dfCellMax < ((CELL *) bufer)[col])
+			dfCellMax = ((CELL *) bufer)[col];
+		}
+	    }
+	    G_percent(row + 1, rows, 2);
+	}
+    }
+
+    /* can the GDAL datatype hold the data range to be exported ? */
+    /* f-flag does not override */
+    if (exact_range_check(export_datatype, dfCellMin, dfCellMax, name)) {
+	G_warning("Raster export results in data loss.");
+	ret = -2;
+    }
+
+    /* a default nodata value was used and NULL cells were present */
+    if (n_nulls && default_nodataval) {
+	if (maptype == CELL_TYPE)
+	    G_important_message(_("Input raster map contains cells with NULL-value (no-data). "
+				 "The value %d will be used to represent no-data values in the input map. "
+				 "You can specify a nodata value with the %s option."),
+				(int)nodataval, nodatakey);
+	else
+	    G_important_message(_("Input raster map contains cells with NULL-value (no-data). "
+				 "The value %f will be used to represent no-data values in the input map. "
+				 "You can specify a nodata value with the %s option."),
+				nodataval, nodatakey);
+    }
+
+    /* the nodata value was present in the exported data */
+    if (nodatavalmatch && n_nulls) {
+	/* default nodataval didn't work */
+	if (default_nodataval) {
+	    G_warning(_("The default nodata value is present in raster"
+			"band <%s> and would lead to data loss. Please specify a "
+			"custom nodata value with the %s parameter."),
+		      name, nodatakey);
+	}
+	/* user-specified nodataval didn't work */
+	else {
+	    G_warning(_("The given nodata value is present in raster"
+			"band <%s> and would lead to data loss. Please specify a "
+			"different nodata value with the %s parameter."),
+		      name, nodatakey);
+	}
+	ret = -1;
+    }
+
+    return ret;
+}
+
+/* actual raster band export
+ * returns 0 on success
+ * -1 on raster data read/write error
+ * */
+int export_band(GDALDatasetH hMEMDS, GDALDataType export_datatype, int band,
+		const char *name, const char *mapset,
+		struct Cell_head *cellhead, RASTER_MAP_TYPE maptype,
+		double nodataval, const char *nodatakey,
+		int suppress_main_colortable, int default_nodataval)
+{
     struct Colors sGrassColors;
     GDALColorTableH hCT;
     int iColor;
@@ -39,6 +244,7 @@
     int fd;
     int cols = cellhead->cols;
     int rows = cellhead->rows;
+    int ret = 0;
 
     /* Open GRASS raster */
     fd = G_open_cell_old(name, mapset);
@@ -113,15 +319,15 @@
 		int nRed, nGreen, nBlue;
 		GDALColorEntry sColor;
 
-		if (G_get_color
-		    (iColor, &nRed, &nGreen, &nBlue, &sGrassColors)) {
+		if (G_get_color(iColor, &nRed, &nGreen, &nBlue,
+				     &sGrassColors)) {
 		    sColor.c1 = nRed;
 		    sColor.c2 = nGreen;
 		    sColor.c3 = nBlue;
 		    sColor.c4 = 255;
 
 		    G_debug(3,
-			    "G_get_color: Y, rcount %d, nRed %d, nGreen %d, nBlue %d",
+			    "Rast_get_c_color: Y, rcount %d, nRed %d, nGreen %d, nBlue %d",
 			    rcount, nRed, nGreen, nBlue);
 		    GDALSetColorEntry(hCT, iColor, &sColor);
 		}
@@ -132,11 +338,13 @@
 		    sColor.c4 = 0;
 
 		    G_debug(3,
-			    "G_get_color: N, rcount %d, nRed %d, nGreen %d, nBlue %d",
+			    "Rast_get_c_color: N, rcount %d, nRed %d, nGreen %d, nBlue %d",
 			    rcount, nRed, nGreen, nBlue);
 		    GDALSetColorEntry(hCT, iColor, &sColor);
 		}
 	    }
+
+	    GDALSetRasterColorTable(hBand, hCT);
 	}
 
 	if (rcount > 0) {
@@ -147,6 +355,8 @@
 	}
 
 	/* Add the rules in reverse order */
+	/* This can cause a GDAL warning with many rules, something like
+	 * Warning 1: Lost metadata writing to GeoTIFF ... too large to fit in tag. */
 	for (i = rcount - 1; i >= 0; i--) {
 	    DCELL val1, val2;
 	    unsigned char r1, g1, b1, r2, g2, b2;
@@ -160,9 +370,6 @@
 		    r2, g2, b2);
 	    GDALSetMetadataItem(hBand, key, value, NULL);
 	}
-
-	if (!suppress_main_colortable)
-	    GDALSetRasterColorTable(hBand, hCT);
     }
 
     /* Create GRASS raster buffer */
@@ -179,16 +386,21 @@
 	return -1;
     }
 
-    /* Copy data form GRASS raster to memory raster */
+    /* Copy data form GRASS raster to GDAL raster */
     int row, col;
-    int n_nulls = 0, nodatavalmatch = 0;
+    int n_nulls = 0;
 
+    /* Better use selected GDAL datatype instead of 
+     * the best match with GRASS raster map types ? */
+
     if (maptype == FCELL_TYPE) {
 
 	/* Source datatype understandable by GDAL */
 	GDALDataType datatype = GDT_Float32;
 	FCELL fnullval = (FCELL) nodataval;
 
+	G_debug(1, "FCELL nodata val: %f", fnullval);
+
 	for (row = 0; row < rows; row++) {
 
 	    if (G_get_raster_row(fd, bufer, row, maptype) < 0) {
@@ -197,7 +409,7 @@
 		return -1;
 	    }
 	    G_get_null_value_row(fd, nulls, row);
-	    for (col = 0; col < cols; col++)
+	    for (col = 0; col < cols; col++) {
 		if (nulls[col]) {
 		    ((FCELL *) bufer)[col] = fnullval;
 		    if (n_nulls == 0) {
@@ -205,9 +417,7 @@
 		    }
 		    n_nulls++;
 		}
-		else if (((FCELL *) bufer)[col] == fnullval) {
-		    nodatavalmatch = 1;
-		}
+	    }
 
 	    if (GDALRasterIO
 		(hBand, GF_Write, 0, row, cols, 1, bufer, cols, 1, datatype,
@@ -223,6 +433,8 @@
 	GDALDataType datatype = GDT_Float64;
 	DCELL dnullval = (DCELL) nodataval;
 
+	G_debug(1, "DCELL nodata val: %f", dnullval);
+
 	for (row = 0; row < rows; row++) {
 
 	    if (G_get_raster_row(fd, bufer, row, maptype) < 0) {
@@ -231,7 +443,7 @@
 		return -1;
 	    }
 	    G_get_null_value_row(fd, nulls, row);
-	    for (col = 0; col < cols; col++)
+	    for (col = 0; col < cols; col++) {
 		if (nulls[col]) {
 		    ((DCELL *) bufer)[col] = dnullval;
 		    if (n_nulls == 0) {
@@ -239,9 +451,7 @@
 		    }
 		    n_nulls++;
 		}
-		else if (((DCELL *) bufer)[col] == dnullval) {
-		    nodatavalmatch = 1;
-		}
+	    }
 
 	    if (GDALRasterIO
 		(hBand, GF_Write, 0, row, cols, 1, bufer, cols, 1, datatype,
@@ -257,6 +467,8 @@
 	GDALDataType datatype = GDT_Int32;
 	CELL inullval = (CELL) nodataval;
 
+	G_debug(1, "CELL nodata val: %d", inullval);
+
 	for (row = 0; row < rows; row++) {
 
 	    if (G_get_raster_row(fd, bufer, row, maptype) < 0) {
@@ -265,7 +477,7 @@
 		return -1;
 	    }
 	    G_get_null_value_row(fd, nulls, row);
-	    for (col = 0; col < cols; col++)
+	    for (col = 0; col < cols; col++) {
 		if (nulls[col]) {
 		    ((CELL *) bufer)[col] = inullval;
 		    if (n_nulls == 0) {
@@ -273,9 +485,7 @@
 		    }
 		    n_nulls++;
 		}
-		else if (((CELL *) bufer)[col] == inullval) {
-		    nodatavalmatch = 1;
-		}
+	    }
 
 	    if (GDALRasterIO
 		(hBand, GF_Write, 0, row, cols, 1, bufer, cols, 1, datatype,
@@ -287,35 +497,104 @@
 	}
     }
 
-    if (n_nulls > 0 && default_nodataval) {
-	if (maptype == CELL_TYPE)
-	    G_important_message(_("Input raster map contains cells with NULL-value (no-data). "
-		       "The value %d was used to represent no-data values in the input map. "
-		       "You can specify a nodata value with the %s option."),
-		      (int)nodataval, nodatakey);
+    return ret;
+}
+
+int exact_range_check(double min, double max, GDALDataType datatype,
+		      const char *name)
+{
+
+    switch (datatype) {
+    case GDT_Byte:
+	if (min < TYPE_BYTE_MIN || max > TYPE_BYTE_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+		      GDALGetDataTypeName(datatype), TYPE_BYTE_MIN,
+		      TYPE_BYTE_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
+	}
 	else
-	    G_important_message(_("Input raster map contains cells with NULL-value (no-data). "
-		       "The value %g was used to represent no-data values in the input map. "
-		       "You can specify a nodata value with the %s option."),
-		      nodataval, nodatakey);
-    }
+	    return 0;
 
-    if (nodatavalmatch && n_nulls) {
-	if (default_nodataval) {  /* default nodataval didn't work */
-	    G_warning(_("The default nodata value is present in raster"
-	    "band <%s> and would lead to data loss. Please specify a "
-	    "different nodata value with the %s parameter."),
-	     name, nodatakey);
+    case GDT_UInt16:
+	if (min < TYPE_UINT16_MIN || max > TYPE_UINT16_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+		      GDALGetDataTypeName(datatype), TYPE_UINT16_MIN,
+		      TYPE_UINT16_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
 	}
-	else {  /* user-specified nodataval didn't work */
-	    G_warning(_("The given nodata value is present in raster"
-	    "band <%s> and would lead to data loss. Please specify a "
-	    "different nodata value with the %s parameter."),
-	     name, nodatakey);
+	else
+	    return 0;
+
+    case GDT_Int16:
+    case GDT_CInt16:
+	if (min < TYPE_INT16_MIN || max > TYPE_INT16_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+		      GDALGetDataTypeName(datatype), TYPE_INT16_MIN,
+		      TYPE_INT16_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
 	}
+	else
+	    return 0;
 
-	return -2;
+    case GDT_Int32:
+    case GDT_CInt32:
+	if (min < TYPE_INT32_MIN || max > TYPE_INT32_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+		      GDALGetDataTypeName(datatype), TYPE_INT32_MIN,
+		      TYPE_INT32_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
+	}
+	else
+	    return 0;
+
+    case GDT_UInt32:
+	if (min < TYPE_UINT32_MIN || max > TYPE_UINT32_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %u - %u"),
+		      GDALGetDataTypeName(datatype), TYPE_UINT32_MIN,
+		      TYPE_UINT32_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
+	}
+	else
+	    return 0;
+
+    case GDT_Float32:
+    case GDT_CFloat32:
+	if (min < TYPE_FLOAT32_MIN || max > TYPE_FLOAT32_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %f - %f"),
+		      GDALGetDataTypeName(datatype), TYPE_FLOAT32_MIN,
+		      TYPE_FLOAT32_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
+	}
+	else
+	    return 0;
+
+    case GDT_Float64:
+    case GDT_CFloat64:
+	/* not possible because DCELL is FLOAT64, not 128bit floating point, but anyway... */
+	if (min < TYPE_FLOAT64_MIN || max > TYPE_FLOAT64_MAX) {
+	    G_warning(_("Selected GDAL datatype does not cover data range."));
+	    G_warning(_("GDAL datatype: %s, range: %f - %f"),
+		      GDALGetDataTypeName(datatype), TYPE_FLOAT64_MIN,
+		      TYPE_FLOAT64_MAX);
+	    G_warning(_("Raster map <%s> range: %f - %f"), name, min, max);
+	    return 1;
+	}
+	else
+	    return 0;
+
+    default:
+	return 0;
     }
-
-    return 0;
 }

Modified: grass/branches/develbranch_6/raster/r.out.gdal/local_proto.h
===================================================================
--- grass/branches/develbranch_6/raster/r.out.gdal/local_proto.h	2009-12-23 12:49:43 UTC (rev 40119)
+++ grass/branches/develbranch_6/raster/r.out.gdal/local_proto.h	2009-12-23 12:51:14 UTC (rev 40120)
@@ -47,19 +47,18 @@
 #define TYPE_FLOAT64_MIN	(-MAXDOUBLE)
 #define TYPE_FLOAT64_MAX	MAXDOUBLE
 #else
-#define TYPE_FLOAT64_MIN	-1.79E308f
-#define TYPE_FLOAT64_MAX	1.79E308f
+#define TYPE_FLOAT64_MIN	-1.79E308
+#define TYPE_FLOAT64_MAX	1.79E308
 #endif
 
-#define GRASS_MAX_COLORS TYPE_UINT16_MAX
+#define GRASS_MAX_COLORS TYPE_UINT16_MAX    /* ok? */
 
-/* main.c */
-int range_check(double, double, GDALDataType, char *);
-int nullvalue_check(double, GDALDataType);
-
 /* export_band.c */
-int export_band(GDALDatasetH, int, const char *, const char *,
-		struct Cell_head *, RASTER_MAP_TYPE, double,
-		const char *, int, int);
+int export_band(GDALDatasetH, GDALDataType, int, const char *, 
+		const char *, struct Cell_head *, RASTER_MAP_TYPE, 
+		double, const char *, int, int);
+int exact_checks(GDALDataType, const char *, const char *,
+                 struct Cell_head *, RASTER_MAP_TYPE, double,
+		 const char *, int);
 
 #endif /* __LOCAL_PROTO_H__ */

Modified: grass/branches/develbranch_6/raster/r.out.gdal/main.c
===================================================================
--- grass/branches/develbranch_6/raster/r.out.gdal/main.c	2009-12-23 12:49:43 UTC (rev 40119)
+++ grass/branches/develbranch_6/raster/r.out.gdal/main.c	2009-12-23 12:51:14 UTC (rev 40120)
@@ -7,7 +7,7 @@
 *               based on GDAL library.
 *               Replaces r.out.gdal.sh script which used the gdal_translate
 *               executable and GDAL grass-format plugin.
-* COPYRIGHT:    (C) 2006-2008 by the GRASS Development Team
+* COPYRIGHT:    (C) 2006-2009 by the GRASS Development Team
 *
 *               This program is free software under the GNU General Public
 *   	    	License (>=v2). Read the file COPYING that comes with GRASS
@@ -32,13 +32,15 @@
 #include "cpl_string.h"
 #include "local_proto.h"
 
+int range_check(double, double, GDALDataType);
+int nodataval_check(double, GDALDataType);
+double set_default_nodata_value(GDALDataType, double, double);
 
-void supported_formats(char **formats)
+void supported_formats(const char **formats)
 {
     /* Code taken from r.in.gdal */
 
     int iDr;
-
     dbString gdal_formats;
 
     db_init_string(&gdal_formats);
@@ -110,13 +112,13 @@
 
     struct Cell_head cellhead;
     struct Ref ref;
-    char *mapset, *gdal_formats = NULL;
-    RASTER_MAP_TYPE maptype;
+    const char *mapset, *gdal_formats = NULL;
+    RASTER_MAP_TYPE maptype, testmaptype;
     int bHaveMinMax;
-    double dfCellMin;
-    double dfCellMax;
+    double dfCellMin, export_min;
+    double dfCellMax, export_max;
     struct FPRange sRange;
-    int retval;
+    int retval = 0, check_range;
 
     G_gisinit(argv[0]);
 
@@ -138,7 +140,7 @@
     flag_f = G_define_flag();
     flag_f->key = 'f';
     flag_f->label = _("Force raster export also if data loss may occur");
-    flag_f->description = _("Overrides saftey checks.");
+    flag_f->description = _("Overrides nodata saftey check.");
 
     input = G_define_standard_option(G_OPT_R_INPUT);
     input->required = NO;
@@ -275,9 +277,8 @@
 			  format->answer);
     }
 
-    /* Determine raster data type */
+    /* Determine GDAL data type */
     GDALDataType datatype = GDT_Unknown;
-    double nodataval;
 
     maptype = CELL_TYPE;
 
@@ -286,141 +287,236 @@
 	if (type->answer[0] == 'B') {
 	    datatype = GDT_Byte;
 	    maptype = CELL_TYPE;
-	    nodataval = TYPE_BYTE_MIN;
 	}
 	else if (type->answer[0] == 'I') {
 	    if (strcmp(type->answer, "Int16") == 0) {
 		datatype = GDT_Int16;
 		maptype = CELL_TYPE;
-		nodataval = TYPE_INT16_MIN;
 	    }
 	    else if (strcmp(type->answer, "Int32") == 0) {
 		datatype = GDT_Int32;
 		maptype = CELL_TYPE;
-		nodataval = TYPE_INT32_MIN;
 	    }
 	}
 	else if (type->answer[0] == 'U') {
 	    if (strcmp(type->answer, "UInt16") == 0) {
 		datatype = GDT_UInt16;
 		maptype = CELL_TYPE;
-		nodataval = TYPE_UINT16_MIN;
 	    }
 	    else if (strcmp(type->answer, "UInt32") == 0) {
 		datatype = GDT_UInt32;
 		maptype = DCELL_TYPE;
-		nodataval = TYPE_UINT32_MIN;
 	    }
 	}
 	else if (type->answer[0] == 'F') {
 	    if (strcmp(type->answer, "Float32") == 0) {
 		datatype = GDT_Float32;
 		maptype = FCELL_TYPE;
-		/* nodataval = TYPE_FLOAT32_MIN; */
-		nodataval = 0.0/0.0;
 	    }
 	    else if (strcmp(type->answer, "Float64") == 0) {
 		datatype = GDT_Float64;
 		maptype = DCELL_TYPE;
-		/* nodataval = TYPE_FLOAT64_MIN; */
-		nodataval = 0.0/0.0;
 	    }
 	}
 	else if (type->answer[0] == 'C') {
 	    if (strcmp(type->answer, "CInt16") == 0) {
 		datatype = GDT_CInt16;
 		maptype = CELL_TYPE;
-		nodataval = TYPE_INT16_MIN;
 	    }
 	    else if (strcmp(type->answer, "CInt32") == 0) {
 		datatype = GDT_CInt32;
 		maptype = CELL_TYPE;
-		nodataval = TYPE_INT32_MIN;
 	    }
 	    else if (strcmp(type->answer, "CFloat32") == 0) {
 		datatype = GDT_CFloat32;
 		maptype = FCELL_TYPE;
-		/* nodataval = TYPE_FLOAT32_MIN; */
-		nodataval = 0.0/0.0;
 	    }
 	    else if (strcmp(type->answer, "CFloat64") == 0) {
 		datatype = GDT_CFloat64;
 		maptype = DCELL_TYPE;
-		/* nodataval = TYPE_FLOAT64_MIN; */
-		nodataval = 0.0/0.0;
 	    }
 	}
     }
 
-    /* If file type not set by user ... */
-    /* Get min/max values. */
-    if (G_read_fp_range(ref.file[0].name, ref.file[0].mapset, &sRange) == -1) {
-	bHaveMinMax = FALSE;
+    /* get min/max values */
+    int band;
+
+    check_range = 0;
+    bHaveMinMax = TRUE;
+    export_min = TYPE_FLOAT64_MIN;
+    export_max = TYPE_FLOAT64_MAX;
+    for (band = 0; band < ref.nfiles; band++) {
+	if (G_read_fp_range
+	    (ref.file[band].name, ref.file[band].mapset, &sRange) == -1) {
+	    bHaveMinMax = FALSE;
+	    G_warning(_("Could not read data range of raster <%s>"),
+		      ref.file[band].name);
+	}
+	else {
+	    G_get_fp_range_min_max(&sRange, &dfCellMin, &dfCellMax);
+	    if (band == 0) {
+		export_min = dfCellMin;
+		export_max = dfCellMax;
+	    }
+	    else {
+		if (export_min > dfCellMin)
+		    export_min = dfCellMin;
+		if (export_max < dfCellMax)
+		    export_max = dfCellMax;
+	    }
+	}
+	G_debug(3, "Range of <%s>: min: %f, max: %f", ref.file[band].name,
+		dfCellMin, dfCellMax);
     }
-    else {
-	bHaveMinMax = TRUE;
-	G_get_fp_range_min_max(&sRange, &dfCellMin, &dfCellMax);
+    if (bHaveMinMax == FALSE) {
+	export_min = TYPE_FLOAT64_MIN;
+	export_max = TYPE_FLOAT64_MAX;
     }
-    G_debug(3, "Range: min: %f, max: %f", dfCellMin, dfCellMax);
+    G_debug(3, "Total range: min: %f, max: %f", export_min, export_max);
 
+    /* GDAL datatype not set by user, determine suitable datatype */
     if (datatype == GDT_Unknown) {
-	/* ... determine raster data type from first GRASS raster in a group */
+	/* Use raster data type from first GRASS raster in a group */
 	maptype = G_raster_map_type(ref.file[0].name, ref.file[0].mapset);
 	if (maptype == FCELL_TYPE) {
 	    datatype = GDT_Float32;
-	    /* nodataval = TYPE_FLOAT32_MIN; */
-	    nodataval = 0.0/0.0;
 	}
 	else if (maptype == DCELL_TYPE) {
 	    datatype = GDT_Float64;
-	    /* nodataval = TYPE_FLOAT64_MIN; */
-	    nodataval = 0.0/0.0;
 	}
 	else {
 	    /* Special tricks for GeoTIFF color table support and such */
-	    if (dfCellMin >= 0 && dfCellMax < 256) {
+	    if (export_min >= TYPE_BYTE_MIN && export_max <= TYPE_BYTE_MAX) {
 		datatype = GDT_Byte;
-		nodataval = TYPE_BYTE_MIN;
 	    }
 	    else {
-		if (dfCellMin >= 0 && dfCellMax < 65536) {
+		if (export_min >= TYPE_UINT16_MIN &&
+		    export_max <= TYPE_UINT16_MAX) {
 		    datatype = GDT_UInt16;
-		    nodataval = TYPE_UINT16_MIN;
 		}
+		else if (export_min >= TYPE_INT16_MIN &&
+			 export_max <= TYPE_INT16_MAX) {
+		    datatype = GDT_Int16;
+		}
 		else {
 		    datatype = GDT_Int32;	/* need to fine tune this more? */
-		    nodataval = TYPE_INT32_MIN;
 		}
 	    }
 	}
     }
 
-    /* default or user-specified nodata-value ? */
+    /* got a GDAL datatype, report to user */
+    G_message(_("Exporting to GDAL data type: %s"),
+	      GDALGetDataTypeName(datatype));
+
+    G_debug(3, "Input map datatype=%s\n",
+	    (maptype == CELL_TYPE ? "CELL" :
+	     (maptype == DCELL_TYPE ? "DCELL" :
+	      (maptype == FCELL_TYPE ? "FCELL" : "??"))));
+
+
+    /* if GDAL datatype set by user, do checks */
+    if (type->answer) {
+
+	/* Check if raster data range is outside of the range of 
+	 * given GDAL datatype, not even overlapping */
+	if (range_check(export_min, export_max, datatype))
+	    G_fatal_error(_("Raster export would result in complete data loss, aborting."));
+
+	/* Precision tests */
+	for (band = 0; band < ref.nfiles; band++) {
+	    testmaptype =
+		G_raster_map_type(ref.file[band].name, ref.file[band].mapset);
+	    /* Exporting floating point rasters to some integer type ? */
+	    if ((testmaptype == FCELL_TYPE || testmaptype == DCELL_TYPE) &&
+		(datatype == GDT_Byte || datatype == GDT_Int16 ||
+		 datatype == GDT_UInt16 || datatype == GDT_Int32 ||
+		 datatype == GDT_UInt32)) {
+		G_warning(_("Precision loss: Raster map <%s> of type %s to be exported as %s. "
+			   "This can be avoided by using %s."),
+			  ref.file[band].name,
+			  (maptype == FCELL_TYPE ? "FCELL" : "DCELL"),
+			  GDALGetDataTypeName(datatype),
+			  (maptype == FCELL_TYPE ? "Float32" : "Float64"));
+		retval = -1;
+	    }
+	    /* Exporting CELL with large values to GDT_Float32 ? Cap is 2^24 */
+	    if (testmaptype == CELL_TYPE && datatype == GDT_Float32 &&
+		(dfCellMin < -16777216 || dfCellMax > 16777216)) {
+		G_warning(_("Precision loss: The range of <%s> can not be "
+			    "accurately preserved with GDAL datatype Float32. "
+			    "This can be avoided by exporting to Int32 or Float64."),
+			  ref.file[band].name);
+		retval = -1;
+	    }
+	    /* Exporting DCELL to GDT_Float32 ? */
+	    if (testmaptype == DCELL_TYPE && datatype == GDT_Float32) {
+		G_warning(_("Precision loss: Float32 can not preserve the "
+			    "DCELL precision of raster <%s>. "
+			    "This can be avoided by using Float64"),
+			  ref.file[band].name);
+		retval = -1;
+	    }
+	}
+	if (retval == -1) {
+	    if (flag_f->answer)
+		G_warning(_("Forcing raster export."));
+	    else
+		G_fatal_error(_("Raster export aborted."));
+	}
+    }
+
+    /* Nodata value */
+    double nodataval;
     int default_nodataval = 1;
 
+    /* User-specified nodata-value ? */
     if (nodataopt->answer) {
 	nodataval = atof(nodataopt->answer);
 	default_nodataval = 0;
-	/* check nodataval here, not in export_band(), because it is specified only once */
-	if (nullvalue_check(nodataval, datatype)) {
+	/* Check if given nodata value can be represented by selected GDAL datatype */
+	if (nodataval_check(nodataval, datatype)) {
+	    G_fatal_error(_("Raster export aborted."));
+	}
+    }
+    /* Set reasonable default nodata value */
+    else {
+	nodataval =
+	    set_default_nodata_value(datatype, export_min, export_max);
+    }
+
+    /* exact range and nodata checks for each band */
+    G_message(_("Checking GDAl data type and nodata value"));
+    for (band = 0; band < ref.nfiles; band++) {
+	if (ref.nfiles > 1) {
+	    G_verbose_message(_("Checking options for raster map <%s> (band %d)..."),
+			      G_fully_qualified_name(ref.file[band].name,
+						     ref.file[band].mapset),
+			      band + 1);
+	}
+
+	retval = exact_checks
+	    (datatype, ref.file[band].name, ref.file[band].mapset,
+	     &cellhead, maptype, nodataval, nodataopt->key,
+	     default_nodataval);
+
+	/* nodata problem */
+	if (retval == -1) {
 	    if (flag_f->answer)
 		G_warning(_("Forcing raster export."));
 	    else
 		G_fatal_error(_("Raster export aborted."));
 	}
+	/* data don't fit into range of GDAL datatype */
+	else if (retval == -2) {
+	    G_fatal_error(_("Raster export aborted."));
+	}
     }
 
-    G_debug(3, "Input map datatype=%s\n",
-	    (maptype == CELL_TYPE ? "CELL" :
-	     (maptype == DCELL_TYPE ? "DCELL" :
-	      (maptype == FCELL_TYPE ? "FCELL" : "??"))));
-    G_message(_("Exporting to GDAL data type: %s"),
-	      GDALGetDataTypeName(datatype));
-
     /* Create dataset for output with target driver or, if needed, with in-memory driver */
     char **papszOptions = NULL;
 
-    /* parse dataset creation options */
+    /* Parse dataset creation options */
     if (createopt->answer) {
 	int i;
 
@@ -478,8 +574,7 @@
     AttachMetadata(hCurrDS, metaopt->answers);
 
     /* Export to GDAL raster */
-    int band;
-
+    G_message(_("Exporting to GDAL raster"));
     for (band = 0; band < ref.nfiles; band++) {
 	if (ref.nfiles > 1) {
 	    G_verbose_message(_("Exporting raster map <%s> (band %d)..."),
@@ -487,44 +582,21 @@
 						     ref.file[band].mapset),
 			      band + 1);
 	}
-	/* do range check here because we don't have GDALDataType in export_band() */
-	if (G_read_fp_range
-	    (ref.file[band].name, ref.file[band].mapset, &sRange) == -1) {
-	    bHaveMinMax = FALSE;
-	}
-	else {
-	    bHaveMinMax = TRUE;
-	    G_get_fp_range_min_max(&sRange, &dfCellMin, &dfCellMax);
-	}
-	G_debug(3, "Range: min: %f, max: %f", dfCellMin, dfCellMax);
-	if (bHaveMinMax == TRUE) {
-	    if (range_check
-		(dfCellMin, dfCellMax, datatype, ref.file[band].name)) {
-		if (flag_f->answer)
-		    G_warning(_("Forcing raster export."));
-		else
-		    G_fatal_error(_("Raster export aborted."));
-		}
-	}
 
-	/* ready to export */
 	retval = export_band
-	    (hCurrDS, band + 1, ref.file[band].name, ref.file[band].mapset,
-	     &cellhead, maptype, nodataval, nodataopt->key, flag_c->answer,
-	     default_nodataval);
+	    (hCurrDS, datatype, band + 1, ref.file[band].name,
+	     ref.file[band].mapset, &cellhead, maptype, nodataval,
+	     nodataopt->key, flag_c->answer, default_nodataval);
+
+	/* read/write error */
 	if (retval == -1) {
 	    G_warning(_("Unable to export raster map <%s>"),
 		      ref.file[band].name);
 	}
-	else if (retval == -2) {
-	    if (flag_f->answer)
-		G_warning(_("Forcing raster export."));
-	    else
-		G_fatal_error(_("Raster export aborted."));
-	}
     }
 
-    /* Finaly create user required raster format from memory raster if in-memory driver was used */
+    /* Finaly create user requested raster format from memory raster 
+     * if in-memory driver was used */
     if (hMEMDS) {
 	hDstDS =
 	    GDALCreateCopy(hDriver, output->answer, hMEMDS, FALSE,
@@ -546,30 +618,30 @@
 }
 
 
-int range_check(double min, double max, GDALDataType datatype, char *name)
+int range_check(double min, double max, GDALDataType datatype)
 {
     /* what accuracy to use to print min max for FLOAT32 and FLOAT64? %g enough? */
 
     switch (datatype) {
     case GDT_Byte:
-	if (min < TYPE_BYTE_MIN || max > TYPE_BYTE_MAX) {
+	if (max < TYPE_BYTE_MIN || min > TYPE_BYTE_MAX) {
 	    G_warning(_("Selected GDAL datatype does not cover data range."));
 	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
 		      GDALGetDataTypeName(datatype), TYPE_BYTE_MIN,
 		      TYPE_BYTE_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
+	    G_warning(_("Range to be exported: %f - %f"), min, max);
 	    return 1;
 	}
 	else
 	    return 0;
 
     case GDT_UInt16:
-	if (min < TYPE_UINT16_MIN || max > TYPE_UINT16_MAX) {
+	if (max < TYPE_UINT16_MIN || min > TYPE_UINT16_MAX) {
 	    G_warning(_("Selected GDAL datatype does not cover data range."));
 	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
 		      GDALGetDataTypeName(datatype), TYPE_UINT16_MIN,
 		      TYPE_UINT16_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
+	    G_warning(_("Range to be exported: %f - %f"), min, max);
 	    return 1;
 	}
 	else
@@ -577,12 +649,12 @@
 
     case GDT_Int16:
     case GDT_CInt16:
-	if (min < TYPE_INT16_MIN || max > TYPE_INT16_MAX) {
+	if (max < TYPE_INT16_MIN || min > TYPE_INT16_MAX) {
 	    G_warning(_("Selected GDAL datatype does not cover data range."));
 	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
 		      GDALGetDataTypeName(datatype), TYPE_INT16_MIN,
 		      TYPE_INT16_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
+	    G_warning(_("Range to be exported: %f - %f"), min, max);
 	    return 1;
 	}
 	else
@@ -590,24 +662,24 @@
 
     case GDT_Int32:
     case GDT_CInt32:
-	if (min < TYPE_INT32_MIN || max > TYPE_INT32_MAX) {
+	if (max < TYPE_INT32_MIN || min > TYPE_INT32_MAX) {
 	    G_warning(_("Selected GDAL datatype does not cover data range."));
 	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
 		      GDALGetDataTypeName(datatype), TYPE_INT32_MIN,
 		      TYPE_INT32_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
+	    G_warning(_("Range to be exported: %f - %f"), min, max);
 	    return 1;
 	}
 	else
 	    return 0;
 
     case GDT_UInt32:
-	if (min < TYPE_UINT32_MIN || max > TYPE_UINT32_MAX) {
+	if (max < TYPE_UINT32_MIN || min > TYPE_UINT32_MAX) {
 	    G_warning(_("Selected GDAL datatype does not cover data range."));
 	    G_warning(_("GDAL datatype: %s, range: %u - %u"),
 		      GDALGetDataTypeName(datatype), TYPE_UINT32_MIN,
 		      TYPE_UINT32_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
+	    G_warning(_("Range to be exported: %f - %f"), min, max);
 	    return 1;
 	}
 	else
@@ -615,12 +687,12 @@
 
     case GDT_Float32:
     case GDT_CFloat32:
-	if (min < TYPE_FLOAT32_MIN || max > TYPE_FLOAT32_MAX) {
+	if (max < TYPE_FLOAT32_MIN || min > TYPE_FLOAT32_MAX) {
 	    G_warning(_("Selected GDAL datatype does not cover data range."));
 	    G_warning(_("GDAL datatype: %s, range: %g - %g"),
 		      GDALGetDataTypeName(datatype), TYPE_FLOAT32_MIN,
 		      TYPE_FLOAT32_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
+	    G_warning(_("Range to be exported: %f - %f"), min, max);
 	    return 1;
 	}
 	else
@@ -628,33 +700,26 @@
 
     case GDT_Float64:
     case GDT_CFloat64:
-	/* not possible because DCELL is FLOAT64, not 128bit floating point, but anyway... */
-	if (min < TYPE_FLOAT64_MIN || max > TYPE_FLOAT64_MAX) {
-	    G_warning(_("Selected GDAL datatype does not cover data range."));
-	    G_warning(_("GDAL datatype: %s, range: %g - %g"),
-		      GDALGetDataTypeName(datatype), TYPE_FLOAT64_MIN,
-		      TYPE_FLOAT64_MAX);
-	    G_warning(_("Raster map <%s> range: %g - %g"), name, min, max);
-	    return 1;
-	}
-	else
-	    return 0;
+	/* not needed because FLOAT64 should always cover the data range */
+	return 0;
 
     default:
 	return 0;
     }
 }
 
-int nullvalue_check(double nodataval, GDALDataType datatype)
+int nodataval_check(double nodataval, GDALDataType datatype)
 {
-    /* what accuracy to use to print nodataval for FLOAT32 and FLOAT64? %g enough? */
 
     switch (datatype) {
     case GDT_Byte:
-	if (nodataval < TYPE_BYTE_MIN || nodataval > TYPE_BYTE_MAX) {
-	    G_warning(_("Specified nodata value %d is not covered by range of selected GDAL datatype."),
-		      (int) nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+    	/* the additional cast to CELL is what happens in export_band()
+	 * accordingly below for the other GDAL types */
+	if (nodataval != (double)(GByte)(CELL) nodataval) {
+	    G_warning(_("Mismatch between metadata nodata value and actual nodata value in exported raster: "
+		       "specified nodata value %f gets converted to %d by selected GDAL datatype."),
+		      nodataval, (GByte)(CELL) nodataval);
+	    G_warning(_("GDAL datatype: %s, valid range: %d - %d"),
 		      GDALGetDataTypeName(datatype), TYPE_BYTE_MIN,
 		      TYPE_BYTE_MAX);
 	    return 1;
@@ -663,10 +728,11 @@
 	    return 0;
 
     case GDT_UInt16:
-	if (nodataval < TYPE_UINT16_MIN || nodataval > TYPE_UINT16_MAX) {
-	    G_warning(_("Specified nodata value %d is not covered by range of selected GDAL datatype."),
-		      (int) nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+	if (nodataval != (double)(GUInt16)(CELL) nodataval) {
+	    G_warning(_("Mismatch between metadata nodata value and actual nodata value in exported raster: "
+		       "specified nodata value %f gets converted to %d by selected GDAL datatype."),
+		      nodataval, (GUInt16)(CELL) nodataval);
+	    G_warning(_("GDAL datatype: %s, valid range: %u - %u"),
 		      GDALGetDataTypeName(datatype), TYPE_UINT16_MIN,
 		      TYPE_UINT16_MAX);
 	    return 1;
@@ -676,10 +742,11 @@
 
     case GDT_Int16:
     case GDT_CInt16:
-	if (nodataval < TYPE_INT16_MIN || nodataval > TYPE_INT16_MAX) {
-	    G_warning(_("Specified nodata value %d is not covered by range of selected GDAL datatype."),
-		      (int) nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
+	if (nodataval != (double)(GInt16)(CELL) nodataval) {
+	    G_warning(_("Mismatch between metadata nodata value and actual nodata value in exported raster: "
+		       "specified nodata value %f gets converted to %d by selected GDAL datatype."),
+		      nodataval, (GInt16)(CELL) nodataval);
+	    G_warning(_("GDAL datatype: %s, valid range: %d - %d"),
 		      GDALGetDataTypeName(datatype), TYPE_INT16_MIN,
 		      TYPE_INT16_MAX);
 	    return 1;
@@ -687,26 +754,29 @@
 	else
 	    return 0;
 
-    case GDT_Int32:
-    case GDT_CInt32:
-	if (nodataval < TYPE_INT32_MIN || nodataval > TYPE_INT32_MAX) {
-	    G_warning(_("Specified nodata value %lld is not covered by range of selected GDAL datatype."),
-		      (long long) nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %d - %d"),
-		      GDALGetDataTypeName(datatype), TYPE_INT32_MIN,
-		      TYPE_INT32_MAX);
+    case GDT_UInt32:
+	if (nodataval != (double)(GUInt32)(DCELL) nodataval) {
+	    G_warning(_("Mismatch between metadata nodata value and actual nodata value in exported raster: "
+		       "specified nodata value %f gets converted to %d by selected GDAL datatype."),
+		      nodataval, (GUInt32)(DCELL) nodataval);
+	    G_warning(_("GDAL datatype: %s, valid range: %u - %u"),
+		      GDALGetDataTypeName(datatype), TYPE_UINT32_MIN,
+		      TYPE_UINT32_MAX);
 	    return 1;
 	}
 	else
 	    return 0;
 
-    case GDT_UInt32:
-	if (nodataval < TYPE_UINT32_MIN || nodataval > TYPE_UINT32_MAX) {
-	    G_warning(_("Specified nodata value %lld is not covered by range of selected GDAL datatype."),
-		      (long long) nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %u - %u"),
-		      GDALGetDataTypeName(datatype), TYPE_UINT32_MIN,
-		      TYPE_UINT32_MAX);
+    case GDT_Int32:
+    case GDT_CInt32:
+    	/* GInt32 is equal to CELL, but that may change in the future */
+	if (nodataval != (double)(GInt32)(CELL) nodataval) {
+	    G_warning(_("Mismatch between metadata nodata value and actual nodata value in exported raster: "
+		       "specified nodata value %f gets converted to %d by selected GDAL datatype."),
+		      nodataval, (GInt32)(CELL) nodataval);
+	    G_warning(_("GDAL datatype: %s, valid range: %d - %d"),
+		      GDALGetDataTypeName(datatype), TYPE_INT32_MIN,
+		      TYPE_INT32_MAX);
 	    return 1;
 	}
 	else
@@ -714,10 +784,11 @@
 
     case GDT_Float32:
     case GDT_CFloat32:
-	if (nodataval < TYPE_FLOAT32_MIN || nodataval > TYPE_FLOAT32_MAX) {
-	    G_warning(_("Specified nodata value %g is not covered by range of selected GDAL datatype."),
-		      nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %g - %g"),
+	if (nodataval != (double)(float) nodataval) {
+	    G_warning(_("Mismatch between metadata nodata value and actual nodata value in exported raster: "
+		       "specified nodata value %f gets converted to %f by selected GDAL datatype."),
+		      nodataval, (float) nodataval);
+	    G_warning(_("GDAL datatype: %s, valid range: %g - %g"),
 		      GDALGetDataTypeName(datatype), TYPE_FLOAT32_MIN,
 		      TYPE_FLOAT32_MAX);
 	    return 1;
@@ -727,18 +798,67 @@
 
     case GDT_Float64:
     case GDT_CFloat64:
-	/* not possible because double is FLOAT64, not 128bit floating point */
-	if (nodataval < TYPE_FLOAT64_MIN || nodataval > TYPE_FLOAT64_MAX) {
-	    G_warning(_("Specified nodata value %g is not covered by range of selected GDAL datatype."),
-		      nodataval);
-	    G_warning(_("GDAL datatype: %s, range: %g - %g"),
-		      GDALGetDataTypeName(datatype), TYPE_FLOAT64_MIN,
-		      TYPE_FLOAT64_MAX);
-	    return 1;
-	}
+	/* not needed because FLOAT64 is equal to double */
+	return 0;
+
+    default:
+	return 0;
+    }
+}
+
+double set_default_nodata_value(GDALDataType datatype, double min, double max)
+{
+    switch (datatype) {
+    case GDT_Byte:
+	if (max < TYPE_BYTE_MAX)
+	    return (double)TYPE_BYTE_MAX;
+	else if (min > TYPE_BYTE_MIN)
+	    return (double)TYPE_BYTE_MIN;
 	else
-	    return 0;
+	    return (double)TYPE_BYTE_MAX;
 
+    case GDT_UInt16:
+	if (max < TYPE_UINT16_MAX)
+	    return (double)TYPE_UINT16_MAX;
+	else if (min > TYPE_UINT16_MIN)
+	    return (double)TYPE_UINT16_MIN;
+	else
+	    return (double)TYPE_UINT16_MAX;
+
+    case GDT_Int16:
+    case GDT_CInt16:
+	if (min > TYPE_INT16_MIN)
+	    return (double)TYPE_INT16_MIN;
+	else if (max < TYPE_INT16_MAX)
+	    return (double)TYPE_INT16_MAX;
+	else
+	    return (double)TYPE_INT16_MIN;
+
+    case GDT_UInt32:
+	if (max < TYPE_UINT32_MAX)
+	    return (double)TYPE_UINT32_MAX;
+	else if (min > TYPE_UINT32_MIN)
+	    return (double)TYPE_UINT32_MIN;
+	else
+	    return (double)TYPE_UINT32_MAX;
+
+    case GDT_Int32:
+    case GDT_CInt32:
+	if (min > TYPE_INT32_MIN)
+	    return (double)TYPE_INT32_MIN;
+	else if (max < TYPE_INT32_MAX)
+	    return (double)TYPE_INT32_MAX;
+	else
+	    return (double)TYPE_INT32_MIN;
+
+    case GDT_Float32:
+    case GDT_CFloat32:
+	return 0.0 / 0.0;
+
+    case GDT_Float64:
+    case GDT_CFloat64:
+	return 0.0 / 0.0;
+
     default:
 	return 0;
     }



More information about the grass-commit mailing list