It seems that C++ programming can be cumbersome for some users. Well, OTB allows you to work also in other programming languages as for instance Python. The Pylab Project which gathers numeric and scientific Python packages such as NumPy, SciPy and Matplotlib is considered as having the same features as some commercial scientific environments. Some people consider that Pylab is even more powerful than these commercial tools, since you get for free the power of a full featured and general purpose programming language as Python.
Recently I needed to quickly plot the temporal profiles of several pixels coming from temporal image series in order to easily assess the impact of a correction algorithm. I have to admit that my first idea was to generate simple ASCII files with the pixel values and then feed these files to Gnuplot. However, I found a more elegant (at least from my point of view) way to do this.
I use the OTB Python bindings in order to read my images (2 GeoTiff files with 49 bands each). Using the “GetPixel” method, I access the pixels I want to plot. Finally I use pylab to plot 2 temporal profiles (before and after correction). Here goes the code.
import sys if( len(sys.argv) != 5): print("Usage: "+sys.argv[0]+" series1 series2 x0 y0") exit() import itk import otbImages import otbIO import pylab # Types ImageType = otbImages.VectorImage[itk.D, 2] ReaderType = otbIO.ImageFileReader[ImageType] # Pixel coordinates index = itk.Index[2]() index[0] = int(sys.argv[3]) index[1] = int(sys.argv[4]) # Import first series readerBefore = ReaderType.New() readerBefore.SetFileName(sys.argv[1]) readerBefore.Update() pixBefore = readerBefore.GetOutput().GetPixel( index ) # Import second series readerAfter = ReaderType.New() readerAfter.SetFileName(sys.argv[2]) readerAfter.Update() pixAfter = readerAfter.GetOutput().GetPixel( index ) # Chek that both series have the same length if( pixBefore.GetSize() != pixAfter.GetSize() ): print("Series have different lengths") exit() # Plot vectorBefore = [pixBefore.GetElement(i) for i in range(pixBefore.GetSize())] vectorAfter = [pixAfter.GetElement(i) for i in range(pixAfter.GetSize())] pylab.xlabel("Date") pylab.ylabel("Reflectance") pylab.plot(vectorBefore, label="Raw") pylab.plot(vectorAfter, label="Clean") pylab.legend(loc='upper left', shadow=True) pylab.show()
And here you can see a result example.
Temporal profiles
Of course, you can imagine to make this script a little bit more interactive in order to choose other pixel coordinates, etc., but I hope you grasp the power of this approach. By using the OTB Python bindings and other freely available open source Python goodies (NumPy, SciPy, Matplotlib, rpy) you get the user friendly approach of Octave, Scilab and their commercial counterparts plus the power of OTB (which in turn makes available ITK, OSSIM, etc.).
Isn’t it great?
Welcome to my favourite language admittedly I do the same with gdal python bindings … Try orange classifiers and learning stuff and hook up otb to it.
Yay for big snakes
give me direction , how i start pyhton otb wrapping .
Great example
I would also like to peform OTB processing tasks using Python as I am no software engineer nor have any experience in C++. I do however have experience using python for geo-processing tasks.
I am forced to use windows and this is the OS for most of our geoprocessing software
I have downloaded the OTB wrapping package windows installer and am using python 2.7. When I run the OTB python “hello world” script there is no recognition of the otb modules. “ImportError: No module named itk”
Not sure how to get python to recognise the OTB modules.
Any help would be greatly appeciated
Dan
Is it possible to load numpy array (image) and use it with other filters?
For example, I have image of instance python.pillow.Image. I can get numpy array version easily. Now if i want to use this python object as an image source just like a gdalIO/JpegIO. I am asking because itk gives numpy array and not sure if possible with OTB?