[ITK-users] Symmetric eigen analysis discrepancy

[Olivier Commowick]

Hi Matt,

Yes, those are very small numbers indeed :) After digging in the code, I found two discrepancies with respect to the vnl symmetric eigensystem, so I took the vnl code as the reference and did the corrections. I am submitting the two patches on gerrit (I didn’t know they would be separated for review, I can merge them if necessary): 

- http://review.source.kitware.com/#/c/20092/ <http://review.source.kitware.com/#/c/20092/>
- http://review.source.kitware.com/#/c/20093/ <http://review.source.kitware.com/#/c/20093/>

That being said, when digging in the code, it is not clear to me why SymmetricEigenAnalysis is not just a wrapper on vnl_symmetric_eigensystem. This would avoid duplicated code, and I do not see any place in the eigensystem code that wouldn’t be thread safe. But again, I am not an expert of that code and so I may be wrong. In any case, may be something to think of for future work :)

Thanks for the help,
Olivier

> On 13 Aug 2015, at 22:59, Matt McCormick <matt.mccormick at kitware.com> wrote:
> 
> Hi Olivier,
> 
> Those are some pretty small numbers ;-).  Thanks for sharing your experiences.
> 
> It may be necessary to build both ITK and vnl/netlib instances with
> debugging symbols and walk through them in a debugger to see where the
> deviation occurs.
> 
> As always, patches are welcome :-)
> 
> Thanks,
> Matt
> 
> On Thu, Aug 13, 2015 at 9:35 AM, Olivier Commowick
> <Olivier.Commowick at inria.fr> wrote:
>> Hello all,
>> 
>> I have been experiencing a problem when using the itk SymmetricEigenAnalysis class in ITK. When using this class, it sometimes produces NaN eigenvalues/vectors. On the contrary, the vnl_symmetric_eigensystem doesn’t seem to have a problem (and neither does matlab). Looking at the documentation, it seems both use QL decomposition but not implemented in the same place. In addition, the vnl one is said to be not thread safe, although it’s the one working for me.
>> 
>> In fact, I have isolated one 3x3 matrix for which I always have the problem.
>> 
>> A =
>> [7.503107128706855667859804009212206062302e-53 1.870107888834986371214720644344813262375e-52 -8.662047753682068190087788895125120421821e-27;
>> 1.870107888834986371214720644344813262375e-52 4.661140319458029838073740773908360767728e-52 -2.15896741973055935210319454393408888029e-26;
>> -8.662047753682068190087788895125120421821e-27 -2.15896741973055935210319454393408888029e-26 1]
>> 
>> It’s a weird one I know, but again vnl (and matlab) works great on it. Would there be anyone who knows the QL decomposition code enough to help me understand and possibly correct the ITK class so that it does the same as the VNL one ?
>> 
>> Thanks a lot in advance.
>> 
>> ---
>> 
>> Olivier Commowick, Ph.D.
>> Research Scientist
>> INRIA Rennes - Bretagne Atlantique, VISAGES Team
>> Campus de Beaulieu
>> 35042 Rennes
>> FRANCE
>> 
>> Phone: +33 2 99 84 25 92
>> Email: Olivier.Commowick at inria.fr
>> Web: http://olivier.commowick.org/
>> 
>> 
>> 
>> 
>> 
>> _____________________________________
>> Powered by www.kitware.com
>> 
>> Visit other Kitware open-source projects at
>> http://www.kitware.com/opensource/opensource.html
>> 
>> Kitware offers ITK Training Courses, for more information visit:
>> http://www.kitware.com/products/protraining.php
>> 
>> Please keep messages on-topic and check the ITK FAQ at:
>> http://www.itk.org/Wiki/ITK_FAQ
>> 
>> Follow this link to subscribe/unsubscribe:
>> http://public.kitware.com/mailman/listinfo/insight-users


---

Olivier Commowick, Ph.D.
Research Scientist
INRIA Rennes - Bretagne Atlantique, VISAGES Team
Campus de Beaulieu
35042 Rennes
FRANCE

Phone: +33 2 99 84 25 92
Email: Olivier.Commowick at inria.fr
Web: http://olivier.commowick.org/





-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://public.kitware.com/pipermail/insight-users/attachments/20150817/364a7a1b/attachment.html>