>Ari Wrote:
>At FOSS4G my conclusion was that my first approach was not good for many 
>reasons, most >importantly because it did not scale to several bands in one 
>operation. So I've trying with the >following ideas/decisions:
>
>The problem is to compute y=f(x1, x2, ...), where y is a new dataset with one 
>band or an existing >dataset, into which a new band is added. 
>x1, x2, ... are existing bands. f is an expression. The goal is to be able to 
>create an expression >object, with which one can write

Hi Ari,

Perhaps it is worthwhile looking back at the email I wrote earlier in the 
process: 

http://permalink.gmane.org/gmane.comp.gis.gdal.devel/42609 

At the time, I suggested using expression objects and to allow functions with 
more than 2 arguments and gave some suggestions on the implementation.

In my experience the key first step is to develop pixel iterators for raster 
bands. Once you have that, raster bands can be wrapped as regular ranges and 
you can make use of 'regular' tools and methods to develop your expression 
objects, e.g. as in the Range v3 library by Eric Niebler.

Hence, my repeated suggestion for GDAL: create iterators for raster bands :)

As you may recall, the code below works just fine (but requires C++11).

With kind regards, Alex


#include <blink/map_algebra/map_algebra.h> 

namespace ma = blink::map_algebra;

int my_function(int w, int x, int y, int z)
{
  return w * x + y * z;
}

int main()
{
  auto a = ma::open_read_only<int>("input_1.tif");
  auto b = ma::open_read_only<int>("input_2.tif");
  
  auto w = ma::create_from_model<double>("output_1.tif", a);
  auto x = ma::create_from_model<double>("output_2.tif", a);
  auto y = ma::create_from_model<double>("output_3.tif", a);
  auto z = ma::create_from_model<double>("output_4.tif", a);

  // Example 1: Using operators
  w = a + 3 * b;

  // Example 2: Using assigning operators
  x = 1;
  x *= a;
  x += 3 * b;
  
  // Example 3: Map algebra using cell-by-cell functions
  y = ma::apply(my_function, 1, a, 3, b);

  // Example 4: Combination
  z = b + 3 * ma::apply(my_function, 1, a, 2, b);

  return 0;
}


_______________________________________________
gdal-dev mailing list
[email protected]
http://lists.osgeo.org/mailman/listinfo/gdal-dev

Reply via email to