BUG: Add ComputeZYX to FixedParameters of AdvancedEuler3DTransform

Following ITK commit InsightSoftwareConsortium/ITK@47c2491 "BUG: Interpretation of the Euler angles ZYX or ZXY was not exported.", Ziv Yaniv, May 9, 2016.

Implementation of `GetFixedParameters`, `SetFixedParameters`, and `SetComputeZYX` mostly copied from "itkEuler3DTransform.h" and "itkEuler3DTransform.hxx", in  https://github.com/InsightSoftwareConsortium/ITK/tree/v5.4.0/Modules/Core/Transform/include

GoogleTest `TransformedPointSameAsITKEuler3D` extended to test ComputeZYX = true. `DefaultElastixFixedParametersAreZeroOrEmpty` test adjusted to allow fixedParameters.size() unequal to NDimension. `CopyDefaultEulerTransformElastix3DFixedParametersToCorrespondingItkTransform` test removed, as 3D EulerTransform is now fully covered by the `CopyDefaultParametersToCorrespondingItkTransform` test.

Common/GTesting/elxTransformIOGTest.cxx

@@ -142,8 +142,11 @@ struct WithDimension
     Expect_default_elastix_FixedParameters_are_all_zero()
     {
       const auto fixedParameters = CheckNew<ElastixTransformType>()->GetFixedParameters();
-      ASSERT_EQ(fixedParameters.size(), NDimension);
-      ASSERT_EQ(fixedParameters, vnl_vector<double>(NDimension, 0.0));
+
+      for (const auto fixedParameter : fixedParameters)
+      {
+        EXPECT_EQ(fixedParameter, 0.0);
+      }
     }
 
 
@@ -777,32 +780,9 @@ GTEST_TEST(TransformIO, CopyDefaultParametersToCorrespondingItkTransform)
 
     WithDimension<2>::WithElastixTransform<elx::SimilarityTransformElastix>::Test_copying_default_parameters<
       itk::Euler2DTransform<double>>(fixed);
+    WithDimension<3>::WithElastixTransform<elx::EulerTransformElastix>::Test_copying_default_parameters<
+      itk::Euler3DTransform<double>>(fixed);
   }
-  WithDimension<3>::WithElastixTransform<elx::EulerTransformElastix>::Test_copying_default_parameters<
-    itk::Euler3DTransform<double>>(false);
-  // See also CopyDefaultEulerTransformElastix3DFixedParametersToCorrespondingItkTransform
-}
-
-
-GTEST_TEST(TransformIO, CopyDefaultEulerTransformElastix3DFixedParametersToCorrespondingItkTransform)
-{
-  const auto elxTransform = CheckNew<elx::EulerTransformElastix<ElastixType<3>>>();
-  const auto compositeTransform = elx::TransformIO::ConvertToItkCompositeTransform(*elxTransform);
-  ASSERT_NE(compositeTransform, nullptr);
-  const auto & transformQueue = compositeTransform->GetTransformQueue();
-  ASSERT_EQ(transformQueue.size(), 1);
-  const auto & itkTransform = transformQueue.front();
-  ASSERT_NE(itkTransform, nullptr);
-
-  const auto elxFixedParameters = elxTransform->GetFixedParameters();
-  itkTransform->SetFixedParameters(elxFixedParameters);
-  const auto itkFixedParameters = itkTransform->GetFixedParameters();
-
-  // Note: ideally itkFixedParameters and elxFixedParameters should be equal!
-  EXPECT_NE(itkFixedParameters, elxFixedParameters);
-
-  ASSERT_GE(itkFixedParameters.size(), elxFixedParameters.size());
-  ASSERT_EQ(itkFixedParameters, vnl_vector<double>(4, 0));
 }
 
 
@@ -918,12 +898,16 @@ GTEST_TEST(Transform, TransformedPointSameAsITKEuler2D)
 
 GTEST_TEST(Transform, TransformedPointSameAsITKEuler3D)
 {
-  elx::DefaultConstruct<itk::Euler3DTransform<double>> itkTransform;
-  itkTransform.SetTranslation(itk::MakeVector(1.0, 2.0, 3.0));
-  itkTransform.SetCenter(itk::MakePoint(3.0, 2.0, 1.0));
-  itkTransform.SetRotation(M_PI_2, M_PI_4, M_PI_4 / 2.0);
+  for (const bool computeZYX : { false, true })
+  {
+    elx::DefaultConstruct<itk::Euler3DTransform<double>> itkTransform;
+    itkTransform.SetTranslation(itk::MakeVector(1.0, 2.0, 3.0));
+    itkTransform.SetCenter(itk::MakePoint(3.0, 2.0, 1.0));
+    itkTransform.SetRotation(M_PI_2, M_PI_4, M_PI_4 / 2.0);
+    itkTransform.SetComputeZYX(computeZYX);
 
-  Expect_elx_TransformPoint_yields_same_point_as_ITK<elx::EulerTransformElastix>(itkTransform);
+    Expect_elx_TransformPoint_yields_same_point_as_ITK<elx::EulerTransformElastix>(itkTransform);
+  }
 }
 
 

Common/Transforms/itkAdvancedEuler3DTransform.h

@@ -85,6 +85,7 @@ class ITK_TEMPLATE_EXPORT AdvancedEuler3DTransform : public AdvancedRigid3DTrans
   itkStaticConstMacro(ParametersDimension, unsigned int, 6);
 
   using typename Superclass::ParametersType;
+  using typename Superclass::FixedParametersType;
   using typename Superclass::NumberOfParametersType;
   using typename Superclass::JacobianType;
   using typename Superclass::ScalarType;
@@ -120,6 +121,12 @@ class ITK_TEMPLATE_EXPORT AdvancedEuler3DTransform : public AdvancedRigid3DTrans
   const ParametersType &
   GetParameters() const override;
 
+  const FixedParametersType &
+  GetFixedParameters() const override;
+
+  void
+  SetFixedParameters(const FixedParametersType & parameters) override;
+
   /** Set the rotational part of the transform. */
   void
   SetRotation(ScalarType angleX, ScalarType angleY, ScalarType angleZ);
@@ -132,8 +139,15 @@ class ITK_TEMPLATE_EXPORT AdvancedEuler3DTransform : public AdvancedRigid3DTrans
   void
   GetJacobian(const InputPointType &, JacobianType &, NonZeroJacobianIndicesType &) const override;
 
-  /** Set/Get the order of the computation. Default ZXY */
-  itkSetMacro(ComputeZYX, bool);
+  /** The Euler angle representation of a rotation is not unique and
+   * depends on the order of rotations. In general there are 12
+   * options. This class supports two of them, ZXY and ZYX. The
+   * default is ZXY. These functions set and get the value which
+   * indicates whether the rotation is ZYX or ZXY.
+   */
+  virtual void
+  SetComputeZYX(const bool flag);
+
   itkGetConstMacro(ComputeZYX, bool);
 
   void

Common/Transforms/itkAdvancedEuler3DTransform.hxx

@@ -45,6 +45,8 @@ AdvancedEuler3DTransform<TScalarType>::AdvancedEuler3DTransform()
 {
   m_ComputeZYX = false;
   m_AngleX = m_AngleY = m_AngleZ = ScalarType{};
+  Superclass::m_FixedParameters.SetSize(SpaceDimension + 1);
+  Superclass::m_FixedParameters.Fill(0.0);
   this->PrecomputeJacobianOfSpatialJacobian();
 }
 
@@ -93,6 +95,59 @@ AdvancedEuler3DTransform<TScalarType>::GetParameters() const -> const Parameters
   return this->m_Parameters;
 }
 
+
+template <class TScalarType>
+auto
+AdvancedEuler3DTransform<TScalarType>::GetFixedParameters() const -> const FixedParametersType &
+{
+  if (const auto numberOfFixedParameters = Superclass::m_FixedParameters.size();
+      numberOfFixedParameters <= SpaceDimension)
+  {
+    itkExceptionMacro("Error getting fixed parameters: number of fixed parameters (" << numberOfFixedParameters
+                                                                                     << ") is less than expected");
+  }
+
+  const auto & center = this->GetCenter();
+
+  for (unsigned int i = 0; i < SpaceDimension; ++i)
+  {
+    Superclass::m_FixedParameters[i] = center[i];
+  }
+
+  Superclass::m_FixedParameters[3] = m_ComputeZYX ? 1.0 : 0.0;
+  return Superclass::m_FixedParameters;
+}
+
+
+template <class TScalarType>
+void
+AdvancedEuler3DTransform<TScalarType>::SetFixedParameters(const FixedParametersType & parameters)
+{
+  if (parameters.size() < InputSpaceDimension)
+  {
+    itkExceptionMacro("Error setting fixed parameters: parameters array size ("
+                      << parameters.size() << ") is less than expected  (InputSpaceDimension = " << InputSpaceDimension
+                      << ')');
+  }
+
+  InputPointType c;
+  for (unsigned int i = 0; i < InputSpaceDimension; ++i)
+  {
+    c[i] = Superclass::m_FixedParameters[i] = parameters[i];
+  }
+  this->SetCenter(c);
+  // conditional is here for backwards compatibility: the
+  // m_ComputeZYX flag was not serialized so it may or may
+  // not be included as part of the fixed parameters
+  if (parameters.Size() == 4)
+  {
+    const auto parameter = parameters[3];
+    Superclass::m_FixedParameters[3] = parameter;
+    this->SetComputeZYX(parameter != 0.0);
+  }
+}
+
+
 // Set Rotational Part
 template <class TScalarType>
 void
@@ -385,6 +440,22 @@ AdvancedEuler3DTransform<TScalarType>::PrecomputeJacobianOfSpatialJacobian()
 }
 
 
+template <class TScalarType>
+void
+AdvancedEuler3DTransform<TScalarType>::SetComputeZYX(const bool flag)
+{
+  if (m_ComputeZYX != flag)
+  {
+    m_ComputeZYX = flag;
+    this->ComputeMatrix();
+    this->ComputeOffset();
+    // The meaning of the parameters has changed so the transform
+    // has been modified even if the parameter values have not.
+    this->Modified();
+  }
+}
+
+
 // Print self
 template <class TScalarType>
 void