<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">
<HTML><HEAD>
<META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=iso-8859-1">
<TITLE></TITLE>
<META content="MSHTML 6.00.2900.3059" name=GENERATOR></HEAD>
<BODY text=#000000 bgColor=#ffffff>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>Hi Dan,</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>I think the problem could be the amount of memory the
filter is consuming in the process. For an
output with PixelType float the output data in the first case is
consuming :</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007></SPAN><SPAN
class=281511214-10052007></SPAN><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>128x128x128 tensors x 6 pixels/tensor x 2 bytes /
pixel = 25.165.824 bytes</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>next case is 8 times the number of pixels so you
have</FONT></SPAN></DIV>
<DIV> </DIV>
<DIV><FONT face=Arial color=#0000ff size=2></FONT>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>256x256x256 tensors x 6 pixels/tensor x 2 bytes /
pixel = 201.326.592 bytes</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>512x512x512 tensors x 6 pixels/tensor x 2 bytes /
pixel = 1.610.612.736 bytes</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>but probably the filter is consuming more memory in the
process itself. So maybe in the second case you have already run out of physical
memory. If you are using double precision then multiply these numbers by 2 so in
the second case only with the output you are consuming 400 Mb. I would check out
this first before considering a problem in the filter
itself.</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>Answering your second questions, for computing gaussian
derivatives probably recursive filters are among the fastest.
Implementation in ITK follows Deriche's work :</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2><A
href="http://citeseer.ist.psu.edu/deriche93recursively.html">http://citeseer.ist.psu.edu/deriche93recursively.html</A></FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>You may find useful information and references
here</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2><A
href="http://en.wikipedia.org/wiki/Scale_space_implementation">http://en.wikipedia.org/wiki/Scale_space_implementation</A></FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>Regards</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2>Ivan Macía</FONT></SPAN></DIV>
<DIV dir=ltr align=left><SPAN class=281511214-10052007><FONT face=Arial
color=#0000ff size=2></FONT></SPAN> </DIV>
<DIV dir=ltr align=left>
<HR tabIndex=-1>
<FONT face=Tahoma size=2><B>De:</B>
insight-users-bounces+imacia=vicomtech.es@itk.org
[mailto:insight-users-bounces+imacia=vicomtech.es@itk.org] <B>En nombre de
</B>Dan Mueller<BR><B>Enviado el:</B> jueves, 10 de mayo de 2007
3:53<BR><B>Para:</B> insight-users<BR><B>Asunto:</B> [Insight-users]
HessianRecursiveGaussian execution speed<BR></FONT><BR></DIV></DIV>
<DIV></DIV><FONT face="Helvetica, Arial, sans-serif"><SMALL>Hi
insight-users,<BR><BR>I am using the <FONT
face="Courier New, Courier, monospace">itk::HessianRecursiveGaussianImageFilter</FONT>
and have noticed it executes quite slow for 'large' 3-D datasets. Here are some
times I recorded:<BR><BR>itkTestHessian TestHessian_128x128x128.mhd 1 1.0
>> Time = 4.28sec<BR></SMALL></FONT><FONT
face="Helvetica, Arial, sans-serif"><SMALL>itkTestHessian
TestHessian_256x256x256.mhd 1 1.0 >> Time = 6.2e+003sec = 1hr
43min<BR></SMALL></FONT><FONT
face="Helvetica, Arial, sans-serif"><SMALL>itkTestHessian
TestHessian_512x512x512.mhd 1 1.0 >> I gave up...</SMALL></FONT><BR><FONT
face="Helvetica, Arial, sans-serif"><SMALL><BR><B>My questions:</B><BR>1. Have
other people experienced similar execution times using this filter, or am I
doing something wrong? I understand the filter performs separable convolution
multiple times in multiple dimensions (ie. it's not a simple computation), but
these times seem impractical for datasets 256x256x256 or bigger...<BR>2. Are
there any techniques (besides relying on Moore's Law) to speed up this filter? I
am interested in the <I>whole</I> image, so using an image function would not
help. </SMALL></FONT><FONT face="Helvetica, Arial, sans-serif"><SMALL>(NOTE: I
compiled my test in <FONT face="Courier New, Courier, monospace">Release</FONT>
mode with optimization for maximum speed ie. <FONT
face="Courier New, Courier, monospace">/O2</FONT>).</SMALL></FONT><BR><FONT
face="Helvetica, Arial, sans-serif"><SMALL><BR>Thanks for any help and/or
ideas.<BR><BR>Cheers, Dan<BR><A class=moz-txt-link-abbreviated
href="mailto:d.mueller@qut.edu.au">d.mueller@qut.edu.au</A><BR><B><BR></B></SMALL></FONT><SMALL><FONT
face="Helvetica, Arial, sans-serif"><B>My
system:</B><BR></FONT></SMALL><SMALL><FONT
face="Helvetica, Arial, sans-serif">Platform: Windows XP SP2<BR>Computer: Intel
Pentium 3.00 GHz, dual core, 1GB physical RAM<BR>Compiler: Microsoft Visual
Studio 2005 (v8.0.50727.762 SP1)<BR>CMake: 2.4.2<BR>ITK: 3.2.0 CVS on Tues 27
Mar 2007</FONT></SMALL><BR><FONT
face="Helvetica, Arial, sans-serif"><SMALL><BR><B>My test code:</B><BR><FONT
face="Courier New, Courier, monospace">/*=========================================================================<BR>
itkTestHessian.cxx<BR>=========================================================================*/<BR><BR>#include
<iomanip><BR>#include "itkImage.h"<BR>#include
"itkNumericTraits.h"<BR>#include "itkCommand.h"<BR>#include
"itkImageFileReader.h"<BR>#include "itkImageFileWriter.h"<BR>#include
"itkHessianRecursiveGaussianImageFilter.h"<BR>#include
"itkTimeProbe.h"<BR><BR>// Declare general types<BR>const unsigned int Dimension
= 3;<BR>typedef float PixelType;<BR>typedef itk::Image< PixelType, Dimension
> ImageType;<BR>typedef itk::ImageFileReader< ImageType >
ReaderType;<BR>typedef itk::ImageFileWriter< ImageType >
WriterType;<BR>typedef itk::HessianRecursiveGaussianImageFilter< ImageType
> HessianFilterType;<BR><BR>// Declare a command to display the progress in a
nice block format<BR>class ProgressCommand : public itk::Command
<BR>{<BR>public:<BR> typedef ProgressCommand
Self;<BR> typedef
itk::Command
Superclass;<BR> typedef itk::SmartPointer<Self>
Pointer;<BR> itkNewMacro( Self );<BR><BR>protected:<BR>
ProgressCommand() { m_LastProgress = -1; };<BR> float
m_LastProgress;<BR><BR>public:<BR><BR> void Execute(itk::Object *caller,
const itk::EventObject & event)<BR> {<BR> Execute(
(const itk::Object *)caller, event);<BR> }<BR><BR> void
Execute(const itk::Object * object, const itk::EventObject &
event)<BR> {<BR> const itk::ProcessObject*
process = dynamic_cast< const itk::ProcessObject* >( object
);<BR><BR> if( ! itk::ProgressEvent().CheckEvent(
&event ) )<BR>
return;<BR><BR> int fprogress =
(process->GetProgress() * 100.0);<BR> int progress =
(int)(process->GetProgress() * 100.0);<BR> if
((int)m_LastProgress == progress)<BR>
return;<BR> if ((int)m_LastProgress != (progress
- 1))<BR> {<BR>
std::cout << std::setfill('0') << std::setw(3) <<
(progress - 1) << " ";<BR>
if (fprogress > 0.0 && (progress - 1)
% 10 == 0)<BR> std::cout
<< std::endl;<BR> }<BR>
if (fprogress > 0.0 && fprogress <= 100.0)<BR>
std::cout << std::setfill('0') <<
std::setw(3) << progress << " ";<BR> if
(fprogress > 0.0 && progress % 10 == 0)<BR>
std::cout << std::endl;<BR>
m_LastProgress = fprogress;<BR> } <BR>};<BR><BR>int
main(int argc, char* argv[])<BR>{<BR> // Display
header<BR> std::cout << "Hessian Test
-------------------------" << std::endl;<BR>
<BR> // Check arguments<BR> if (argc <
4)<BR> {<BR> std::cerr
<< "Usage: " << argv[0] <<<BR>
" InputFilename"
<<<BR>
" NormalizeAcrossScale" <<<BR>
" Sigma"
<<<BR>
" \n";<BR> return
EXIT_FAILURE;<BR> }<BR><BR> // Setup the
algorithm parameters<BR> unsigned int argn =
1;<BR> char*
InputFilename =
argv[argn++];<BR> bool NormalizeAcrossScale =
(bool)atoi( argv[argn++] );<BR> double
Sigma
= atof( argv[argn++] );<BR><BR> // Display
parameters<BR> std::cout << "InputFilename:" <<
InputFilename << std::endl;<BR> std::cout <<
"NormalizeAcrossScale:" << NormalizeAcrossScale <<
std::endl;<BR> std::cout << "Sigma:" << Sigma
<< std::endl;<BR> std::cout <<
"------------------------------------" << std::endl;<BR>
<BR> try<BR> {<BR>
// Read input image<BR>
std::cout << "Reading input..." << std::endl;<BR>
ReaderType::Pointer reader =
ReaderType::New();<BR>
reader->SetFileName( InputFilename );<BR>
reader->Update();<BR><BR>
// Setup Hessian<BR>
std::cout << "Computing Hessian..." <<
std::endl;<BR> HessianFilterType::Pointer
filterHessian = HessianFilterType::New();<BR>
filterHessian->SetInput( reader->GetOutput()
);<BR> filterHessian->SetSigma( Sigma
);<BR>
filterHessian->SetNormalizeAcrossScale( NormalizeAcrossScale
);<BR> ProgressCommand::Pointer
observerHessian = ProgressCommand::New();<BR>
filterHessian->AddObserver( itk::ProgressEvent(),
observerHessian );<BR>
<BR> // Compute and time
hessian<BR> itk::TimeProbe
time;<BR>
time.Start();<BR> filterHessian->Update(
);<BR> time.Stop();<BR>
std::cout << std::setprecision(3) << "Time: "
<< time.GetMeanTime() << std::endl;<BR>
}<BR><BR> catch (itk::ExceptionObject &
err)<BR> { <BR> std::cout
<< "ExceptionObject caught !" << std::endl; <BR>
std::cout << err << std::endl;
<BR> return
EXIT_FAILURE;<BR> }<BR><BR>
//Return<BR> return
EXIT_SUCCESS;<BR>}</FONT><BR></SMALL></FONT><BR>
<P><FONT size=2>No virus found in this incoming message.<BR>Checked by AVG Free
Edition.<BR>Version: 7.5.467 / Virus Database: 269.6.6/794 - Release Date:
08/05/2007 14:23<BR></FONT></P>
<P><FONT face=Arial size=2></FONT></P></BODY></HTML>
<BR>
<P><FONT SIZE=2>No virus found in this outgoing message.<BR>
Checked by AVG Free Edition.<BR>
Version: 7.5.467 / Virus Database: 269.6.6/795 - Release Date: 09/05/2007 15:07<BR>
</FONT> </P>