GDALAutoCreateWarpedVRT() on datasets without a defined NO_DATA value initializes grids to 0.0 instead of to chosen NO_DATA value

[latwood]

Found this while doing angleFromNorth stuff, I was wondering why wind_sd_to_uv() and wind_uv_to_sd() needed the NO_DATA check for the ascii grids, but not for the wxModel data. So why did the wxModel data place 0.0 values instead of NO_DATA values around the edges of the warped grids? And I found that it was because the pszWarpOptions used to try to set the NO_DATA value for the warp was being completely ignored.

The following style of code does NOT work, at least for the gdal version that we are currently using:

psWarpOptions->papszWarpOptions = CSLSetNameValue( psWarpOptions->papszWarpOptions, "INIT_DEST", "NO_DATA" );

psWarpOptions->padfSrcNoDataReal = (double*) CPLMalloc( sizeof( double ) * nBandCount );
psWarpOptions->padfSrcNoDataImag = (double*) CPLMalloc( sizeof( double ) * nBandCount );
psWarpOptions->padfDstNoDataReal = (double*) CPLMalloc( sizeof( double ) * nBandCount );
psWarpOptions->padfDstNoDataImag = (double*) CPLMalloc( sizeof( double ) * nBandCount );
for( i = 0; i < nBandCount; i++ )
{
    psWarpOptions->padfSrcNoDataReal[i] = dfNoData;
    psWarpOptions->padfSrcNoDataImag[i] = dfNoData;
    psWarpOptions->padfDstNoDataReal[i] = dfNoData;
    psWarpOptions->padfDstNoDataImag[i] = dfNoData;
}

Neither does going through each band, pre- or post- warp, trying to set the NO_DATA value for the given band
GDALSetRasterNoDataValue(hBand, dfNoData);

Tbh I'm surprised that GDALSetRasterNoDataValue(hBand, dfNoData); doesn't throw an error when used on a gdal dataset opened as GA_ReadOnly. Regardless, doing so doesn't alter the NO_DATA value to the pre-warp dataset, doing so DOES alter the NO_DATA value for the post-warp dataset but the raw NO_DATA values introduced by the warp are still 0.0 in value.

I think what is going on, is that "INIT_DEST"="NO_DATA" tries to force the warp to use the NO_DATA value as the initial value for the new warped grid, but NO_DATA is not set as a value for the source grid, so 0.0 is used as the NO_DATA value instead. And so the warped grid is initialized with values of 0.0 instead of using the dfNoData value as the initial value.

I did find a work around, "INIT_DEST"="<NO_DATA_VALUE>", just put the raw NO_DATA value that you want into that spot, rather than the "NO_DATA" text. So basically the following block of code:

GDALDatasetH hBand = GDALGetRasterBand( hSrcDS, band );
int bSuccess = false;
double dfNoData = GDALGetRasterNoDataValue( hBand, &bSuccess );
if( bSuccess == false )
{
    dfNoData = -9999.0;
}

psWarpOptions->papszWarpOptions = CSLSetNameValue( psWarpOptions->papszWarpOptions, "INIT_DEST", "NO_DATA" );
if( bSuccess == false )  // if GDALGetRasterNoDataValue( hBand, &bSuccess ) fails to return that a NO_DATA value is in the source dataset
{
    psWarpOptions->papszWarpOptions = CSLSetNameValue( psWarpOptions->papszWarpOptions, "INIT_DEST", boost::lexical_cast<std::string>(dfNoData).c_str() );
}

I've already started using this method successfully for the ascii grid warping in griddedInitialization.cpp and ninja.cpp, during calls to gdal_util.cpp GdalWarp().

I tried the method out in nomads_wx_init.cpp, but having proper NO_DATA values downstream from that point still needs to be properly checked. It seemed to work fine, it ran to completion and the weather model kmz outputs now drop the band of zero values around the edges, but it still needs to be better checked.

This needs to be verified to be a proper method, then all instances of GDALAutoCreateWarpedVRT() in the code base will need to be checked and updated, and all downstream code will need to be checked and updated to more properly handle actually set NO_DATA values if applicable.