Category Archives: OTB

OTB++11

Modern C++

I have just recently had a look at the last C++ ISO standard (known as C++11). I was not aware of all the fancy features that have been included in the language.

Some of the most interesting features are:

  • auto,
  • range-based for loops,
  • general constant expressions,
  • uniform initialisation,
  • type aliases,
  • decltype,
  • variadic templates,
  • type safe threading,
  • user-defined literals,
  • move semantics

There are also new containers and algorithms in the STL.

During the end of the year break, I spent some time reading and watching on line presentations about the new C++. I recommend watching Bjarne Stroustrup’s keynote at Going Native 2012 (link) and Scott Meyer’s talk at Going Native 2013 where he even anticipates some C++14 features (link).

Stroustrup introduces his talk by saying:

I use C++11 freely. Examples include auto, general constant expressions, uniform initialisation, type aliases, type safe threading, and user-defined literals. C++11 features are only just starting to appear in production compilers, so some of my suggestions are conjecture. Developing a “modern style,” however, is essential if we don’t want to maintain newly-written 1970s and 1980s style code in 2020.

After listening at him, I realised that the most easy to use features alone could make C++ programming much more fun and simple. Since I have spent some time in the last 2 years flirting with functional and dynamic typed (or at least with type inference) languages (Scheme, Clojure, Python, SML), I wanted to see how using C++11 affects the style of OTB code.

As usual, GCC has been the first compiler to make available these new features (closely followed by Clang). C++11 is available in recent versions of GCC using -std=c++11. So I recompiled OTB with this flag and started modifying some of the Software Guide examples using C++11 features.

All in all, I can say that these features make OTB code shorter, clearer and maybe safer. I report some examples here.

Examples

auto

The auto keyword allows to declare variables without explicitly stating their type and puts this burden on the compiler. This is not dynamic typing, this is type inference and nothing changes in terms of what the compiler generates. So you can write:

auto x = sqrt(2.0);

and the type of x will be deduced from the return value type of sqrt. You can also do this:

auto x=0;
auto y=0.0;

and x will be of type int and y will be of type double. If you want to force other types, you have to use the classical declarations.

How to benefit from auto when using OTB? A common pattern in OTB code is to instantiate a smart pointer using the New() class method. Sometimes, we just want to create one object of a given type and we do this:

typedef otb::ImageFileReader<ImageType> ReaderType;
ReaderType::Pointer reader = ReaderType::New();

typedef otb::ImageFileWriter<ImageType> WriterType;
WriterType::Pointer writer = WriterType::New();

in order to create a reader and a writer. That is, since type names are long we define a type using typedef and then we use it twice, once for Pointer and once for New.

Using auto, we can reduce the boilerplate like this:

auto reader = otb::ImageFileReader<ImageType>::New();
auto writer = otb::ImageFileWriter<ImageType>::New();

The same applies for filters, for which usually we only use one instance:

typedef itk::RescaleIntensityImageFilter<ImageType,
                                         OutputImageType> RescalerType;
RescalerType::Pointer rescaler = RescalerType::New();

becomes

auto rescaler = itk::RescaleIntensityImageFilter<ImageType,
                                                 OutputImageType>::New();

which is once again much shorter.

Another example in a for loop:

for (PolygonListType::Iterator it = polygonList->Begin();
     it != polygonList->End(); ++it)
  {
  DataNodeType::Pointer newPolygon = DataNodeType::New();
  newPolygon->SetPolygonExteriorRing(it.Get());
  tree->Add(newPolygon, multiPolygon);
  }

becomes

for (auto it = polygonList->Begin();
     it != polygonList->End(); ++it)
  {
  auto newPolygon = DataNodeType::New();
  newPolygon->SetPolygonExteriorRing(it.Get());
  tree->Add(newPolygon, multiPolygon);
  }

using

The using keyword allows to define template aliases. Whenever we write this:

typedef otb::Image<unsigned int, 2> ImageType;
typedef otb::VectorImage<unsigned int, 2> VectorImageType;
typedef otb::Image<double, 2> DoubleImageType;
typedef otb::VectorImage<double, 2> DoubleVectorImageType;

we could do this:

using ImageType = otb::Image<unsigned int, 2>;
using VectorImageType = otb::VectorImage<unsigned int, 2>;
using DoubleImageType = otb::Image<double, 2>;
using DoubleVectorImageType = otb::VectorImage<double, 2>;

For people like me who have always found counterintuitive the order of the types in a typedef, this is much clear. I know that the code is not shorter, but for a C++ beginner, this should be much clearer.

One cool thing that using allows and that typedef doesn’t is partial template specialisation1:

template<typename T> 
using Image2D = otb::Image<T, 2>;
auto im = Image2D<int>::New();

decltype

Sometimes, we don’t want to mentally track the type names, but we know that to objects should have the same type. decltype allows to declare the type of one object by using the name of another object:

ImageType::IndexType start = {0, 0};
decltype(start) start2 = {0, 0};

In OTB, at first, I didn’t see a way to use decltype on smart pointers and templates, since we usually manipulate ObjectType::Pointer, but the template parameters are simply ObjectType. Fortunately, ITK designers declared an ObjectType in the SmartPointer class which is an alias for the type of the pointed object.

I found an example where this is useful (the orthofusion example in the Tutorial chapter of the Software Guide).

The original code defines a type for a projection and creates an object:

typedef otb::GenericMapProjection<otb::TransformDirection::INVERSE>
        InverseProjectionType;

InverseProjectionType::Pointer utmMapProjection =
                               InverseProjectionType::New();

Later in the program, we can use the defined type as a template parameter for a filter:

typedef otb::OrthoRectificationFilter<ImageType, DoubleImageType,
                                      InverseProjectionType> 
             OrthoRectifFilterType;
OrthoRectifFilterType::Pointer orthoRectifPAN =
                               OrthoRectifFilterType::New();

If we created the projection without defining the type like this:

auto utmMapProjection =
     otb::GenericMapProjection<otb::TransformDirection::INVERSE>::New();

how can we then instantiate the filter with a type which has not been defined? Like this:

auto orthoRectifPAN =
      otb::OrthoRectificationFilter<ImageType,
                                    DoubleImageType,
                                    decltype(utmMapProjection)::ObjectType>::New();

I am not sure whether this is a good solution in production code which has to be maintained long time after it was written, but it is for sure an interesting solution when prototyping processing chains.

Lambdas

Lambdas, or anonymous functions, are one of the most exciting things in C++11 when one has used them in LISP-like languages or even in Python. They make fast prototyping much easier mainly when making a heavy use of the STL algorithms because they can replace functors or plain functions.

Usually, STL algorithms can be specialised by passing them functions as parameters. If one wants to store state, a functor or function object (a class which defines the () operator) is used.

The inconvenient of functions and functors is that they are defined away from where they are used, so they are not ideal when they are only used once. Functors are worse than functions since they need much boilerplate for only one useful method.

Lambdas are the ideal solution for this, since they are defined locally. They can store state by capturing local variables either by reference or by value.

OTB iterators are not fully compatible with the STL, so one can think that lambdas would not be very useful with OTB specific code. But actually, there are many filters which can be implemented by specialising one of the N-ary functor filters. These filters are useful when one wants to implement a filtering operation which can’t be obtained with the otb::BandMathFilter, as for instance when a local neighbourhood is needed.

One example of this is the change detection framework. A long time ago when I wrote this example I had to use a lot of code to achieve a simple thing like computing the mean difference around a pixel.

First of all, I needed a functor to implement the processing in the neighbourhood:

template<class TInput1, class TInput2, class TOutput>
class MyChangeDetector{
public:
  MyChangeDetector() {}
  ~MyChangeDetector() {}
  inline TOutput operator ()(const TInput1& itA,
                             const TInput2& itB)
  {
    TOutput result = 0.0;

    for (unsigned long pos = 0; pos < itA.Size(); ++pos)
      {

      result += static_cast<TOutput>(itA.GetPixel(pos) - itB.GetPixel(pos));

      }
    return static_cast<TOutput>(result / itA.Size());
  }
};

That is, a class definition with empty constructor and destructor. Then I needed to inherit from the binary functor which performs a neighbourhood filtering and specialise it using my functor.

template <class TInputImage1, class TInputImage2, class TOutputImage>
class ITK_EXPORT MyChangeDetectorImageFilter :
  public otb::BinaryFunctorNeighborhoodImageFilter<
      TInputImage1, TInputImage2, TOutputImage,
      MyChangeDetector<
          typename itk::ConstNeighborhoodIterator<TInputImage1>,
          typename itk::ConstNeighborhoodIterator<TInputImage2>,
          typename TOutputImage::PixelType> >
{
public:
  typedef MyChangeDetectorImageFilter Self;

  typedef typename otb::BinaryFunctorNeighborhoodImageFilter<
      TInputImage1, TInputImage2, TOutputImage,
      MyChangeDetector<
          typename itk::ConstNeighborhoodIterator<TInputImage1>,
          typename itk::ConstNeighborhoodIterator<TInputImage2>,
          typename TOutputImage::PixelType>
      >  Superclass;

  typedef itk::SmartPointer<Self>       Pointer;
  typedef itk::SmartPointer<const Self> ConstPointer;

  itkNewMacro(Self);

protected:
  MyChangeDetectorImageFilter() {}
  virtual ~MyChangeDetectorImageFilter() {}

private:
  MyChangeDetectorImageFilter(const Self &);
  void operator =(const Self&);

};

So again, a class definition for nothing. I could maybe just have defined a type, but the code using it would not have been as simple as this:

typedef MyChangeDetectorImageFilter<InputImageType1, InputImageType2,
    ChangeImageType>      FilterType;
FilterType::Pointer   filter = FilterType::New();

Even if I had defined a type alias for the filter in order to avoid to inherit from the existing filter, the functor definition is an overkill and it is located far from where it is used.

Using lambdas, the code becomes this:

using CPT = ChangeImageType::PixelType;
using CNIType1 = itk::ConstNeighborhoodIterator<InputImageType1>;
using CNIType2 = itk::ConstNeighborhoodIterator<InputImageType2>;

That is, 3 lines of type aliases for a simpler notation. Then the lambda:

std::function<CPT (CNIType1, CNIType2)>
  cdt = [](CNIType1 itA, CNIType2 itB){
  CPT result{0};
  for (auto pos = 0; pos < itA.Size(); ++pos)
    result += static_cast<CPT>(itA.GetPixel(pos) - itB.GetPixel(pos));
  return static_cast<CPT>(result / itA.Size());
};

And finally the filter instantiation:

auto filter = otb::BinaryFunctorNeighborhoodImageFilter<InputImageType1,
                                                        InputImageType2,
                                                        ChangeImageType,
                                                        decltype(cdt)>::New();

If we had wanted to have a function with state, objects in the scope where the lambda is defined can be captured by value or by reference within the [], so we are very close to what a functor may achieve with much less code.

The main advantage here is that all the code is in the same location thanks to the possibility of creating a function anywhere in the sources. If I understand it well, lambdas are inlined, which means that the code is faster than a function call.

Conclusion

I have only started to use this kind of code, but I think that once you start using auto and lambdas it must be very difficult to go back to C++98/03.

The only drawback with using C++11 is dropping compatibility with old compilers, but this is not an issue for most of the things I do. So I can start writing modern OTB code instead of poo (plain old OTB)!

Footnotes:

1

With GCC 4.7 this only seems to work out of a scoped block.

Hack OTB on the Clojure REPL

The WHY and the WHAT

Lately, I have been playing with several configurations in order to be able to fully exploit the ORFEO Toolbox library from Clojure. As explained in previous posts, I am using OTB through the Java bindings.

Recently, the folks of the OTB team have made very cool things – the application wrappers – which enrich even further the ways of using OTB from third party applications.

The idea of the application wrappers is that when you code a processing chain, you may want to use it on the command line, with a GUI or even as a class within another application. What Julien and the others have done is to build an application engine which given a pipeline and its description (inputs, outputs and parameters) is able to generate a CLI interface and a Qt-based GUI interface.

But once you are there, there is no major problem (well, if your name is Julien, at least) to generate a SWIG interface which will then allow to build bindings for the different languages supported by SWIG.

This may seem a chicken sink – usine à gaz in French – and actually it is if you don’t pay attention (like me). The application wrappers use advanced techniques which require a recent version of GCC. At least recent from a Debian point of view (GCC 4.5 or later).

On the other hand, the bindings for the library (the OTB-Wrapping package) use CableSwig, which in turn uses its own SWIG interfaces.

So if you want to be able to combine the application wrappers with the library bindings, odds are that you may end up with a little mess in terms of compiler versions, etc.

And this is exactly what happened to me when I had to build GCC from source in Debian Squeeze (no GCC 4.5 in the backports yet). So I started making some experiments with virtual machines in order to understand what was happening.

Eventually I found that the simplest configuration was the best to sort things out.

This is the HOW-TO for having a working installation of OTB with application wrappers, OTB-Wrapping and access everything from a Clojure REPL.

 

 

The HOW-TO

Installing Arch Linux on a virtual machine

I have chosen Arch Linux because with few steps you can have a minimal Linux install. Just grab the net-install ISO and create a new virtual machine on VirtualBox. In terms of packages just get the minimal set (the default selection, don’t add anything else).

After the install, log in as root, and perform a system update with:

pacman -Suy

If you used the net-install image, there should be no need for this step, but in case you used an ISO with the full system, it’s better to make sure that your system is up to date.

Then, create a user account and set up a password if you want to:

useradd -m jordi
passwd jordi

Again, this is not mandatory, but things are cleaner this way, and you may want to keep on using this virtual machine. Don’t forget to add your newly created user to the sudoers file:

pacman -S sudo
vi /etc/sudoers

You can now log out and log in with the user name you created. Oh, yes, I forgot to point out that you don’t have Xorg, or Gnome or anything graphical. Even the mouse is useless, but since you are a real hacker, you don’t mind.

 

 

Install the packages you need

Now the fun starts. You are going to keep things minimal and only install the things you really need.

sudo pacman -S mercurial make cmake gcc gdal openjdk6 \
               mesa swig cvs bison links unzip rlwrap

Yes, I know, Mesa (OpenGL) is useless without X, but we want to check that everything builds OK and OTB uses OpenGL for displaying images.

  • Mercurial is needed to get the OTB sources, cvs is needed for getting CableSwig
  • Make, cmake and gcc are the tools for building from sources
  • Swig is needed for the binding generation, and bison is needed by CableSwig
  • OpenJDK is our Java version of choice
  • Links will be used to grab the Clojure distro, unzip to extract the corresponding jar file and rlwrap is used by the Clojure REPL.

And that’s all!

 

 

Get the source code for OTB and friends

I prefer using the development version of OTB, so I can complain when things don’t work.

hg clone http://hg.orfeo-toolbox.org/OTB
hg clone http://hg.orfeo-toolbox.org/OTB-Wrapping

CableSwig is only available through CVS.

cvs -d :pserver:anonymous@pubilc.kitware.com:/cvsroot/CableSwig \
        co CableSwig

 

 

Start building everything

We will compile everything on a separate directory:

mkdir builds
cd builds/

 

OTB and OTB-Wrapping

We create a directory for the OTB build and configure with CMake

mkdir OTB
cd OTB
ccmake ../../OTB

Don’t forget to set to ON the application build and the Java wrappers. Then just make (literally):

make

By the way, I have noticed that the compilation of the application engine can blow up your gcc if you don’t allocate enough RAM for your virtual machine.

At this point, you should be able to use the Java application wrappers. But we want also the library bindings so we gon on.

We can now build CableSwig which will be needed by OTB-Wrapping. Same procedure as before:

cd ../
mkdir CableSwig
cd CableSwig/
ccmake ../../CableSwig/
make

And now, OTB-Wrapping. Same story:

cd ../
mkdir OTB-Wrapping
cd OTB-Wrapping/
ccmake ../../OTB-Wrapping/

In the cmake configuration, I choose only to build Java, but even in this case a Python interpreter is needed. I think that CableSwig needs it to generate XML code to represent the C++ class hierarchies. If you did not install Python explicitly in Arch, you will have by default a 2.7 version. This is OK. If you decided to install Python with pacman, you will have both, Python 2.7 and 3.2 and the default Python executable will point to the latter. In this case, don’t forget set the PYTHON_EXECUTABLE in CMake to /usr/bin/python2.

Then, just make and cd to your home directory.

make
cd

And you are done.

Well not really. Right now, you can do Java, but what’s the point? You might as well use the C++ version, right?

 

 

 

Land of Lisp

Since Lisp is the best language out there, and OTB is the best remote sensing image processing software (no reference here, but trust me), we’ll do OTB in Lisp.

 

PWOB

You may want to get some examples that I have gathered on Bitbucket.

mkdir Dev
cd Dev/
hg clone http://bitbucket.org/inglada/pwob

PWOB stands for Playing With OTB Bindings. I have only put there 3 languages which run on the JVM for reasons I stated in a previous post. You will of course avoid the plain Java ones. I have mixed feelings about Scala. I definetly love Clojure since it is a Lisp.

 

 

Get Clojure

The cool thing about this Lisp implementation is that it is contained into a jar file. You can get it with the text-based web browser links:

links http://clojure.org

Go to the download page and grab the latest release. It is a zip file which contains, among other things the needed jar. You can unzip the file:

mkdir src
mv clojure-1.3.0.zip src/
cd src/
unzip clojure-1.3.0.zip

Copy the jar file to the user .clojure dir:

cd clojure-1.3.0
mkdir ~/.clojure
mv clojure-1.3.0.jar ~/.clojure/

Make a sym link so we have a clojure.jar:

ln -s /home/inglada/.clojure/clojure-1.3.0.jar /home/inglada/.clojure/clojure.jar

And clean up useless things

cd ..
rm -rf clojure-1.3.0*

 

 

Final steps before hacking

Some final minor steps are needed before the fun starts. You may want to create a file for the REPL to store completions:

touch ~/.clj_completions

In order for Clojure to find all the jars and shared libraries, you have to define some environment variables. You may choose to set them into your .bashrc file:

export LD_LIBRARY_PATH="~/builds/OTB-Wrapping/lib/:~/builds/OTB/bin/"
export ITK_AUTOLOAD_PATH="~/builds/OTB/bin/"

PWOB provides a script to run a Clojure REPL with everything set up:

cd ~/Dev/pwob/Clojure/src/Pwob
./otb-repl.sh

Now you should see something like this:

Clojure 1.3.0
user=>

Welcome to Lisp!

If you want to use the applications you can for instance do:

(import '(org.otb.application Registry))
(def available-applications (Registry/GetAvailableApplications))

In the PWOB source tree you will find other examples.

 

 

 

 

Final remarks

Using the REPL is fun, but you will soon need to store the lines of code, test things, debug, etc. In this case, the best choice is to use SLIME, the Superior Lisp Interaction Mode for Emacs.

There are many tutorials on the net on how to set it up for Clojure. Search for it using DuckDuckGo.

In the PWOB tree (classpath.clj) you will find some hints on how to set it up for OTB and Clojure. A simpler config for Emacs is to use the inferior lisp mode, for which I have also written a config (otb-clojure-config.el).

I may write a post some day about that.

Have fun!

Lambda Omega Lambda

In my current research work, I am investigating the possibility of integrating domain expert knowledge together with image processing techniques (segmentation, classification) in order to provide accurate land cover maps from remote sensing image time series.

When we look at this in terms of tools (software) there are a set of requirements which have to be met in order to be able to go from toy problems to operational processing chains. These operational chains are needed by scientists for carbon and water cycle studies, climate change analysis, etc.

In the coming years, the image data needed for this kind of near-real-time mapping will be available through Earth observation missions like ESA’s Sentinel program.

Briefly, the constraints that the software tools for this task have are the following.

  1. Interactive environment for exploratory analysis
  2. Concurrent/parallel processing
  3. Availability of state of the art, validated image processing algorithms
  4. Symbolic, semantic information processing

This long post / short article describes some of the thoughts I have and some of the conclusions I have come to after a rather long analysis of the tools available.

Need for a repl

The repl acronym stands for Read-Evaluate-Print-Loop and I think has its origins in the first interactive Lisp systems. Nowadays, we call this an interactive interpreter. People used to numerical and/or scientific computing will think about Matlab, IDL, Octave, Scilab and what not. Computer scientists will see here the Python/Ruby/Tcl/etc. interpreters.

Many people use nowadays Pyhton+Numpy or Scipy (now gathered in Pylab) in order to have something similar to Matlab/IDL/Scilab/Octave, but with a general purpose programming language.

What a repl allows is to interactively explore the problem at hand without going through the cycle of writing the code, compiling, running etc.

The simple script that I wrote for the OTB blog using OTB’s Python bindings could be typed on the interpreter and bring together OTB and Pylab.

 

 

Need for concurrent programming

Without going into the details and differences between parallel and concurrent programming, it seems clear that Moore’s law can only continue to hold through multi-core architectures.

In terms of low-level (pixel-wise) processing, OTB provides an interesting solution (multi-threaded execution of filters by dividing images into chunks). This approach can be generalized to GPUs.

However, sometimes an algorithm needs to operate on the whole image because the image splitting affects the results. This is typically the case for Markovian approaches to filtering or methods for image segmentation.

For this cases, one way to speed up things is to process several images in parallel (if the memory footprint allows that!). On way of trying to maximize the use of all available computing cores in Python is using the multiprocessing module which allows to deal with a pool of threads. One example would be as follows:

However, this does not allow for easy inter-thread communication which is not needed in the above example, but can be very useful if the different processes are working on the same image: imagine a multi-agent system where classifiers, algorithms for biophysical parameter extraction, data assimilation techniques, etc. work together to produce an accurate land cover map. They may want to communicate in order to share information. As far as I understand, Python has some limitations due to the global interpreter lock. Some languages as for instance Erlang offer appropriate concurrency primitives for this. I have played a little bit with them in my exploration of 7 languages in 7 weeks. Unfortunately, there are no OTB bindings for Erlang. Scala has copied the Erlang actors, but I didn’t really got into the Scala thing.

 

 

Need for OTB access

 

This one here shouldn’t need much explanation. Efficient remote sensing image processing needs OTB. Period. I am sorry. I am rather biased on that! I like C++ and I have no problem in using it. But there is no repl, one needs several lines of typedef before being able to use anything. This is the price to pay in order to have a good static checking of the types before running the problem. And it’s damned fast! We have Python bindings which allows us to have a clean syntax, like in the pipeline example of the OTB tutorials. However, the lack of easy concurrency is a bad point for Python. Also, the lack of Artificial Intelligence frameworks for Python is an anti-feature. Java has them, but Java has no repl and look at its syntax. It’s worse than C++. You have all these mangled names which were clean in C++ and Python and become things like otbImageFileReaderIUS2. Scala, thanks to its interoperability with Java (Scala runs in the JVM), can use OTB bindings. Actually, we have a cleaner syntax than Java’s:There is still the problem of the mangled names, but with some pattern matching or case classes, this should disappear.

So Scala seems a good candidate. It has:

  1. Python-like syntax (although statically typed)
  2. Concurrency primitives
  3. A repl

Unfortunately, Scala is not a Lisp. Bear with me.

 

 

A Lisp would be nice

I want to build an expert system, it would be nice to have something for remote sensing similar to Wolfram Alpha. We could call it \Omega Toolbox and keep the OTB name (or close).
Why Lisp? Well I am not able to explain that here, but you can read P. Norvig’s or P. Graham’s essays on the topic.

If you have a look at books like PAIP or AIMA, or systems like LISA, CLIPS or JESS, they are either written in Lisp or the offer a Lisp-like DSLs.

I am aware of implementations of the AIMA code in Python, and even P. Norvig himself has reasons to have migrated from Lisp to Python, but as stated above, Python seems to be out of the game for me.

The code is data philosophy of Lisp is, as far as I understand it, together with the repl tool, one of the main assets for AI programming. Another aspect which is also important is the functional programming paradigm used in Lisp (even though other programming paradigms are also available in Lisp).

Concurrency is the main reason for the upheaval of functional languages in recent years (Haskell, for instance). Even though I (still) don’t see the need for pure functional programming for my applications, lambda calculus is elegant and interesting.

Maybe \lambda Toolbox should be a more appropriate name?

 

 

Clojure

If we recap the needs:

  1. A repl (no C++, no Java)
  2. Concurrency (no Python)
  3. OTB bindings available (no Erlang, no Haskell, no Ruby)
  4. Lisp (none of the above)

there is one single candidate: Clojure.

Clojure is a Lisp dialect which runs on the JVM and has nice concurrency features like inmutability and STM and agents.

And by the way, OTB bindings work like a charm:

Admittedly, the syntax is less beautiful than Python’s or Scala’s, but (because!) it’s a Lisp. And it’s a better Java than Java. So you have all the Java libs available, and even Clojure specific repositories like Clojars. A particular interesting project is Incanter which provides is a Clojure-based, R-like platform for statistical computing and graphics. Have a look at this presentation to get an overview of what you can do at the Clojure repl with that. If we bear in mind that in Lisp code is data and that Lisp macros are mega-powerful, one could imagine writing things like:

(make-otb-pipeline reader gradientFilter thresholdFilter writer)

Or even emulating the C++ template syntax to avoid using the mangled names of the OTB classes in the Java bindings (using macros and keywords):

(def filter (RescaleIntensityImageFilter :itk (Image. : otb :Float 2)
                                              (Image. : otb :UnsignedChar 2)))

instead of

(def filter (itkRescaleIntensityImageFilterIF2IUC2.))

I have already found a cool syntax for using the Setters and Getters of a filter using the doto macro (see line 19 in the example below):

Conclusion

I am going to push further the investigation of the use of Clojure because is seems to fit my needs:

  1. Has an interactive interpreter
  2. Access to OTB (through the Java bindings)
  3. Concurrency primitives (agents, STM, etc.)
  4. It’s a lisp, so I can easily port existing rule-based expert systems.

Given the fact that this would be the sum of many cool features, I think I should call it \Sigma Toolbox, but I don’t like the name.

The mix of \lambda calculus and \Omega Toolbox, should be called \lambda \Omega \lambda, which is LoL in Greek.

 

 

Geospatial analysis and Object-Based image analysis

I was searching in the web about the use of PostGIS data bases for object based image analysis (OBIA) and Google sent me to the OTB Wiki to a page that I wrote 6 months ago (first hit for “postgis object based image analysis”).
It seems that this is a subject for which no much work is ongoing. Julien Michel will be presenting some results about how to put together OBIA and geospatial analysis (+ active learning and some other cool things) in a workshop in Paris next May.

Changing the processing paradigm?

No I am not going to write about compressive sensing. I am sorry.

The paradigm I want to write about is the one currently used for mapping applications from remote sensing imagery. This way of turning images into maps can be resumed as ortho-analysis-GIS, that is, you get the image and convert it into a map projection, then you analyze it (segmentation, classification, etc.) in order to produce vector layers and finally you import these vector layers into a GIS in order to perform geospatial analysis (eventually fusing with existing vector maps) and you produce the final map.

This works OK, but is not the most efficient way of working. If you look at the final map, how much information really came from the images? How may pixels were really useful for producing this information? The answer is usually “not much”.

Now look at the computation time needed by all the processing steps in the map production chain. Can you guess what is the most expensive step? With current high resolution images, the ortho-rectification step is the most time consuming.

One solution to this could be to ortho-rectify only the interesting area for the application. The drawback of this approach is that usually, you need to process the image (detect the pertinent features, changes, etc.) before you know where is the interesting information.

In this case, the solution, would be to process the image before ortho-rectification. There is one main problem to this : many modern software tools for segmentation and classification are not very good at processing huge images. Even if they were good at it, you would need to process the whole image before ortho-rectification.

The thing is that often, the existing maps tell you where the interesting things are likely to be found, but since your maps are “ortho”, you still need to ortho-rectify your image before processing.

Also, the geospatial reasoning step is made at the end of the processing, inside the GIS tool, which usually knows very little about image processing and so.

So it seems that the paradigm cited above (which could also be named ERDAS-Definiens-ArcGIS, for example), although useful has real drawbacks for efficiency. And I am not talking about import/export and format issues.

In order to be really efficient, we would need a tool which would allow us to send the existing shapefile or KML maps on top of the image in sensor geometry, perform some geospatial reasoning up there, segment and classify only the areas of interest (still in sensor geometry), produce vector data and finally send only the useful information down to the map projection.

Hey, I finally nearly wrote about compressive sensing, didn’t I?

To finnish: don’t tell anybody, but it seems there is a free software out there which is able to do what I have just wrote. Well the PostGIS interface is not yet ready, but it is on its way.

The Wiki Power

Last week I could experiment the power of cooperative work on the Internet.

In the frame of the development of the ORFEO Toolbox (OTB), I was doing some bibliographic research. I wanted to find as many radiometric indices as possible that I could compute using optical multispectral remote sensing data. My problem was mainly focussed on Pleiades-like data, but I thought that it would be interesting to widen the search. After all, these indices are going to be coded in OTB, which has a rather heterogeneous user base.

Unfortunately, I have no easy access to application oriented remote sensing journals, so I was facing a tedious work.

I started with the most well known indices as NDVI and I wrote some descriptions for them in the OTB Wiki. As I was on the wiki, I thought I could ask some other people to complete the documentation, so I sent an e-mail to several mailing lists where I knew that people knew about radiometry issues.

To be honest, at the beginning I thought that it would be useless. I had a good surprise, when the next morning I saw that somebody had created a new account on the wiki and was formatting and enriching (adding some descriptions and bibliographic references) to the list of indices I had started. The day after that 2 other people joined the effort and added new indices and useful information.

The list is accessible here. At the moment of this writing, the list seems to be stable and contains 14 vegetation indices, 8 water indices, 3 soil indices and 2 built-up indices. Most of them where immediately coded in OTB.

This is a real win-win approach: you give us the formulas, we give you the code.

Pragmatic programmers

This is a blog on remote sensing image processing, not a blog on programming or computers. However, if you want to do efficient remote sensing image processing, you are often led to do some programming yourself.

Programming is a thing that I like very much. High quality, bug free, efficient programs is something I like to use. Since I am using many of my own programs, I try to use good approaches for this.

Of course, ORFEO Toolbox is my library of choice. It relies on many interesting concepts as object orientation, generic programming, streaming, multi-threading, etc. These concepts invite the programmer to write good code, but they may not be enough. Many tips and techniques have to be applied in order to be happy with your code at the end of the day.

I have found 2 interesting bibliographic references from the Pragmatic Programmer which are worth reading:

1. Practices of an Agile Developer
2. Ship it!: A Practical Guide to Successful Software Projects

These 2 books are very easy (and fun!) to read and give many hints on how to proceed to produce good code.

ORFEO Toolbox (OTB) on YouTube

Emmanuel Christophe has produced a very interesting video about ORFEO Toolbox.

The animation has been generated using an open source tool called Code Swarm http://vis.cs.ucdavis.edu/~ogawa/codeswarm/. You can watch it here:


This video shows OTB’s history in terms of change in source code. Each dot represents a file gravitating around the developer who edited it.

It is interesting to note the increase of activity during the weeks before a release and the decrease during summer holidays (the development team lives in France …).

Anybody wants to join the adventure?

Summer School in Very High Resolution Remote Sensing

You may be interested in the upcoming 2008 Summer School on Very High
Resolution Remote Sensing which will be held in Grenoble, France, over a
week (September 8-12).

I will be giving a talk about “New sensors, new missions and new challenges on very high resolution remote sensing imagery”.

It is a pleasure for me having been invited by the organizing committee led by Jocelyn Chanussot to participate with international experts as:

  • Lorenzo Bruzzone, University of Trento, Italy
  • Paolo Gamba, University of Pavia, Italy
  • Antonio Plaza, University of Extremadura, Spain
  • Kostas Papathanassiou, DLR, Germany
  • Irena Hajnsek, DLR, Germany

CNES will take care of organizing lab sessions using the ORFEO Toolbox
(an open source library for remote sensing image processing algorithms developed by CNES).

You can find more information about the program and the registration procedure here.