[Insight-users] how to ensure multi threading and best iterators

[devieill at irit.fr]

Hi all,
I am trying to issue some trivial bench between blitz++ and itk since
I wrote some code (optimization stuff : deconvolution RL-TV with non
linear Conjugate Gradient method) and that my program is rather slow for
26M 3D images..
So my question is  : assuming i have the following declarations and
directives in a main.cxx file (see after), how can I be sure that what is
inside the for statement will be executed with all the cores my computer
have ?
Although some cryptic compile options (e.g. -fununroll-loop
-ftreevectorize) exist on gcc 4.4 I wonder if itk uses directly the
multi-threading or if there exist a specific directive to enable it.

On the other side I wonder if what I wrote is the fastest way to achieve a
it3 = it1*log(1 + it2) computation, but to my knowledge there is no "array
operations" in itk (maybe I am wrong or they are disguised).

I would appreciate any comments/help regarding the multi-threading problem
and the fast traversal of the 3D images.

with kind regards,
de Vieilleville François


main.cxx
..
typedef itk::Image<double, 3> WorkingImageType;
typedef itk::ImageRegionIterator< WorkingImageType > ImageIteratorType;
..
int main () {
..
WorkingImageType::Pointer img1;
WorkingImageType::Pointer img2;
WorkingImageType::Pointer img3;
..
// initialisation & fill of the 3 3D images
..
ImageIteratorType it1( img1, img1->GetLargestPossibleRegion());
ImageIteratorType it2( img2, img2->GetLargestPossibleRegion());
ImageIteratorType it3( img3, img3->GetLargestPossibleRegion());
..
for ( it1.GoToBegin(), it2.GoToBegin(), it3.GoToBegin(); !it3.IsAtEnd();
++it1, ++it2, ++it3)
{
      it3.Set( it1.Get()*log(1.0+it2.Get()));
}
..
}