[Gdal-dev] gdalwarp + nosrcdata
Robert Osfield
robert at openscenegraph.com
Wed Mar 10 18:16:11 EST 2004
Hi,
I have been using the GDAL Warp API successful for a while without filtering
out no source data. I now am using a dataset that I'd like to manage the no
source data, and the API appears to support this, but as yet I haven't been
able to get it to work. I'm setting up both the padfSrcNoDataReal/Imag and
padfDstNoDataReal/Img, but none of these settings appear to make any
difference to the final resultant image.
I haven't found anything in the docs or warp tutorial to clarify all the
details, but looking at the gdalwap app source it looks like my
implementation shouldn't be too far off the mark.
I have tried out the gdalwap app as well, and still see the same problem with
no source data making its way into the final warp image.
I'm not sure what to expect from gdalwarp, so a few pointers of what to expect
to be able to do, and how to go about it would be great. I've pasted details
on what I've tried just in case this helps pinpoint the what steps I'm
missing.
Thanks for any help that you provide,
Robert.
---------------------------------------------------------------------------------------------------------------
Here's the gdalinfo report for the source data:
Driver: DTED/DTED Elevation Raster
Size is 3601, 3601
Coordinate System is:
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.257223563]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433]]
Origin = (104.999861,-9.999861)
Pixel Size = (0.00027778,-0.00027778)
Metadata:
DTED_VerticalAccuracy_UHL= NA
DTED_VerticalAccuracy_ACC= NA
DTED_SecurityCode_UHL=U
DTED_SecurityCode_DSI=U
DTED_UniqueRef_UHL=
DTED_UniqueRef_DSI=000000000000000
DTED_DataEdition=01
DTED_MatchMergeVersion=A
DTED_MaintenanceDate=0000
DTED_MatchMergeDate=0000
DTED_MaintenanceDescription=0000
DTED_Producer=
DTED_VerticalDatum=MSL
DTED_DigitizingSystem= ARC/INFO
DTED_CompilationDate=0210
DTED_HorizontalAccuracy= NA
DTED_RelHorizontalAccuracy= NA
DTED_RelVerticalAccuracy= NA
Corner Coordinates:
Upper Left ( 104.9998611, -9.9998611) (104d59'59.50"E, 9d59'59.50"S)
Lower Left ( 104.9998611, -11.0001389) (104d59'59.50"E, 11d 0'0.50"S)
Upper Right ( 106.0001389, -9.9998611) (106d 0'0.50"E, 9d59'59.50"S)
Lower Right ( 106.0001389, -11.0001389) (106d 0'0.50"E, 11d 0'0.50"S)
Center ( 105.5000000, -10.5000000) (105d30'0.00"E, 10d30'0.00"S)
Band 1 Block=1x3601 Type=Int16, ColorInterp=Undefined
Computed Min/Max=0.000,334.000
NoData Value=-32767
---------------------------------------------------------------------------------------------------------------
The gdalwarp options that I used to convert:
gdalwarp -t_srs '+proj=utm +zone=47 +datum=WGS84' -srcnodata -32767
s11e105.dt2 warped.tif
---------------------------------------------------------------------------------------------------------------
The gdalinfo of the final warped.tif:
Driver: GTiff/GeoTIFF
Size is 3657, 3693
Coordinate System is:
PROJCS["WGS 84 / UTM zone 47N",
GEOGCS["WGS 84",
DATUM["WGS_1984",
SPHEROID["WGS 84",6378137,298.2572235629972,
AUTHORITY["EPSG","7030"]],
AUTHORITY["EPSG","6326"]],
PRIMEM["Greenwich",0],
UNIT["degree",0.0174532925199433],
AUTHORITY["EPSG","4326"]],
PROJECTION["Transverse_Mercator"],
PARAMETER["latitude_of_origin",0],
PARAMETER["central_meridian",99],
PARAMETER["scale_factor",0.9996],
PARAMETER["false_easting",500000],
PARAMETER["false_northing",0],
UNIT["metre",1,
AUTHORITY["EPSG","9001"]],
AUTHORITY["EPSG","32647"]]
Origin = (1156568.438741,-1111402.317524)
Pixel Size = (30.74708996,-30.74708996)
Corner Coordinates:
Upper Left ( 1156568.439,-1111402.318) (104d58'49.94"E, 10d 0'0.76"S)
Lower Left ( 1156568.439,-1224951.321) (105d 0'1.07"E, 11d 1'17.53"S)
Upper Right ( 1269010.547,-1111402.318) (105d59'59.47"E, 9d58'48.96"S)
Lower Right ( 1269010.547,-1224951.321) (106d 1'22.49"E, 10d59'58.24"S)
Center ( 1212789.493,-1168176.819) (105d30'3.22"E, 10d30'2.88"S)
Band 1 Block=3657x1 Type=Int16, ColorInterp=Gray
Computed Min/Max=-32767.000,320.000
---------------------------------------------------------------------------------------------------------------
In my own code I've tried using a destination no data value of 0.0, but its
always overwritten by the -32767 no data value.
The relevant code is (note bits #if def'd out are experiements at trying to
get things to function):
// check to see if no value values exist in source datasets.
int numNoDataValues = 0;
for(i = 0; i < _sourceData->_gdalDataSet->GetRasterCount(); i++ )
{
int success = 0;
GDALRasterBand* band = _sourceData->_gdalDataSet->GetRasterBand(i+1);
band->GetNoDataValue(&success);
if (success) ++numNoDataValues;
}
if (numNoDataValues)
{
// no data values exist, so populate the no data arrays.
psWO->padfSrcNoDataReal = (double*)
CPLMalloc(psWO->nBandCount*sizeof(double));
psWO->padfSrcNoDataImag = (double*)
CPLMalloc(psWO->nBandCount*sizeof(double));
psWO->padfDstNoDataReal = (double*)
CPLMalloc(psWO->nBandCount*sizeof(double));
psWO->padfDstNoDataImag = (double*)
CPLMalloc(psWO->nBandCount*sizeof(double));
for(i = 0; i < _sourceData->_gdalDataSet->GetRasterCount(); i++ )
{
int success = 0;
GDALRasterBand* band = _sourceData->_gdalDataSet->GetRasterBand(i
+1);
double noDataValue = band->GetNoDataValue(&success);
double new_noDataValue = 0.0;
if (success)
{
std::cout<<"\tassinging no data value "<<noDataValue<<" to
band "<<i+1<<std::endl;
psWO->padfSrcNoDataReal[i] = noDataValue;
psWO->padfSrcNoDataImag[i] = 0.0;
psWO->padfDstNoDataReal[i] = new_noDataValue;
psWO->padfDstNoDataImag[i] = 0.0;
GDALRasterBandH band = GDALGetRasterBand(hDstDS,i+1);
GDALSetRasterNoDataValue( band, new_noDataValue);
}
}
#if 0
psWO->papszWarpOptions = (char**)CPLMalloc(2*sizeof(char*));
psWO->papszWarpOptions[0] = strdup("INIT_DEST=NO_DATA");
psWO->papszWarpOptions[1] = 0;
#endif
#if 0
psWO->pfnSrcValidityMaskFunc = (GDALMaskFunc)GDALWarpNoDataMasker;
#endif
}
---------------------------------------------------------------------------------------------------------------
More information about the Gdal-dev
mailing list