/*========================================================================= * * Copyright NumFOCUS * * Licensed under the Apache License, Version 2.0 (the "License"); * you may not use this file except in compliance with the License. * You may obtain a copy of the License at * * https://www.apache.org/licenses/LICENSE-2.0.txt * * Unless required by applicable law or agreed to in writing, software * distributed under the License is distributed on an "AS IS" BASIS, * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. * See the License for the specific language governing permissions and * limitations under the License. * *=========================================================================*/ #ifndef itkCenteredTransformInitializer_hxx #define itkCenteredTransformInitializer_hxx namespace itk { template void CenteredTransformInitializer::InitializeTransform() { // Sanity check if (!m_FixedImage) { itkExceptionMacro("Fixed Image has not been set"); } if (!m_MovingImage) { itkExceptionMacro("Moving Image has not been set"); } if (!m_Transform) { itkExceptionMacro("Transform has not been set"); } // If images come from filters, then update those filters. m_FixedImage->UpdateSource(); m_MovingImage->UpdateSource(); InputPointType rotationCenter; OutputVectorType translationVector; if (m_UseMoments) { m_FixedCalculator->SetImage(m_FixedImage); m_FixedCalculator->Compute(); m_MovingCalculator->SetImage(m_MovingImage); m_MovingCalculator->Compute(); typename FixedImageCalculatorType::VectorType fixedCenter = m_FixedCalculator->GetCenterOfGravity(); typename MovingImageCalculatorType::VectorType movingCenter = m_MovingCalculator->GetCenterOfGravity(); for (unsigned int i = 0; i < InputSpaceDimension; ++i) { rotationCenter[i] = fixedCenter[i]; translationVector[i] = movingCenter[i] - fixedCenter[i]; } } else { // Here use the geometrical center of each image. const typename FixedImageType::RegionType & fixedRegion = m_FixedImage->GetLargestPossibleRegion(); const typename FixedImageType::IndexType & fixedIndex = fixedRegion.GetIndex(); const typename FixedImageType::SizeType & fixedSize = fixedRegion.GetSize(); InputPointType centerFixedPoint; using CoordRepType = typename InputPointType::ValueType; using ContinuousIndexType = ContinuousIndex; using ContinuousIndexValueType = typename ContinuousIndexType::ValueType; ContinuousIndexType centerFixedIndex; for (unsigned int k = 0; k < InputSpaceDimension; ++k) { centerFixedIndex[k] = static_cast(fixedIndex[k]) + static_cast(fixedSize[k] - 1) / 2.0; } m_FixedImage->TransformContinuousIndexToPhysicalPoint(centerFixedIndex, centerFixedPoint); const typename MovingImageType::RegionType & movingRegion = m_MovingImage->GetLargestPossibleRegion(); const typename MovingImageType::IndexType & movingIndex = movingRegion.GetIndex(); const typename MovingImageType::SizeType & movingSize = movingRegion.GetSize(); InputPointType centerMovingPoint; ContinuousIndexType centerMovingIndex; for (unsigned int m = 0; m < InputSpaceDimension; ++m) { centerMovingIndex[m] = static_cast(movingIndex[m]) + static_cast(movingSize[m] - 1) / 2.0; } m_MovingImage->TransformContinuousIndexToPhysicalPoint(centerMovingIndex, centerMovingPoint); for (unsigned int i = 0; i < InputSpaceDimension; ++i) { rotationCenter[i] = centerFixedPoint[i]; translationVector[i] = centerMovingPoint[i] - centerFixedPoint[i]; } } m_Transform->SetCenter(rotationCenter); m_Transform->SetTranslation(translationVector); } template void CenteredTransformInitializer::PrintSelf(std::ostream & os, Indent indent) const { Superclass::PrintSelf(os, indent); itkPrintSelfObjectMacro(Transform); itkPrintSelfObjectMacro(FixedImage); itkPrintSelfObjectMacro(MovingImage); os << indent << "UseMoments: " << (m_UseMoments ? "On" : "Off") << std::endl; itkPrintSelfObjectMacro(MovingCalculator); itkPrintSelfObjectMacro(FixedCalculator); } } // namespace itk #endif /* itkCenteredTransformInitializer_hxx */