> -----Original Message----- > From: gdal-dev [mailto:gdal-dev-boun...@lists.osgeo.org] On Behalf Of > Ari Jolma > Sent: 08 April 2016 07:05 > To: gdal-dev@lists.osgeo.org > Subject: Re: [gdal-dev] Map algebra > > 08.04.2016, 02:21, alex kirjoitti: > >> That is one aspect, yes. However, so far I've considered only > >> operations on one or two bands at time. Operations on more than two > >> can of course be dissolved into a series of operations with max two at a > >> time. > > That will often be possible, but not always. For instance how would > > you dissolve the following function? (without resorting to complete > > ugliness like c * a + !c * b) > > > > int conditional(bool c, int a, int b) { return c ? a : b; } > > As an assignment that would be > > d = c ? a : b; > > If all variables are spatial then that's four bands. To save space > we'll do in-place modifications to bands. If we do not want to change > any of the bands a, b and c, we'll begin by making d a copy of b: > > d = b; > > and then use a conditional > > d = a if c; >
Ha, that was easier than I thought :) But also quite inefficient. And you'll still need to negotiate three different block sizes. [snip] > >> The approach is also strongly supported by the basic design of GDAL > >> and the fundamental algorithm is presented in GDAL documentation (I > >> gave the link in earlier email). > > I didn't see that email and just trawled through the whole thread, > > can you give > it again? > > http://www.gdal.org/classGDALRasterBand.html#ad80cecc562fd48f5783677de > 625360ac > > it is a simple algorithm but demonstrates the two level approach. OK, I also use the ReadBlock and WriteBlock functions instead of the RasterIO. But I still don't see how the principle of iterating over blocks and within blocks is an approach that can be used for two or three bands with different block sizes. > > > >> The solution I already have seems to work well for the one and two > >> band case even with 3x3 focal region and as far as I can see it > >> could be extended to cases involving more bands. > > I'll add the disclaimer about the curse of dimensionality I mention > above for which a solution should be found first. It could be a simple > one like copy all blocks to same datatype first. Can you describe the principle of your method? And where your curse of dimensionality comes in. Perhaps somebody has a solution. [snip] > > It really depends on the block size, which is something GDAL core and > drivers handle. > >>> 3. Zonal operations: operations that calculate indicators for > >>> groups > >>> (zones) of > >> pixels. > >>> This seems more trivial than the other two components of map > >>> algebra and I > >> haven't given it much attention, but I think should be supported. > >> > >> One complication, which I just today worked on, seems to stem from > >> one requirement that I have. I want the API be free of templates > >> and thus there is a need for generic classes for, e.g., numbers > >> (GDAL supports bands where data is various length integers or floats). > >> Thus, for example a simple zonal method, which returns which zones > >> share a boundary returns a rather complex data structure. > > I understand zonal methods to be operating on pixels within zones > > and not on > the relationship between zones. What kind of methods do you have in mind? > > I admit that it is perhaps not a "pure" zonal method. For example, > think about the question: "What are the neighboring countries of each > country?". Sure, that is useful functionality, but I was trying to find out what is the rationale and purpose of your Map Algebra library. I think I misunderstood your original announcement email. I thought you intended to offer generic support for Map Algebra type analysis. Instead it seems that you are offering a specific set of Map Algebra functions (as I can now see in the RFC). _______________________________________________ gdal-dev mailing list gdal-dev@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/gdal-dev