[Insight-users] Memory management

Bradley Lowekamp blowekamp at mail.nih.gov
Mon Mar 21 12:03:24 EDT 2011


Hello Kana,

I would strongly not recommend not taking this approach. It is rather inefficient, and version 3.20 has a bug related to this method.

I have an alternative implementation with exploits the separability of the convolution by a derivative. Simply, I implemented a discrete Hessian filter where a Gaussian convolution can occur as a prior filter. I hopefully will be able to get this to a sharable state this week.

Brad



On Mar 21, 2011, at 10:35 AM, Arunachalam Kana wrote:

> Hi Vera,
> Thank you for the information that the recursive Gaussian cannot be streamed.
>  
> Hi Masslawi,
> I have not tried the release version , i will surely try it.
>  
> I decided to use discretehessiancalculation for calculation of hessian matrix. I thought of increasing the speed by multithreading in my own filter.
> My filter takes 1 input image and gives 1 output image. So i calculate the hessian for each voxel.
> The problem: My program breaks due to the following error.
> ERROR MESSAGE:
> Windows has triggered a breakpoint in iAnalyseGUI.exe.
> This may be due to a corruption of the heap, which indicates a bug in iAnalyseGUI.exe or any of the DLLs it has loaded.
> This may also be due to the user pressing F12 while iAnalyseGUI.exe has focus.
> The output window may have more diagnostic information.
>  
> The breakage happens when the 3rd thread starts, and in DiscreteHessianGaussianImageFunction.txx line 116.
> I need some help in solving the heap problem.
>  
> And finally the question comes to whether is it possible to multithread the DiscreteHessianGaussianImageFunction.txx ?????
> As i have to work with large size data, i really need it, so if i can find the code somewhere please direct me or i would need some help to
> Multithread the DiscreteHessianGaussianImageFunction.txx.
>  
> //code
> template<class TInputImageType,  class TOutputImageType>
> void MYFilter<TInputImageType,     TOutputImageType>
> ::ThreadedGenerateData(const OutputImageRegionType& outputRegionForThread, int threadId)
> {
>   typename OutputImageType::Pointer output = this->GetOutput();
>   typename InputImageType::ConstPointer input = this->GetInput();
>  
>                 //set hessian function
> m_HessianFunction = HessianFunctionType::New();
>                 m_HessianFunction->SetInputImage(m_GradientMagnitudeImage );
>                 m_HessianFunction->SetSigma( 1 );
>                 m_HessianFunction->SetMaximumError( 0.01 );
>                 m_HessianFunction->SetMaximumKernelWidth( 32 );
>                 m_HessianFunction->SetNormalizeAcrossScale( 1 );
>                 m_HessianFunction->SetUseImageSpacing( false );
>                 m_HessianFunction->Initialize( );
>                
>                 //initiate iamge iterators
>                 InputImageRegionIteratorType inIter(m_InputImage, outputRegionForThread);
>                 OutputImageRegionIteratorType outIter(m_OutputImage, outputRegionForThread);
>  
>                 inIter.GoToBegin();
>                 while ( !inIter.IsAtEnd() )
>                 {
>                                InputIndexType currentindex =inIter.GetIndex();
>                                InputPixelType currentpixel = inIter.Get();
>                                               
>                                //initialise the variables
>                                HessianFunctionType::TensorType   hess_matrix;
>                                hess_matrix = m_HessianFunction->EvaluateAtIndex( currentindex );
>                               
>                                ++inIter;                            
>                 }//while
> }
>  
> Thank you,
> Regards,
> Kana
>  
>  
> From: Dawood Masslawi [mailto:masslawi at gmail.com] 
> Sent: 19 March 2011 14:48
> To: Arunachalam Kana
> Cc: insight-users at itk.org
> Subject: RE: Memory management
>  
> Hi Kana,
>  
> Have you tried to build in release mode? It can make a big difference.
>  
> Good luck,
>  
> Dawood
>  
>  
> >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
>  
> >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 
>  
> > Hi ITK Users,
> > 
> > I need some advice on memory management and speed.
> > 
> > My system details: 64 bit windowsXP system and 4Gb RAM.
> > 
> > My goal: To obtain eigen vectors for 100Mb unsigned short image. Eigen
> > vectors not for the whole image but for certain regions (approx occupies 40%
> > of image).
> > 
> > 
> > 
> > I have to do some preprocessing too.
> > 
> > The pipeline is :
> > 
> > Original image -> discretegaussian -> gradientmagnitude ->
> > HessianRecursiveGaussianImageFilter->symmetric eigenvalue analysis ->
> > finally 3 vector image
> > 
> > For the above pipeline for a 100Mb image, i am running out of RAM.
> > 
> > After reading itkmails, i decided to use streamfilter to stream part of the
> > image so that the RAM usage is low. I only tested half way and the RAM was
> > overloaded.
> > 
> > 
> > 
> > I am using streamimagefilter to reduce the memory usage. After calculation
> > of gradientmagnitude i have used 1.08Gb RAM. After this i need 1.2Gb RAM
> > 
> > For hessianimageoutput, so by using streamfilter i thought i would be
> > reaching 2.5Gm RAM but i reach approx 6Gb(virtual memory is used).
> > 
> > 1.       Am i doing something wrong with the stream filter ? (code given
> > below)
> > 
> > 2.       I tried another option: As i do not need hessian for the whole
> > image, I used discretehessianfunction for hessian calculation but it is 30
> > times slower that the recursivehessian. Is discretegaussianimagefunction
> > multithreaded?
> > 
> > 3.       Is there any other way to achieve speed vs RAM compromised
> > solution?
> > 
> > 4.       The 100Mb data set is test data. The real data is 10Gb for which
> > i will use 64 bit linux system (opensuse 11.2) . will the pipeline be
> > executed in linux too?
> > 
> > Below is my code and after application of filter show the RAM usage:
> > 
> > 
> > 
> >   char *infilename = "StreamTest.mhd";
> > 
> >   typedef itk::CovariantVector<float,3> VectorPixelType;
> > 
> >   typedef itk::Vector< VectorPixelType, 3 > EV_PixelType;
> > 
> > 
> > 
> >   //Image type
> > 
> >   typedef itk::Image<unsigned short, 3>              InputImageType;
> > 
> >   typedef itk::Image<float, 3>                       FloatImageType;
> > 
> >   typedef
> > itk::Image<EV_PixelType,3>
> > EVImageType;
> > 
> > **** RAM usage = 674Mb****
> > 
> >   //reader initialisation and reading file
> > 
> >   typedef itk::ImageFileReader<InputImageType>
> > ReaderType;
> > 
> >   ReaderType::Pointer reader = ReaderType::New();
> > 
> >   reader->SetFileName( infilename );
> > 
> >   reader->Update();
> > 
> >   InputImageType::Pointer inImage = reader->GetOutput();
> > 
> > 
> > 
> > **** RAM usage = 783Mb****
> > 
> >   typedef itk::DiscreteGaussianImageFilter<InputImageType, FloatImageType>
> > RGIFType;
> > 
> >   RGIFType::Pointer gaussfilter = RGIFType::New();
> > 
> >   gaussfilter->SetInput(inImage);
> > 
> >   gaussfilter->SetVariance(4.0);
> > 
> >   gaussfilter->SetMaximumError(0.01);
> > 
> > 
> > 
> >   typedef itk::CastImageFilter<FloatImageType,InputImageType> CIFType;
> > 
> >   CIFType::Pointer castfilter = CIFType::New();
> > 
> >   castfilter->SetInput(gaussfilter->GetOutput());
> > 
> > 
> > 
> >   typedef itk::GradientMagnitudeImageFilter<InputImageType,FloatImageType>
> > GMIFType;
> > 
> >   GMIFType::Pointer gmfilter = GMIFType::New();
> > 
> >   gmfilter->SetInput(castfilter->GetOutput());
> > 
> >   gmfilter->Update();
> > 
> > 
> > 
> >   FloatImageType::Pointer gmImage = gmfilter->GetOutput();
> > 
> > **** RAM usage = 1290Mb****
> > 
> > 
> > 
> >   gaussfilter->UnRegister();
> > 
> >   castfilter->UnRegister();
> > 
> >   **** RAM usage = 1080Mb****
> > 
> > 
> > 
> >   typedef itk::HessianRecursiveGaussianImageFilter<FloatImageType>
> > HessianRecursiveGaussianFilterType;
> > 
> >   typedef HessianRecursiveGaussianFilterType::OutputImageType
> > HessianImageType;
> > 
> >   HessianRecursiveGaussianFilterType::Pointer hessianfilter =
> > HessianRecursiveGaussianFilterType::New();
> > 
> >   hessianfilter->SetInput(gmImage);
> > 
> >   hessianfilter->SetSigma(1);
> > 
> > 
> > 
> >   typedef itk::StreamingImageFilter<HessianImageType,HessianImageType>
> > StreamerType;
> > 
> >   StreamerType::Pointer streamer = StreamerType::New();
> > 
> >   streamer->SetInput( hessianfilter->GetOutput() );
> > 
> >   streamer->SetNumberOfStreamDivisions( 20 );
> > 
> >   streamer->Update();
> > 
> > **** RAM usage = 3600Mb**** after few seconds it came to **** RAM usage =
> > 6120Mb****
> > 
> > 
> > 
> > Thank you in advance.
> > 
> > Regards,
> > 
> > Kana Arunachalam Kannappan
> > 
> > Research Associate
> > 
> > FH OÖ Forschungs & Entwicklungs GmbH
> > 
> > Stelzhamer Strasse 23,
> > 
> > 4600 Wels,
> > 
> > Austria.
> > 
> > Phone: +43 (0)7242 72811 -4420
> > 
> > kana.arunachalam at fh-wels.at
> > 
> > www.fh-ooe.at; www.3dct.at
>  
> <ATT00001..txt>

========================================================
Bradley Lowekamp  
Lockheed Martin Contractor for
Office of High Performance Computing and Communications
National Library of Medicine 
blowekamp at mail.nih.gov


-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.itk.org/pipermail/insight-users/attachments/20110321/afc4770e/attachment.htm>


More information about the Insight-users mailing list