<html><head><meta http-equiv="Content-Type" content="text/html charset=iso-8859-1"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; -webkit-line-break: after-white-space; ">Hello lien,<div><br></div><div>Are you comparing 4.2.1 to 3.20? Or the v3 vs v4 registration in 4.2.1?</div><div><br></div><div>Also how does v3 compare to v4 when only a single thread is used?</div><div><br></div><div>Thanks for the additional information.</div><div><br></div><div>Brad</div><div><br><div><div>On Oct 26, 2012, at 9:12 AM, lien lee &lt;<a href="mailto:lienlee@gmail.com">lienlee@gmail.com</a>&gt; wrote:</div><br class="Apple-interchange-newline"><blockquote type="cite">Thanks for the informative reply, Brian. I am going to double check, and will let you know what I could get.<br><br><br><div class="gmail_quote">2012/10/25 brian avants <span dir="ltr">&lt;<a href="mailto:stnava@gmail.com" target="_blank">stnava@gmail.com</a>&gt;</span><br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">hi lien<br>
<br>
thanks for trying out v4.<br>
<br>
i would not give up on it yet as there are many ways to use the v4 framework.<br>
<br>
i don't know off hand which is the most computationally efficient<br>
approach --- to some extent we rely on users to help us understand<br>
these questions as the number of active developers, right now, is<br>
relatively few.<br>
<br>
typically, if we have questions about efficiency , we use a profiler<br>
to help understand sources of additional computation.<br>
<br>
michael stauffer has done the majority of profiling - and he is cc'd here.<br>
<br>
from what i recall, the v4 and v3 metrics are roughly similar in terms<br>
of speed and per iteration cost for computing the metric,<br>
<br>
assuming that you set up the metric in a similar manner.<br>
<br>
so , without seeing your v3 code, &nbsp;we cannot understand what<br>
differences might exist in your v4 and v3 comparison. &nbsp;one obvious<br>
question is per iteration cost of assessing the metric and how many<br>
points one is using in the call<br>
<br>
vRigidRegistration-&gt;SetMetricSamplingPercentage(0.1);<br>
<br>
relative to your v3 setup.<br>
<br>
the "extra" computation from the virtual domain and in the composite<br>
transform should be relatively little as the identity transform does<br>
not require any actual computation.<br>
<br>
just a few thoughts ---- perhaps others will have additional feedback.<br>
<br>
brian<br>
<div><div class="h5"><br>
<br>
<br>
<br>
On Thu, Oct 25, 2012 at 5:20 PM, lien lee &lt;<a href="mailto:lienlee@gmail.com">lienlee@gmail.com</a>&gt; wrote:<br>
&gt; Hi all,<br>
&gt;<br>
&gt; As a starting point with v4, a simple rigid transform image registration (as<br>
&gt; attached at the bottom of this message) was created and tested on a pair of<br>
&gt; 256x256x187 and 256x256x229 images. &nbsp;It takes about 90s.<br>
&gt;<br>
&gt; I had an same registration based on v3, and tested against the same data,<br>
&gt; but it just took about 12s on the same machine with almost the same matching<br>
&gt; result.<br>
&gt;<br>
&gt; As a newbie, I am not sure whether I did the right thing and I am trying to<br>
&gt; understand more about v4. By debugging into the v4 code, I noticed that:<br>
&gt; 1. the new virtual domain (same as the fixed image in my case) introduces<br>
&gt; one more layer which needs some extra computation.<br>
&gt; 2. the transformation on a point was done through two Transformxxx()<br>
&gt; operations by two transform instances in<br>
&gt; ImageToImageMetricv4::m_CompositeTransform, although one of which is<br>
&gt; actually an identity transform.<br>
&gt; and, I am guessing maybe they are reasons for more computing time, but, I am<br>
&gt; not so sure.<br>
&gt;<br>
&gt; Of course, I can just stick to v3, but, I am just curious whether there are<br>
&gt; some ways that I can avoid those extra computations with v4.<br>
&gt;<br>
&gt;<br>
&gt; //=== Start of the code ===================================<br>
&gt; //<br>
&gt; bool<br>
&gt; RigidTransform(itk::CompositeTransform&lt;double,3&gt;::Pointer vComposite,<br>
&gt; &nbsp; &nbsp; ImageType const&amp; vFixImage, ImageType const&amp; vMovImage)<br>
&gt; {<br>
&gt; &nbsp; &nbsp; //- The Euler transform is a rotation and translation about a center, so<br>
&gt; we<br>
&gt; &nbsp; &nbsp; // &nbsp; &nbsp;need to find the rotation center.<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; typedef itk::Euler3DTransform&lt;double&gt; RigidTransformType;<br>
&gt; &nbsp; &nbsp; RigidTransformType::Pointer vRigid = RigidTransformType::New();<br>
&gt; &nbsp; &nbsp; typedef itk::CenteredTransformInitializer&lt; &nbsp; &nbsp;RigidTransformType,<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; ImageType,<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; ImageType &gt;<br>
&gt; InitializerType;<br>
&gt; &nbsp; &nbsp; InitializerType::Pointer Initializer = InitializerType::New();<br>
&gt; &nbsp; &nbsp; Initializer-&gt;SetTransform(vRigid);<br>
&gt; &nbsp; &nbsp; Initializer-&gt;SetFixedImage(&amp;vFixImage);<br>
&gt; &nbsp; &nbsp; Initializer-&gt;SetMovingImage(&amp;vMovImage);<br>
&gt; &nbsp; &nbsp; Initializer-&gt;GeometryOn();<br>
&gt; &nbsp; &nbsp; Initializer-&gt;InitializeTransform();<br>
&gt;<br>
&gt; &nbsp; &nbsp; vComposite-&gt;AddTransform(vRigid);<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- Metric<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; typedef itk::MattesMutualInformationImageToImageMetricv4&lt;ImageType,<br>
&gt; ImageType&gt; MetricType;<br>
&gt; &nbsp; &nbsp; MetricType::Pointer vMetric = MetricType::New();<br>
&gt; &nbsp; &nbsp; vMetric-&gt;SetNumberOfHistogramBins(32);<br>
&gt; &nbsp; &nbsp; vMetric-&gt;SetUseFixedImageGradientFilter(false);<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- Optimizer<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; typedef itk::GradientDescentOptimizerv4 OptimizerType;<br>
&gt; &nbsp; &nbsp; OptimizerType::Pointer vOptimizer = OptimizerType::New();<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetNumberOfIterations( 50 );<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetDoEstimateLearningRateOnce( true );<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetMinimumConvergenceValue( 1e-6 );<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetConvergenceWindowSize( 5 );<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetMaximumStepSizeInPhysicalUnits( 0.5 );<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- Scale estimator<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; itk::OptimizerParameterScalesEstimator::Pointer vScalesEstimator;<br>
&gt; &nbsp; &nbsp; typedef itk::RegistrationParameterScalesFromJacobian&lt;MetricType&gt;<br>
&gt; JacobianScalesEstimatorType;<br>
&gt; &nbsp; &nbsp; {<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; JacobianScalesEstimatorType::Pointer vJacobianScalesEstimator<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; = JacobianScalesEstimatorType::New();<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; vJacobianScalesEstimator-&gt;SetMetric(vMetric);<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; vJacobianScalesEstimator-&gt;SetTransformForward(true);<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; vScalesEstimator = vJacobianScalesEstimator;<br>
&gt; &nbsp; &nbsp; }<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetScalesEstimator(vScalesEstimator);<br>
&gt; &nbsp; &nbsp; vOptimizer-&gt;SetDoEstimateScales(true);<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- The RegistrationMethod class coordinates the registration operation.<br>
&gt; &nbsp; &nbsp; // &nbsp; &nbsp;It needs all the pieces that come together to perform the<br>
&gt; registration<br>
&gt; &nbsp; &nbsp; // &nbsp; &nbsp;operation.<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; typedef itk::ImageRegistrationMethodv4&lt;ImageType, ImageType,<br>
&gt; itk::Euler3DTransform&lt;double&gt;&gt; RigidRegistrationType;<br>
&gt; &nbsp; &nbsp; RigidRegistrationType::Pointer vRigidRegistration =<br>
&gt; RigidRegistrationType::New();<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetOptimizer(vOptimizer);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetFixedImage(&amp;vFixImage);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetMovingImage(&amp;vMovImage);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetMovingInitialTransform(vRigid);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetNumberOfLevels(3);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetMetric(vMetric);<br>
&gt;<br>
&gt; vRigidRegistration-&gt;SetMetricSamplingStrategy(RigidRegistrationType::RANDOM);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetMetricSamplingPercentage(0.1);<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- Shrink the virtual domain by specified factors for each level.<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; RigidRegistrationType::ShrinkFactorsArrayType vRigidShrinkFactors;<br>
&gt; &nbsp; &nbsp; vRigidShrinkFactors.SetSize( 3 );<br>
&gt; &nbsp; &nbsp; vRigidShrinkFactors[0] = 4;<br>
&gt; &nbsp; &nbsp; vRigidShrinkFactors[1] = 2;<br>
&gt; &nbsp; &nbsp; vRigidShrinkFactors[2] = 1;<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetShrinkFactorsPerLevel( vRigidShrinkFactors );<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- Smoothing sigmas array<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; RigidRegistrationType::SmoothingSigmasArrayType vRigidSmoothingSigmas;<br>
&gt; &nbsp; &nbsp; vRigidSmoothingSigmas.SetSize(3);<br>
&gt; &nbsp; &nbsp; vRigidSmoothingSigmas.Fill(0);<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;SetSmoothingSigmasPerLevel(vRigidSmoothingSigmas);<br>
&gt;<br>
&gt; &nbsp; &nbsp; //- Observer<br>
&gt; &nbsp; &nbsp; //<br>
&gt; &nbsp; &nbsp; typedef CommandIterationUpdate&lt; RigidRegistrationType &gt; CommandType;<br>
&gt; &nbsp; &nbsp; CommandType::Pointer observer = CommandType::New();<br>
&gt; &nbsp; &nbsp; vRigidRegistration-&gt;AddObserver( itk::InitializeEvent(), observer );<br>
&gt;<br>
&gt; &nbsp; &nbsp; try<br>
&gt; &nbsp; &nbsp; {<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; std::cout &lt;&lt; "Starting rigid registration..." &lt;&lt; std::endl;<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; vRigidRegistration-&gt;Update();<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; std::cout &lt;&lt; "Rigid parameters after registration: " &lt;&lt; std::endl<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &nbsp; &lt;&lt; vOptimizer-&gt;GetCurrentPosition() &lt;&lt; std::endl;<br>
&gt; &nbsp; &nbsp; }<br>
&gt; &nbsp; &nbsp; catch( itk::ExceptionObject &amp;e )<br>
&gt; &nbsp; &nbsp; {<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; std::cerr &lt;&lt; "Exception caught: " &lt;&lt; e &lt;&lt; std::endl;<br>
&gt; &nbsp; &nbsp; &nbsp; &nbsp; return false;<br>
&gt; &nbsp; &nbsp; }<br>
&gt;<br>
&gt; vComposite-&gt;AddTransform(const_cast&lt;RigidTransformType*&gt;(vRigidRegistration-&gt;GetOutput()-&gt;Get()));<br>
&gt; &nbsp; &nbsp; return true;<br>
&gt; }<br>
&gt; //<br>
&gt; //=== End of the code =====================================<br>
&gt;<br>
&gt;<br>
</div></div>&gt; _____________________________________<br>
&gt; Powered by <a href="http://www.kitware.com/" target="_blank">www.kitware.com</a><br>
&gt;<br>
&gt; Visit other Kitware open-source projects at<br>
&gt; <a href="http://www.kitware.com/opensource/opensource.html" target="_blank">http://www.kitware.com/opensource/opensource.html</a><br>
&gt;<br>
&gt; Kitware offers ITK Training Courses, for more information visit:<br>
&gt; <a href="http://www.kitware.com/products/protraining.php" target="_blank">http://www.kitware.com/products/protraining.php</a><br>
&gt;<br>
&gt; Please keep messages on-topic and check the ITK FAQ at:<br>
&gt; <a href="http://www.itk.org/Wiki/ITK_FAQ" target="_blank">http://www.itk.org/Wiki/ITK_FAQ</a><br>
&gt;<br>
&gt; Follow this link to subscribe/unsubscribe:<br>
&gt; <a href="http://www.itk.org/mailman/listinfo/insight-users" target="_blank">http://www.itk.org/mailman/listinfo/insight-users</a><br>
&gt;<br>
</blockquote></div><br>
_____________________________________<br>Powered by <a href="http://www.kitware.com">www.kitware.com</a><br><br>Visit other Kitware open-source projects at<br><a href="http://www.kitware.com/opensource/opensource.html">http://www.kitware.com/opensource/opensource.html</a><br><br>Kitware offers ITK Training Courses, for more information visit:<br>http://www.kitware.com/products/protraining.php<br><br>Please keep messages on-topic and check the ITK FAQ at:<br>http://www.itk.org/Wiki/ITK_FAQ<br><br>Follow this link to subscribe/unsubscribe:<br>http://www.itk.org/mailman/listinfo/insight-users<br></blockquote></div><br></div></body></html>