For this post, we are glad to have Thierry as a guest blogger. He talks about an optimisation example he recently encountered that is a good indirect example of ‘inlining‘.
As my first post here, I will first describe a standard “auto-contrast” method with some optimizing hints. I will then try a more powerful one, and hit a situation where the buit-in Matlab function (quantile) can be beaten by a simple Matlab script.
Suppose you want to “auto-contrast” an image.
A simple way to do this is to scale the value of your matrix from 0 to 255 (or [0..2^16-1] for a 16 bit image). You can make it with this code :
imageMin=min(startImage(:)); imageMax=max(startImage(:)); scaledImage = uint8(255*(startImage-imageMin)/(imageMax-imageMin))
Please note the “min(var(:))” where “(:)” linearize the 2D (or more) matrix, avoiding the slower solution “min(min(var))”.
Of course, this code makes more sense if you start from a matrix of double and want to store it in an image file (TIFF, png…). I am sure many of you know that the “imagesc” Matlab function is doing it for you, in the background.
With such a formula, you are sensitive to “hot pixels” or dust in optic which are out of scale compared to the rest of your data. To further improve the result, you therefore probably want to “saturate” both the highest (hot) pixels values as well as the lowest. This is what the “clims” parameter of imagesc is for. In ImageJ, there is a toot that “enhances contrast with XX% saturation“.
To perform the same calculation in Matlab, we just have to modify the min and max value in the previous formula by the XX% lowest and highest values. We might get these value with this algorithm :
sortedValues=sort(myImage(:)); N=length(sortedValues); satMin=sortedValues(floor(N*saturation)); satMax=sortedValues(floor(N*(1-saturation))); [and now, you use the previous autocontrast/scaling formula with this "satMin" and satMax values]
saturation is here a number between 0 and 1. I often choose a saturation of 0.001, which means that 0.1% of all my pixels will have the 0 value and the same amount will have 255.
Please note again the “(:)” to linearize myImage, allowing to sort all the value at the same time.
The above code happened to be what I produced as my first (after some debugging) guess. Using the profiler, I got a noticeable long calculation time on the “sort(” line. So I searched for a dedicated Matlab function (which in principle should have been the first thing to do). And guess what?
There is a “quantile” function which gives the same result.
First thing to be proud of: I reproduce a Matlab function with no bug!
This is the now the equivalent code using that quantile function.
I also found that the quantile function has a strange scaling with various input parameters, in term of computational speed, as shown below.
tic; for i=1:220,quantile(I(:),.9);end; toc
elapsed time : 11sec
tic; for i=1:220,quantile(I(:),[.1 .9]);end; toc
elapsed time : 19 sec
My small algorithm with “sort” does the job in 6 seconds, no matter how many values I want to extract (sorting the image is the bottleneck; extracting the “n-th” value is instantaneous).
To conclude, I used the profiler to find the line which takes the longest to compute. Usually, one can improve the speed by a factor (between 2 and 10) by removing bottlenecks. It was not possible here as even the builtin Matlab function was slower. The lesson to this is the following : if the slowest thing in your code is a basic Matlab function like matrix multiplication, convolution or (in this example) sort, there is usually nothing to do.
I’m sure you’ve already encountered such a situation, feel free to write it down here as a comment. As you can see, by writing simpler function, you can even beat Matlab provided function with a Matlab script of your own.