[Mrtrix-discussion] Evaluate SH in MRview
Peter Neher
p.neher at dkfz-heidelberg.de
Mon Jan 14 06:46:25 PST 2013
Hi Donald,
thanks for your help. We will try to figure out the differences between
our and your code. If we find the issue I will let you know.
Best,
Peter
On 01/07/2013 02:44 AM, Donald Tournier wrote:
> Hi Peter,
>
> OK, at least this confirms that MRtrix is internally consistent. As to
> what the discrepancy might be due to, I'm not sure. The "csdeconv
> -help" command should already provide a good idea of how the
> coefficients are stored. The only other thing I can do is point out
> the exact code used in MRtrix: indexing of the SH coefficients is done
> as per line 41 in src/dwi/SH.h:
>
> 41 inline int index (int l, int m) { return (l*(l+1)/2 + m); }
>
> This basically means pack all /m/ terms of each even harmonic degree
> /l/ sequentially, in order of increasing /l/. As to the definition of
> the basis functions, the simplest form of it is at line 148 in
> src/dwi/SH.cpp:
>
> 148 float value (float azimuth, float elevation, int l, int m)
> 149 {
> 150 elevation = gsl_sf_legendre_sphPlm (l, abs(m),
> cos(elevation));
> 151 if (!m) return (elevation);
> 152 if (m > 0) return (elevation * cos (m*azimuth));
> 153 return (elevation * sin (-m*azimuth));
> 154 }
>
> See the GSL documentation
> <http://www.gnu.org/software/gsl/manual/html_node/Associated-Legendre-Polynomials-and-Spherical-Harmonics.html>
> for the definition of gsl_sf_legendre_sphPlm(). You'd have to see
> exactly what convention the other packages use to figure out where the
> differences are...
>
> Ariel: I hope this is the information you were after. Otherwise, I'm
> not sure what you mean by "the actual numbers that comprise the basis
> set"...
>
> Hope this helps.
> Cheers,
>
> Donald.
>
>
> On 4 January 2013 22:54, Peter Neher <p.neher at dkfz-heidelberg.de
> <mailto:p.neher at dkfz-heidelberg.de>> wrote:
>
> H Donald,
>
> I checked the CSD result with "disp_profile" and the result looks
> just fine (attachment). Only if I try to display the SH
> coefficients in MITK I get this strange rotation-like error. If I
> try it the other way round, meaning if I try to visualize the
> coefficients reconstructed with MITK and also FSL using
> "disp_profile" I also get such a rotation, so there seems to be
> some discrepancy in the handling of those coefficients. MITK and
> FSL don't perform a CSD but a CSA-Q-ball reconstruction, but that
> still the results should look somehow comparable. Unfortunately I
> have no idea how the SH coefficients could be handled differently
> to explain these discrepancies.
>
> Best
> Peter
>
>
> On 01/03/2013 11:55 PM, Donald Tournier wrote:
>> Hi Peter,
>>
>> As René mentioned, the SH basis in MRtrix is not orthonormal -
>> something that I hope will be fixed in a future major release.
>> That said, I'm not sure this is the problem you're having, it
>> wouldn't cause a pure rotation.
>>
>> I'm struggling to understand the figure you sent through. The ODF
>> doesn't look like it was displayed using MRtrix - what software
>> was used for that, and what SH basis did it assume? What does the
>> ODF look like when you display it using MRtrix's disp_profile
>> command? Also, which directions are the peaks actually pointing
>> along? The effect of the non-orthonormality of the basis will be
>> different for peaks aligned along z than for peaks in the x-y
>> plane, causing a differential scaling in this case. Also, you
>> said the directions estimated by find_SH_peaks are correct, but
>> they clearly don't match the ODF lobes as displayed. If so, and
>> the ODF was not displayed using MRtrix, then I would suspect that
>> whatever routine you use for display is at fault. Using
>> disp_profile would at least verify that MRtrix is consistent with
>> itself (i.e. the ODF lobes it displays should point along the
>> directions estimated by find_SH_peaks).
>>
>> Finally, If you're going to change the code in MRtrix, you'll
>> need to know that there are a number of functions that compute SH
>> functions, each coded up slightly differently (for performance
>> reasons) - that said, they are all in the same file, so it's not
>> difficult to modify them all.
>>
>> Hope that helps,
>>
>> Donald.
>>
>>
>> On 4 January 2013 01:22, Peter Neher <p.neher at dkfz-heidelberg.de
>> <mailto:p.neher at dkfz-heidelberg.de>> wrote:
>>
>> Hi René,
>>
>> thank you for your reply! Unfortunately that did not fix my
>> problem. I attached an image of the ODF I get with my method
>> (negative Lobes are dark blue). The red lines indicate the
>> peaks extracted from the SH coefficient file using
>> "find_SH_peaks". These peaks are correct but the ODF is
>> somehow rotated. The ODF should actually represent a simple
>> 90° crossing of two fibers corresponding to the two peaks.
>>
>> Here a code snippet of my SH basis calculation. Maybe someone
>> can directly identify the issue:
>>
>> j=0;
>> for (int l=0; l<=ShOrder; l=l+2)
>> for (m=-l; m<=l; m++)
>> {
>> // evaluation of legendre polynome
>> float mag =
>> sqrt((double)(2*l+1)/(4.0*M_PI)*factorial<double>(l-abs(m))/factorial<double>(l+abs(m)))
>> * legendre_p<double>(l,abs(m),cos(sphCoords(0)));
>>
>> if (m<0)
>> m_ShBasis(j) = sqrt(2.0)*mag*cos(-m*sphCoords(1));
>> else if (m==0)
>> m_ShBasis(j) = mag;
>> else
>> m_ShBasis(j) = pow(-1.0,
>> m)*sqrt(2.0)*mag*sin(m*sphCoords(1));
>>
>> j++;
>> }
>>
>> I also tried it without the sqrt(2.0) and pow(-1.0, m) terms,
>> but that does not change anything.
>>
>> Best,
>> Peter
>>
>> On 01/03/2013 01:01 PM, René Besseling wrote:
>>> Hi Peter,
>>> I ran into this before; it has to do with the MRtrix SH
>>> basis being orthogonal but not orthonormal, see the
>>> following quote from an e-mail from Donald.
>>> Best regards and happy new year to you too,
>>>
>>> René
>>>
>>> Quote Donald:
>>> "To answer your question: the matlab code does use an
>>> orthonormal basis, but MRtrix unfortunately does not. Main
>>> reason is that I didn't bother checking that the basis was
>>> orthonormal in the MRtrix code until relatively recently,
>>> and it's now too late to change it without making all
>>> current data sets ambiguous (i.e. how do you know whether a
>>> data set uses the orthonormal basis or not?). I have since
>>> changed the matlab code though, which is why there is a
>>> difference.
>>>
>>> It's pretty easy to convert between them: the m=0 terms are
>>> not affected, while the m!=0 terms are all scaled by a
>>> factor of sqrt(2) in MRtrix compared with Matlab. Easy to
>>> fix, but very annoying I have to admit."
>>>
>>>
>>> On Thu, Jan 3, 2013 at 12:54 PM, Peter Neher
>>> <p.neher at dkfz-heidelberg.de
>>> <mailto:p.neher at dkfz-heidelberg.de>> wrote:
>>>
>>> Hi everyone and happy new year!
>>>
>>> I am trying to import the SH coefficient file (output of
>>> csdeconv) into my own program but the ODFs are not
>>> displayed correctly. The same method worked fine for the
>>> FSL coefficient files so I guess you are calculating the
>>> SH basis in a different way (the file format seems to be
>>> the same as the FSL file format). MRview renders the
>>> ODFs correctly, so I wanted to compare may code to the
>>> MRview code. Can you tell me the location in the source
>>> code where I can find your calculation and evaluation of
>>> the SH basis?
>>>
>>> Best,
>>> Peter
>>> _______________________________________________
>>> Mrtrix-discussion mailing list
>>> Mrtrix-discussion at www.nitrc.org
>>> <mailto:Mrtrix-discussion at www.nitrc.org>
>>> http://www.nitrc.org/mailman/listinfo/mrtrix-discussion
>>>
>>>
>>
>> --
>> Dipl.-Inform. Peter Neher
>> German Cancer Research Center
>> (DKFZ, Deutsches Krebsforschungszentrum in der Helmholtz-Gemeinschaft, Stiftung des öffentlichen Rechts)
>> Division of Medical and Biological Informatics
>> Im Neuenheimer Feld 280, D-69120 Heidelberg
>>
>> Phone: 49-(0)6221-42-3552, Fax: 49-(0)6221-42-2345
>> E-Mail:p.neher at dkfz-heidelberg.de <mailto:p.neher at dkfz-heidelberg.de>, Web:www.dkfz.de <http://www.dkfz.de>
>>
>> The information contained in this message may be confidential and legally protected under applicable law. The message is intended solely for the addressee(s). If you are not the intended recipient, you are hereby notified that any use, forwarding, dissemination, or reproduction of this message is strictly prohibited and may be unlawful. If you are not the intended recipient, please contact the sender by return e-mail and destroy all copies of the original message.
>>
>>
>> _______________________________________________
>> Mrtrix-discussion mailing list
>> Mrtrix-discussion at www.nitrc.org
>> <mailto:Mrtrix-discussion at www.nitrc.org>
>> http://www.nitrc.org/mailman/listinfo/mrtrix-discussion
>>
>>
>>
>>
>> --
>> *Dr Jacques-Donald Tournier
>> *
>> Research Fellow
>>
>> The Florey Institute of Neuroscience and Mental Health
>> Melbourne Brain Centre - Austin Campus
>> 245 Burgundy Street
>> Heidelberg Vic 3084
>> Ph: +61 3 9035 7033
>> Fax: +61 3 9035 7307
>> www.florey.edu.au <http://www.florey.edu.au>
>>
>
> --
> Dipl.-Inform. Peter Neher
> German Cancer Research Center
> (DKFZ, Deutsches Krebsforschungszentrum in der Helmholtz-Gemeinschaft, Stiftung des öffentlichen Rechts)
> Division of Medical and Biological Informatics
> Im Neuenheimer Feld 280, D-69120 Heidelberg
>
> Phone: 49-(0)6221-42-3552, Fax: 49-(0)6221-42-2345
> E-Mail:p.neher at dkfz-heidelberg.de <mailto:p.neher at dkfz-heidelberg.de>, Web:www.dkfz.de <http://www.dkfz.de>
>
> The information contained in this message may be confidential and legally protected under applicable law. The message is intended solely for the addressee(s). If you are not the intended recipient, you are hereby notified that any use, forwarding, dissemination, or reproduction of this message is strictly prohibited and may be unlawful. If you are not the intended recipient, please contact the sender by return e-mail and destroy all copies of the original message.
>
>
>
>
> --
> *Dr Jacques-Donald Tournier
> *
> Research Fellow
>
> The Florey Institute of Neuroscience and Mental Health
> Melbourne Brain Centre - Austin Campus
> 245 Burgundy Street
> Heidelberg Vic 3084
> Ph: +61 3 9035 7033
> Fax: +61 3 9035 7307
> www.florey.edu.au <http://www.florey.edu.au>
>
--
Dipl.-Inform. Peter Neher
German Cancer Research Center
(DKFZ, Deutsches Krebsforschungszentrum in der Helmholtz-Gemeinschaft, Stiftung des öffentlichen Rechts)
Division of Medical and Biological Informatics
Im Neuenheimer Feld 280, D-69120 Heidelberg
Phone: 49-(0)6221-42-3552, Fax: 49-(0)6221-42-2345
E-Mail: p.neher at dkfz-heidelberg.de, Web: www.dkfz.de
The information contained in this message may be confidential and legally protected under applicable law. The message is intended solely for the addressee(s). If you are not the intended recipient, you are hereby notified that any use, forwarding, dissemination, or reproduction of this message is strictly prohibited and may be unlawful. If you are not the intended recipient, please contact the sender by return e-mail and destroy all copies of the original message.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://www.nitrc.org/pipermail/mrtrix-discussion/attachments/20130114/5ab7e666/attachment-0001.html
More information about the Mrtrix-discussion
mailing list