diff --git a/autotest/gdrivers/aaigrid.py b/autotest/gdrivers/aaigrid.py index 0025b288ca31..551d1f8cc077 100755 --- a/autotest/gdrivers/aaigrid.py +++ b/autotest/gdrivers/aaigrid.py @@ -250,7 +250,7 @@ def test_aaigrid_10(): except OSError: pass - + ############################################################################### # Test SIGNIFICANT_DIGITS creation option (same as DECIMAL_PRECISION test) @@ -351,7 +351,52 @@ def test_aaigrid_15(): gdal.Unlink('/vsimem/aaigrid_15.asc') + +############################################################################### +# Test support for D12 AAIGRID with null value (#5095) + + +def test_aaigrid_null(): + + gdal.FileFromMemBuffer('/vsimem/test_aaigrid_null.asc', """ncols 4 +nrows 1 +xllcorner 0 +yllcorner -1 +cellsize 1 +NODATA_value null +null 1.5 null 3.5 +""") + + ds = gdal.Open('/vsimem/test_aaigrid_null.asc') + assert ds.GetRasterBand(1).DataType == gdal.GDT_Float32 + assert ds.GetRasterBand(1).GetNoDataValue() < -1e38 + assert ds.GetRasterBand(1).ComputeRasterMinMax() == (1.5,3.5) + ds = None + + gdal.Unlink('/vsimem/test_aaigrid_null.asc') + + ############################################################################### +# Test support for D12 AAIGRID with null value and force AAIGRID_DATATYPE=Float64 (#5095) + +def test_aaigrid_null_float64(): + + gdal.FileFromMemBuffer('/vsimem/test_aaigrid_null.asc', """ncols 4 +nrows 1 +xllcorner 0 +yllcorner -1 +cellsize 1 +NODATA_value null +null 1.5 null 3.5 +""") + + with gdaltest.config_option('AAIGRID_DATATYPE', 'Float64'): + ds = gdal.Open('/vsimem/test_aaigrid_null.asc') + assert ds.GetRasterBand(1).DataType == gdal.GDT_Float64 + assert ds.GetRasterBand(1).GetNoDataValue() < -1e308 + assert ds.GetRasterBand(1).ComputeRasterMinMax() == (1.5,3.5) + ds = None + gdal.Unlink('/vsimem/test_aaigrid_null.asc') diff --git a/frmts/aaigrid/aaigriddataset.cpp b/frmts/aaigrid/aaigriddataset.cpp index eb33626938e8..ec719d1a354b 100644 --- a/frmts/aaigrid/aaigriddataset.cpp +++ b/frmts/aaigrid/aaigriddataset.cpp @@ -194,11 +194,25 @@ CPLErr AAIGRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff, if( pImage != nullptr ) { + // "null" seems to be specific of D12 software + // See https://github.com/OSGeo/gdal/issues/5095 if( eDataType == GDT_Float64 ) - reinterpret_cast(pImage)[iPixel] = CPLAtofM(szToken); + { + if( strcmp(szToken, "null") == 0 ) + reinterpret_cast(pImage)[iPixel] = + -std::numeric_limits::max(); + else + reinterpret_cast(pImage)[iPixel] = CPLAtofM(szToken); + } else if( eDataType == GDT_Float32 ) - reinterpret_cast(pImage)[iPixel] = - DoubleToFloatClamp(CPLAtofM(szToken)); + { + if( strcmp(szToken, "null") == 0 ) + reinterpret_cast(pImage)[iPixel] = + -std::numeric_limits::max(); + else + reinterpret_cast(pImage)[iPixel] = + DoubleToFloatClamp(CPLAtofM(szToken)); + } else reinterpret_cast(pImage)[iPixel] = static_cast(atoi(szToken)); @@ -553,24 +567,42 @@ int AAIGDataset::ParseHeader(const char *pszHeader, const char *pszDataType) const char *pszNoData = papszTokens[i + 1]; bNoDataSet = true; - dfNoDataValue = CPLAtofM(pszNoData); - if( pszDataType == nullptr && - (strchr(pszNoData, '.') != nullptr || - strchr(pszNoData, ',') != nullptr || - std::numeric_limits::min() > dfNoDataValue || - dfNoDataValue > std::numeric_limits::max()) ) + if( strcmp(pszNoData, "null") == 0 ) { - eDataType = GDT_Float32; - if( !CPLIsInf(dfNoDataValue) && - (fabs(dfNoDataValue) < std::numeric_limits::min() || - fabs(dfNoDataValue) > std::numeric_limits::max()) ) + // "null" seems to be specific of D12 software + // See https://github.com/OSGeo/gdal/issues/5095 + if( pszDataType == nullptr || eDataType == GDT_Float32 ) { + dfNoDataValue = -std::numeric_limits::max(); + eDataType = GDT_Float32; + } + else + { + dfNoDataValue = -std::numeric_limits::max(); eDataType = GDT_Float64; } } - if( eDataType == GDT_Float32 ) + else { - dfNoDataValue = MapNoDataToFloat(dfNoDataValue); + dfNoDataValue = CPLAtofM(pszNoData); + if( pszDataType == nullptr && + (strchr(pszNoData, '.') != nullptr || + strchr(pszNoData, ',') != nullptr || + std::numeric_limits::min() > dfNoDataValue || + dfNoDataValue > std::numeric_limits::max()) ) + { + eDataType = GDT_Float32; + if( !CPLIsInf(dfNoDataValue) && + (fabs(dfNoDataValue) < std::numeric_limits::min() || + fabs(dfNoDataValue) > std::numeric_limits::max()) ) + { + eDataType = GDT_Float64; + } + } + if( eDataType == GDT_Float32 ) + { + dfNoDataValue = MapNoDataToFloat(dfNoDataValue); + } } } @@ -1000,7 +1032,11 @@ GDALDataset *AAIGDataset::CommonOpen( GDALOpenInfo *poOpenInfo, poOpenInfo->pabyHeader[i - 1] == '\r' || poOpenInfo->pabyHeader[i - 2] == '\r' ) { - if( !isalpha(poOpenInfo->pabyHeader[i]) && + if( (!isalpha(poOpenInfo->pabyHeader[i]) || + // null seems to be specific of D12 software + // See https://github.com/OSGeo/gdal/issues/5095 + (i + 5 < poOpenInfo->nHeaderBytes && + memcmp(poOpenInfo->pabyHeader + i, "null ", 5) == 0)) && poOpenInfo->pabyHeader[i] != '\n' && poOpenInfo->pabyHeader[i] != '\r') {