[Insight-users] ScaleSkewVersor3DTransform: how to initialize it?

[Margarete Herbst]

hello itkusers,

I would like to initialize a ScaleSkewVersor3DTransform with
center-of-rotation(x,y,z), rotation(x,y,z), translation(x,y,z) and
scaling(x,y,z) (but no skewing) so that an AffineTransform with the same
'external' parameters (rotation, translation, scaling) results in the same
transformation.
I tried this in a test-program (attached below). in the program I try to
initialize an Euler3DTransform, an AffineTransform and a
ScaleSkewVersor3DTransform with the same 'external' parameters and then look
if the resulting transformation is equivalent (it is clear that an
Euler3DTransform will never be equivalent if at least one scaling factor is
different from 1.0). I test the resulting transformation by transforming a
test point and comparing the results.

Results:
- for rigid parameters (scaling = 1.0, 1.0, 1.0) all transformations'
results are equal
- for arbitrary scaling without rotation (0.0, 0.0, 0.0) AffineTransform and
ScaleSkewVersor3DTransform are equal
- BUT: if scaling is arbitrary (e.g. 0.5,1.0,1.0) and rotation is different
from zero (e.g. 0,90deg,0) the results are different!!!

to my knowledge the AffineTransform transforms a point as follows (where
rotation and scaling are w.r.t. center of rotation):
rotate-x, rotate-y, rotate-z, scale-x, scale-y, scale-z, translate-x,
translate-y, translate-z

how should I understand the scaling-vector in ScaleSkewVersor3DTransform?

I would be grateful for comments on my little program.
thank you!

regards, marge


 
#include <iostream>

#include "itkAffineTransform.h"
#include "itkEuler3DTransform.h"
#include "itkScaleSkewVersor3DTransform.h"

using namespace std;

int main(int argc, char *argv[])
{

  typedef double ScalarType;
  typedef itk::Euler3DTransform<ScalarType> EulerTransformType;
  typedef itk::AffineTransform<ScalarType> AffineTransformType;
  typedef itk::ScaleSkewVersor3DTransform<ScalarType>
ScaleSkewTransformType;

  const double degToRad = acos(-1) / 180.;

  double rotx, roty, rotz;
  double translx, transly, translz;
  double centerx, centery, centerz;
  double scalex, scaley, scalez;

  EulerTransformType::InputPointType testPoint;

  // 'external' parameters:
  rotx = 0 * degToRad;
  roty = 90 * degToRad;
  rotz = 0 * degToRad;

  translx = 0;
  transly = 0;
  translz = 0;

  centerx = 0;
  centery = 0;
  centerz = 0;

  scalex = .5;
  scaley = 1.;
  scalez = 1.;

  // 'test point':
  testPoint[0] = 15;
  testPoint[1] = 25;
  testPoint[2] = 35;


  cout << endl << "Test parameters: " << endl;
  cout << " center of rotation: " << centerx << ", " << centery << ", " <<
centerz << endl;
  cout << " rotation: " << rotx << " (" << rotx / degToRad << "°), " << roty
<< " (" << roty / degToRad << "°), " << rotz << " (" << rotz / degToRad <<
"°) " << endl;
  cout << " translation: " << translx << ", " << transly << ", " << translz
<< endl;
  cout << " scaling (affine): " << scalex << ", " << scaley << ", " <<
scalez << endl;
  


  cout << endl << endl << "EULER3D" << endl;

  EulerTransformType::Pointer euler = EulerTransformType::New();
  EulerTransformType::OutputVectorType eulerTransl;
  EulerTransformType::InputPointType eulerCenter;

  euler->SetComputeZYX(true);

  eulerCenter[0] = centerx;
  eulerCenter[1] = centery;
  eulerCenter[2] = centerz;
  euler->SetCenter(eulerCenter);

  euler->SetRotation(rotx, roty, rotz);

  eulerTransl[0] = translx;
  eulerTransl[1] = transly;
  eulerTransl[2] = translz;
  euler->SetTranslation(eulerTransl);

  cout << "matrix: " << endl << euler->GetMatrix() << endl;
  cout << "translation: " << endl << euler->GetTranslation() << endl;
  cout << "parameters: " << endl << euler->GetParameters() << endl;
  cout << "test-point: " << testPoint << " ==> " <<
euler->TransformPoint(testPoint) << endl;



  cout << endl << endl << "AFFINE" << endl;

  AffineTransformType::Pointer affine = AffineTransformType::New();
  AffineTransformType::OutputVectorType affineTransl;
  AffineTransformType::InputPointType affineCenter;
  AffineTransformType::OutputVectorType rotAxis;
  AffineTransformType::OutputVectorType scaleVector;

  affineCenter[0] = centerx;
  affineCenter[1] = centery;
  affineCenter[2] = centerz;
  affine->SetCenter(affineCenter);  

  // scaling outgoing from center!
  scaleVector[0] = scalex;
  scaleVector[1] = scaley;
  scaleVector[2] = scalez;
  affine->Scale(scaleVector, false);

  // z-axis:
  rotAxis[0] = 0;
  rotAxis[1] = 0;
  rotAxis[2] = 1;
  affine->Rotate3D(rotAxis, rotz, true);
  // y-axis:
  rotAxis[0] = 0;
  rotAxis[1] = 1;
  rotAxis[2] = 0;
  affine->Rotate3D(rotAxis, roty, true);
  // x-axis:
  rotAxis[0] = 1;
  rotAxis[1] = 0;
  rotAxis[2] = 0;
  affine->Rotate3D(rotAxis, rotx, true);

  affineTransl[0] = translx;
  affineTransl[1] = transly;
  affineTransl[2] = translz;
  affine->SetTranslation(affineTransl);  

  cout << "matrix: " << endl << affine->GetMatrix() << endl;
  cout << "translation: " << endl << affine->GetTranslation() << endl;  
  cout << "parameters: " << endl << affine->GetParameters() << endl;
  cout << "test-point: " << testPoint << " ==> " <<
affine->TransformPoint(testPoint) << endl;



  cout << endl << endl << "SCALE/SKEW" << endl;

  ScaleSkewTransformType::Pointer scaleSkew = ScaleSkewTransformType::New();
  ScaleSkewTransformType::OutputVectorType scaleSkewTransl;
  ScaleSkewTransformType::InputPointType scaleSkewCenter;
  ScaleSkewTransformType::VersorType scaleSkewRotVersorX;
  ScaleSkewTransformType::VersorType scaleSkewRotVersorY;
  ScaleSkewTransformType::VersorType scaleSkewRotVersorZ;
  ScaleSkewTransformType::ScaleVectorType scaleSkewScaleVector;

  scaleSkewCenter[0] = centerx;
  scaleSkewCenter[1] = centery;
  scaleSkewCenter[2] = centerz;
  scaleSkew->SetCenter(scaleSkewCenter);  

  // rotation:
  scaleSkewRotVersorX.SetIdentity();
  scaleSkewRotVersorX.SetRotationAroundX(rotx);
  scaleSkewRotVersorY.SetIdentity();
  scaleSkewRotVersorY.SetRotationAroundY(roty);
  scaleSkewRotVersorZ.SetIdentity();
  scaleSkewRotVersorZ.SetRotationAroundZ(rotz);
  // -> compose the versors (not commutative!) (rotz -> roty -> rotx):
  scaleSkewRotVersorZ *= scaleSkewRotVersorY;
  scaleSkewRotVersorZ *= scaleSkewRotVersorX;
  scaleSkew->SetRotation(scaleSkewRotVersorZ);

  // scaling outgoing from center!
  scaleSkewScaleVector[0] = scalex;
  scaleSkewScaleVector[1] = scaley;
  scaleSkewScaleVector[2] = scalez;
  scaleSkew->SetScale(scaleSkewScaleVector);

  scaleSkewTransl[0] = translx;
  scaleSkewTransl[1] = transly;
  scaleSkewTransl[2] = translz;
  scaleSkew->SetTranslation(scaleSkewTransl);  

  cout << "matrix: " << endl << scaleSkew->GetMatrix() << endl;
  cout << "translation: " << endl << scaleSkew->GetTranslation() << endl;  
  cout << "parameters: " << endl << scaleSkew->GetParameters() << endl;
  cout << "test-point: " << testPoint << " ==> " <<
scaleSkew->TransformPoint(testPoint) << endl;



  int precision = 10000; // 4 

  cout << endl << endl << "*** TEST-RESULT: EULER3D vs. AFFINE (precision="
<< precision << "): ";
  if (round(euler->TransformPoint(testPoint)[0] * precision) ==
round(affine->TransformPoint(testPoint)[0] * precision) && 
      round(euler->TransformPoint(testPoint)[1] * precision) ==
round(affine->TransformPoint(testPoint)[1] * precision) &&
      round(euler->TransformPoint(testPoint)[2] * precision) ==
round(affine->TransformPoint(testPoint)[2] * precision))
    cout << "SUCCESSFUL";
  else
    cout << "FAILURE";
  cout << endl << endl;

  cout << endl << endl << "*** TEST-RESULT: AFFINE vs. SCALE/SKEW
(precision=" << precision << "): ";
  if (round(scaleSkew->TransformPoint(testPoint)[0] * precision) ==
round(affine->TransformPoint(testPoint)[0] * precision) && 
      round(scaleSkew->TransformPoint(testPoint)[1] * precision) ==
round(affine->TransformPoint(testPoint)[1] * precision) &&
      round(scaleSkew->TransformPoint(testPoint)[2] * precision) ==
round(affine->TransformPoint(testPoint)[2] * precision))
    cout << "SUCCESSFUL";
  else
    cout << "FAILURE";
  cout << endl << endl;


  return EXIT_SUCCESS;
}


-- 
View this message in context: http://www.nabble.com/ScaleSkewVersor3DTransform%3A-how-to-initialize-it--tp20224894p20224894.html
Sent from the ITK - Users mailing list archive at Nabble.com.