Hi Michiel,
This is actually the reason why using Euler angles is a bad idea.
The lack of a standard for the order of
application of the individual rotations.
plus their vulnerability to the gimbal lock problem:
http://en.wikipedia.org/wiki/Gimbal_lock
http://en.wikipedia.org/wiki/Gimbal_lock#Loss_of_a_degree_of_freedom_with_Euler_angles
---
I would strongly encourage you to use itk::Versors instead.
Versors (Unit Quaternions) are an unambiguous representation of
rotations.
Regards,
Luis
------------------------------------------------------------------------
On Tue, Mar 30, 2010 at 9:11 AM, michiel mentink
<michael.mentink at st-hughs.ox.ac.uk> wrote:
> in
> /itkEuler3DTransform.txx
> the rotation is applied either in X Y Z order, or Z X Y order (standard),
> where Y rotation is applied first.
>
> Indeed, if I try different orders of rotation, I get different outcomes. I
> suppose this is due to round-off
> errors.
> Why is Z X Y chosen as the standard way of multiplication?
>
>
> Z Y X
> image->GetDirection():
> 0.999695 0.0171452 0.0177517
> -0.0174497 0.999701 0.0171452
> -0.0174524 -0.0174497 0.999695
>
> X Y Z
> image->GetDirection():
> 0.999695 0.0174497 0.0174524
> -0.0177543 0.99969 0.0174497
> -0.0171425 -0.0177543 0.999695
>
> Z X Y
> image->GetDirection():
> 0.99969 0.0174497 0.0177543
> -0.0177543 0.999695 0.0171425
> -0.0174497 -0.0174524 0.999695
>
> cheers,
>
> Michael
>
>
>
>
>
> On Tue, Mar 30, 2010 at 12:18 PM, michiel mentink
> <michael.mentink at st-hughs.ox.ac.uk> wrote:
>>
>> ah just after posting I found it (figures!)
>>
>> complete and correct code now:
>>
>> const ImageType::DirectionType & direction = image->GetDirection()
>> std::cout << "direction: " << std::endl << direction << std::endl <<
>> std::endl;
>>
>> float angleX, angleY, angleZ;
>> angleX = angleY = angleZ = 1 * vnl_math::pi / 180.0; //
>> in this case: rotation = 1 degree
>>
>> const double cx = vcl_cos(angleX);
>> const double sx = vcl_sin(angleX);
>>
>> typedef itk::Matrix<double,3,3> Matrix;
>> Matrix RotationX;
>> Matrix FinalRotation = image->GetDirection();
>>
>> std::cout << "sin: " << sx << " Cos: " << cx << std::endl;
>>
>> RotationX[0][0] = 1; RotationX[0][1] = 0; RotationX[0][2] = 0;
>> RotationX[1][0] = 0; RotationX[1][1] = cx; RotationX[1][2] = sx;
>> RotationX[2][0] = 0; RotationX[2][1] = -sx; RotationX[2][2] = cx;
>>
>> FinalRotation = direction*RotationX;
>>
>> std::cout << "image->GetDirection(): " << std::endl <<
>> image->GetDirection() << std::endl;
>> std::cout << "RotationX: " << std::endl << RotationX << std::endl;
>> std::cout << "FinalRotation: " << std::endl << FinalRotation <<
>> std::endl;
>>
>> image->SetDirection(FinalRotation);
>>
>> std::cout << "image->GetDirection(): " << std::endl <<
>> image->GetDirection() << std::endl;
>>
>>
>> cheers, Michael
>>
>>
>>
>>
>>
>> On Tue, Mar 30, 2010 at 12:10 PM, michiel mentink
>> <michael.mentink at st-hughs.ox.ac.uk> wrote:
>>>
>>> Hello Frederic,
>>>
>>> thanks for your suggestion.
>>>
>>> although my code was working, I tried your suggestion. Unfortunately, it
>>> produces
>>>
>>> error: no match for ‘operator*’ in
>>> ‘image.itk::SmartPointer<TObjectType>::operator-> [with TObjectType =
>>> itk::Image<float, 3u>]()->itk::ImageBase<VImageDimension>::GetDirection
>>> [with unsigned int VImageDimension = 3u] * RotationX’
>>>
>>>
>>> Anyway, I forgot to mention to use:
>>>
>>> float angleX, angleY, angleZ;
>>> angleX = angleY = angleZ = 1 * (3.14/180);
>>>
>>> (multiply by pi and divide by 180 degrees, because ITK internally works
>>> with radians instead of degrees)
>>>
>>> Which leads me to the question: vnl has a pi constant, and I remember
>>> vaguely having seen it somewhere as vnl::PI or something.
>>> Does anyone know how to convince vnl to hand me pi constant?
>>>
>>> cheers, Michael
>>>
>>>
>>>
>>> On Tue, Mar 30, 2010 at 11:55 AM, Frederic Perez <fredericpcx at gmail.com>
>>> wrote:
>>>>
>>>> Hello Michiel,
>>>>
>>>> perhaps you could use a const Matrix object after all, since it looks to
>>>> me that FinalRotation is first built with image->GetDirection() but this
>>>> value is not actually used, and that the signature of itk::Image's is
>>>> SetDirection(const DirectionType direction).
>>>>
>>>> So here you are, my quickly written proposal (caution, I haven't
>>>> compiled it):
>>>>
>>>> float angleX, angleY, angleZ;
>>>> angleX = angleY = angleZ = 5;
>>>>
>>>> const double cx = vcl_cos(angleX);
>>>> const double sx = vcl_sin(angleX);
>>>>
>>>> typedef itk::Matrix<double,3,3> Matrix;
>>>> Matrix RotationX;
>>>> // Matrix FinalRotation = image->GetDirection(); -- Commented now
>>>>
>>>> RotationX[0][0] = 1; RotationX[0][1] = 0; RotationX[0][2] = 0;
>>>> RotationX[1][0] = 0; RotationX[1][1] = cx; RotationX[1][2] = sx;
>>>> RotationX[2][0] = 0; RotationX[2][1] = -sx; RotationX[2][2] = cx;
>>>>
>>>> const Matrix FinalRotation = direction*RotationX;
>>>>
>>>> std::cout << "image->GetDirection(): " << std::endl <<
>>>> image->GetDirection() << std::endl;
>>>> std::cout << "RotationX: " << std::endl << RotationX << std::endl;
>>>> std::cout << "FinalRotation: " << std::endl << FinalRotation <<
>>>> std::endl;
>>>>
>>>> image->SetDirection(FinalRotation);
>>>>
>>>> std::cout << "image->GetDirection(): " << std::endl <<
>>>> image->GetDirection() << std::endl;
>>>>
>>>> Cheers,
>>>>
>>>> Frederic
>>>>
>>>
>>
>
>