[PATCH] Fix handling of NODATA values (for DTED files)

[tbonfort]

Reporter: rouault Date: 2007/07/26 - 23:43 Currently, when doing WCS on DTED files (whose NODATA value is -32767), you may experience problems when using INTERPOLATION=AVERAGE or INTERPOLATION=BILINEAR. That's because the NODATA value is not recognized as such and the resampling algorithms are using -32767 as a normal value, leading to 'fun' values like -16145 when computing averages. My proposal is to force the noData value as offsite.red when the source is a one band image.

The following addition in msResampleGDALToMap fixes this :

    if (srcImage->format->bands == 1)
    {
      void* hBand1 = GDALGetRasterBand( hDS, 1 );
      noData = msGetGDALNoDataValue( layer, hBand1, &bGotNoData );
      if (bGotNoData)
      {
        layer->offsite.red = noData;
      }
    }

But, this fix is not really sufficient to address the correct handling of DTED files. In aeraunatical applications using DTED files, the meaning of the NODATA value is really 'unknown', 'there's maybe a big mountain here but no data collected'. So, you often prefer the propagation of the uncertainty of the source to the result, rather than a clever interpolation that may be dangerous for the final user. So, if for some processing of an application, you need to resample a DTED Level 2 file to a lower resolution, like near the resolution of a DTED1 file, you want -32767 to be the result of the resampling as soon as one of source grid elements is -32767. With the small patch above, those source grid elements are just ignored for the computation.[[BR]]

So, I've introduced a new PROCESSING directive ('STRICT_NODATA=YES') that enforces the propagation of the noData value. Of course, when using STRICT_NODATA=yes, you should also set 'LOAD_FULL_RES_IMAGE=YES' to be sure that all the source grid elements will be read. And of course, the default resampling method, NEAREST, doesn't make much sense... Debug warnings are emitted when those 2 conditions are not met.

PS : I've a longer patch that adds another INTERPOLATION method, 'MAXIMUM' useful for aeronautical applications too (safety altitude computations for example). As its name suggests, it takes the maximum value of the source grid elements (except a noData value ;-)). As it depends of this current patch, I've thought that I should delay its posting.

PS2 : It doesn't include the small patch I've previously sent for the reprojection/tileindex issue (http://trac.osgeo.org/mapserver/ticket/2180)

[tbonfort]

Author: sdlime Date: 2007/08/01 - 05:35 Frank: What do you think? Can we fix in 5.0 or wait...

Steve

[tbonfort]

Author: warmerdam Date: 2007/08/01 - 16:07 Steve,

My concern with this patch is that it incorporates some very specific kinds of processing for particular applications. It feels like unwarranted complication. So I'm ambivalent.

I think picking up the nodata value from GDAL for the resampling code might be reasonable, but I wonder:

• Shouldn't we be handling offsite/nodata stuff in mapdrawgdal.c instead? It would be disturbing (I think) if nodata got used as the offsite value only if the map gets resampled as opposed to the normal up/down sampling logic.
• Doesn't this code wipe out a user supplied offsite value?

I think a lot of this bug report is new functionality and shouldn't be added for 5.0 as we are past the new functionality cutoff. There may be a portion that qualifies as a bug fix, but I'm not all that comfortable with it yet.

As this deeply affects the general raster rendering engine, I'd rather take this bug if you don't mind. But I'm not promising to act on it for 5.0.

[tbonfort]

Author: sdlime Date: 2007/08/01 - 17:04 I figured it belongs with you anyway to will reasign. I'm ok with waiting until 5.2...

Steve

[tbonfort]

Author: warmerdam Date: 2010/03/25 - 20:38 Fix userid typos.

[tbonfort]

Author: warmerdam Date: 2010/03/26 - 06:05 I considered implementing the offsite mechanism but decided that it was a bit wrong since it only applied to raw raster data run through the resampler, and directly read data got no nodata processing.

So instead I have implemented preliminary support for a nodata bit array to indicate validity in the raw arrays of an imageObj. This is used to test pixels for validity and the logic is done during initial loading to a raw mode imageObj.

Potentially the mechanism might be extended to non-rawmode image data at some point.

The changes are in trunk only (e30dc1eefc831bf30f3bb88c0399ac483337e94b (r10017)) and should be considered very experimental.

[tbonfort]

Author: warmerdam Date: 2010/03/26 - 15:04 I have added modest testing (bda63e624b19ee064e5fe14e8e6f36321d77020b (r10019)).

I am going to close this ticket even though I have ignores Even's suggestion of strict nodata handling where any pixel influenced by nodata is emitted as nodata. I still think this is pretty specialized logic and I'm not inclined to address it in mapserver general release.