Iteration over a raster coverage - 52North/IlwisCore GitHub Wiki

This text describes the various ways how to use the pixel-iterator both in an explicit way as in an implicit way.

Contruction

IRasterCoverage raster("file:///d:/data/someraster.img");
PixelIterator iterator(raster); // iterator works over the whole coverage

PixelIterator iterator(raster, BoundingBox(Pixel(100,200), Pixel(750, 1200)); // over a window

UPGeometry geom(GeometryHelper::fromWkt("Polygon((12.34 34.22, 17.22 36.78,...))"));
PixelIterator iterator(raster, geom); // the internals of a polygon(will respect holes and islands)

UPGeometry geom(GeometryHelper::fromWkt("Multipoint((12.56 33.90),(34.13 45.99),...)"));
PixelIterator iterator(raster, geom); // only iterate over the points

Usage

As explained in the referred text, the iterator behaves like a normal STL container iterator. So normal STL algorithms like

IRasterCoverage raster("file:///d:/data/someraster.img");
std::fill(begin(raster), end(raster), 0); // whole map with zeroes

//move to the first pixel with a value bigger than 200
PixelIterator iter = std::find_if(begin(raster), end(raster), [](double v){ return v > 200;});

should work just fine. Note that the examples work over the whole raster. Of course one could define for example a 'windowed' iterator and use that

IRasterCoverage raster("file:///d:/data/someraster.img");
PixelIterator iterator(raster, BoundingBox(Pixel(100,200), Pixel(750, 1200));
std::fill(iter, iter.end(), 0); // only a window filled with zeroes

Apart from using the stl to do some (simple) processing we can of course also implement our own algorithms. First of all, there is an implementation of the free methods end() and begin() ( consistent with c++11) for rastercoverages. So the next will work

IRasterCoverage raster("file:///d:/data/someraster.img");
PixelIterator beginIter = begin(raster);
PixelIterator endIter = end(raster);

This also implies that the implicit iterator like in the loop below, will work. Suppose i want to know how many pixels have a value over 100. Below is a crude way of doing this

int count = 0;
IRasterCoverage myraster("file:///d:/data/ethiopia-20110322-0845.img");
PixelIterator iter(myraster);
while(iter != end(myraster)){
    if ( *iterator > 100)
       ++count;
}

One could also do (using the implicit iterator)

unsigned int count = 0;
for(auto value : raster)
   if ( value > 100) ++count;

which is a bit more compact. Values are returned by reference, so one can change them change them. E.g. all value below 100 to 0.

for( auto& value : myraster)
  if ( value < 100) value = 0;

The default 'flow' of a iterator is XYZ (see [also] (https://github.com/52North/IlwisCore/wiki/pixeliterator#description) But what if you need a different flow. Say ZXY, basically band interleaved per pixel. It can usefull when aggregating over spatial temporal data or classifying. So lets do a simpel aggregation.

 IRasterCoverage raster("file:///d:/data/spatialtemporal.img");
 PixelIterator iterIn(myraster);
 iterIn.flow(PixelIterator::fZXY);
 
 IRasterCoverage aggregate;
 OperationHelperRaster::initialize(raster,aggregate, itDOMAIN | itGEOREF | itCOORDSYSTEM);
 PixelIterator iterOut(aggregate);
 double sum = 0; 

 for_each(iterIn, iterIn.end(),[&](double v){
    if (iterIn.xChanged()){ // flow goes zxy, so when a z column is done, x will change
       *iterOut = sum / raster->size().zsize();
       sum = 0 ; 
       ++iterOut;
    }
    sum += v;
 });
 *iterIn = sum; // adding the last column
⚠️ **GitHub.com Fallback** ⚠️