[gdal-dev] GDAL raster processing library?

Even Rouault even.rouault at spatialys.com
Tue Jul 12 06:55:41 PDT 2016


James,

There's some intersection with Ari's proposal :
https://trac.osgeo.org/gdal/wiki/rfc62_raster_algebra . At least regarding the 
overall purposes, since technically this is quite different.

Actually what you propose is very close to the existing VRT pixel functions of 
derived bands. See "Using Derived Bands" in 
http://www.gdal.org/gdal_vrttut.html . In the last days, we've merged 
Antonio's work regarding a predefined set of pixel functions.
Perhaps some extension to allow passing user parameters to the pixel func 
could be useful. It is possible to use pixel functions from Python as shown in 
https://github.com/OSGeo/gdal/blob/trunk/autotest/gcore/testnonboundtoswig.py#L303 
although this is a bit ugly as it uses ctypes and not SWIG. But should be 
possible through SWIG by introducing proper types similarly to what is done 
for the progress functions or error handler functions.

Even

Le mardi 12 juillet 2016 14:40:27, jramm a écrit :
> I wonder if there is a use case/interest in a small library within GDAL to
> facilitate generic raster processing.
> 
> My idea would be to have a user-extensible framework to run processing
> functions on whole rasters, with a growing library of common-such functions
> within GDAL.
> 
> Something along the lines of this:
> 
> /*******************************************************************/
> typedef CPLErr (*GDALRasterProcessFn)(double *padfInArray, double
> *padfOutArray,
>                 int nWindowXSize, int nWindowYSize, void *pData, double
> *pdfNoDataValue);
> /**
> * \brief Definition of a raster processing function.
> *
> * A GDALRasterProcessFn accepts an array of data as input, applies custom
> logic and writes the output to padfOutArray.
> * Such a function can be passed to GDALRunRasterProcess to apply custom
> processing to a GDALDataset in chunks and create
> * a new GDALDataset.
> *
> * @param padfInArray The input array of data.
> *
> * @param padfOutArray The output array of data. On first call (via
> GDALRunRasterProcess) this will be an empty, initialised array,
> *    which should be populated with the result of calculations on
> padfInArray. In subsequent calls it will contain the result of the
> *    previous window.
> *
> * @param nWindowXSize the actual x size (width) of the read window.
> *
> * @param nWindowYSize the actual y size (height) of the read window. The
> length of padfInArray == padfOutArray == nWindowXSize * nWindowYSize
> *
> * @param pData Process-specific data. This data is passed straight through
> to the GDALRasterProcessFn and may contain e.g user defined parameters.
> *     The GDALRasterProcessFn definition would define the structure/type of
> such data.
> *
> * @param pdfNoDataValue The no data value of the dataset
> */
> 
> CPLErr GDALRunRasterProcess(GDALRasterProcessFn processFn, GDALDataset
> *poSrcDataset,
>                         GDALDataset *poDstDataset, void *pData, int
> *pnWindowXSize,
>                         int *pnWindowYSize, int *pnPixelBuffer);
> /**
> * \brief Apply a raster processing function to each sub-window of a raster.
> *
> * The input raster dataset is read in chunks of nWindowXSize * nWindowYSize
> and each chunk is passed to the processing
> * function. The output array from the function is written to the
> destination dataset.
> * An optional 'pixel buffer' can be specified to allow overlaps between
> successive windows. This is useful for
> * some algorithms, e.g. blob extraction, watershed/stream flow analysis,
> convolution etc.
> * Process specific data can be passed (e.g. configuration parameters). This
> data is simply passed straight through to the processing
> * function on each call.
> *
> * @param processFn A GDALRasterProcessFn to apply to each sub window of the
> raster.
> *
> * @param poSrcDataset The source raster dataset from which pixel values are
> read
> *
> * @param poDstDataset The destination raster dataset to which pixel values
> are written. Must support RasterIO in write mode.
> *
> * @param pData Process-specific data. This is passed straight through to
> the GDALRasterProcessFn on each call.
> *
> * @param pnWindowXSize The desired width of each read window. If NULL it
> defaults to the 'natural' block size of the raster
> *
> * @param pnWindowYSize The desired height of each read window. If NULL it
> defaults to the 'natural' block size.
> *
> * @param pnPixelBuffer A pixel buffer to apply to the read window. The read
> window is expanded by pnPixelBuffer pixels in all directions such that
> *    each window overlaps by pnPixelBuffer pixels. Ignored if NULL or 0
> *
> * @return a CPLErr enum indicating whether the function succesfully ran.
> */
> 
> 
> CPLErr GDALRunRasterProcesses(GDALRasterProcessFn *paProcessFn, int
> nProcesses, GDALDataset *poSrcDataset,
>                         GDALDataset *poDstDataset, void **paData, int
> *pnWindowXSize,
>                         int *pnWindowYSize, int *pnPixelBuffer);
> 
> /**
> * \brief Apply multiple raster processing functions to each sub-window of a
> raster
> *
> * For each window, the functions defined by the paProcessFn array are
> called in turn, with the array output of the previous function forming the
> input * to the next function. This allows processing 'toolchains' to be
> built without having to create intermediate datasets, which can be less
> efficient in time and space.
> *
> *
> * @param paProcessFn An array of GDALRasterProcessFn to apply to each sub
> window of the raster
> *
> * @param nProcesses The size of paProcessFn
> *
> * @param poSrcDataset The source raster dataset from which pixel values are
> read
> *
> * @param poDstDataset The destination raster dataset to which pixel values
> are written. Must support RasterIO in write mode.
> *
> * @param paData an array of process-specific data objects of size
> nProcesses. Each data object will be passed to the corresponding
> GDALRasterProcessFn
> *
> * @param pnWindowXSize The desired width of each read window. If NULL it
> defaults to the 'natural' block size of the raster
> *
> * @param pnWindowYSize The desired height of each read window. If NULL it
> defaults to the 'natural' block size.
> *
> * @param pnPixelBuffer A pixel buffer to apply to the read window. The read
> window is expanded by pnPixelBuffer pixels in all directions such that
> *    each window overlaps by pnPixelBuffer pixels. If NULL, it is ignored.
> *
> * @return a CPLErr enum indicating whether the function succesfully ran.
> */
> /*******************************************************************/
> 
> So GDALRunRasterProcess would take care of I/O and looping through the
> dataset in chunks, which can be configure by the user. Each chunk is
> submitted to the processing function, which does 'stuff' and populates the
> output array. GDALRunRasterProcess would write it to the output dataset
> object.
> 
> An example of a processing function could be something like this:
> 
> /*******************************************************************/
> typedef struct ReclassArgs {
>     int nReclassArgs;
>     double *padfMinBound;
>     double *padfMaxBound;
>     double *padfBinValue;
> 
> } ReclassArgs;
> 
> 
> CPLErr Reclass(double *padfInArray, double *padfOutArray,
>                int nWindowXSize, int nWindowYSize, void *pData,
>                double *pdfNoDataValue)
> {
>     ReclassArgs *poRargs = static_cast<ReclassArgs *>(pData);
>     double pixel;
>     for (int i = 0; i < nWindowXSize * nWindowYSize - 1; ++i){
>         pixel = padfInArray[i];
>         for (int k = 0; k < poRargs->nReclassArgs; ++k){
>                 if ((pixel > poRargs->padfMinBound[k]) && (pixel <=
> poRargs->padfMaxBound[k])){
>                     padfOutArray[i] = poRargs->padfBinValue[k];
>                     break;
>                 }
>         }
>     }
>     return CE_None;
> }
> /*******************************************************************/
> 
> To use it you would:
> 
> /*******************************************************************/
> // Open a source dataset
> GDALDataset *poSrc = GDALOpen(...)
> 
> // Create some output dataset with same dimensions etc..
> GDALDataset *poDst = ....
> 
> // Setup the reclass args
> ReclassArgs *poRargs = ...
> 
> // Call the process
> CPLErr retval = GDALRunRasterProcess(Reclass, poSrc, poDst, (void
> *)poRargs, NULL, NULL, NULL);
> /*******************************************************************/
> 
> GDAL could feature its own library of common functions...reclass, truncate,
> blob extraction etc etc and a user of GDAL could write their own and have
> GDAL process it.
> 
> It would be useful to expose enough to python via SWIG so that a
> GDALRasterProcessFn could be defined in python - this would then allow
> rapid development in python, but push the expensive work of looping over
> dataset chunks to C++.
> I dont know enough about SWIG yet to see whether this would be feasible.
> 
> However, even in C++, it allows the user to only think about their
> algorithm and leave the boilerplate to GDAL.
> 
> I imagine all this being implemented in a 'sub library', something like
> gdpl.h (geospatial/gdal dataset processing library), but sitting in the
> GDAL source tree (then I can happily expose all the useful GDAL goodness
> without being worried about users having to include gdal headers etc, and
> avoid any work in having to abstract it away)
> 
> This is a first stab at defining something like this and is probably not
> yet comprehensive enough.
> 
> Any thoughts?
> 
> 
> 
> 
> --
> View this message in context:
> http://osgeo-org.1560.x6.nabble.com/GDAL-raster-processing-library-tp52759
> 48.html Sent from the GDAL - Dev mailing list archive at Nabble.com.
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/gdal-dev

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com


More information about the gdal-dev mailing list