itkCurvilinearArraySpecialCoordinatesImage.h

Go to the documentation of this file.

/*============================================================================

The Medical Imaging Interaction Toolkit (MITK)

Copyright (c) German Cancer Research Center (DKFZ)
All rights reserved.

Use of this source code is governed by a 3-clause BSD license that can be
found in the LICENSE file.

============================================================================*/
/*=========================================================================
 *
 *  Copyright Insight Software Consortium
 *
 *  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
 *
 *         http://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 itkCurvilinearArraySpecialCoordinatesImage_h
#define itkCurvilinearArraySpecialCoordinatesImage_h

#include "itkSpecialCoordinatesImage.h"
#include "itkPoint.h"
#include "vnl/vnl_math.h"
#include "itkNeighborhoodAccessorFunctor.h"

namespace itk
{
template< typename TPixel, unsigned int VDimension >
class CurvilinearArraySpecialCoordinatesImage:
  public SpecialCoordinatesImage< TPixel, VDimension >
{
public:
  typedef CurvilinearArraySpecialCoordinatesImage       Self;
  typedef SpecialCoordinatesImage< TPixel, VDimension > Superclass;
  typedef SmartPointer< Self >                          Pointer;
  typedef SmartPointer< const Self >                    ConstPointer;
  typedef WeakPointer< const Self >                     ConstWeakPointer;

  itkNewMacro(Self);

  itkTypeMacro(CurvilinearArraySpecialCoordinatesImage, SpecialCoordinatesImage);

  typedef TPixel PixelType;

  typedef TPixel ValueType;

  typedef TPixel InternalPixelType;

  typedef typename Superclass::IOPixelType IOPixelType;

  typedef DefaultPixelAccessor< PixelType > AccessorType;

  typedef DefaultPixelAccessorFunctor< Self > AccessorFunctorType;

  typedef NeighborhoodAccessorFunctor< Self > NeighborhoodAccessorFunctorType;

  itkStaticConstMacro(ImageDimension, unsigned int, VDimension);

  typedef typename Superclass::IndexType      IndexType;
  typedef typename Superclass::IndexValueType IndexValueType;

  typedef typename Superclass::OffsetType    OffsetType;

  typedef typename Superclass::SizeType      SizeType;
  typedef typename Superclass::SizeValueType SizeValueType;

  typedef ImportImageContainer< SizeValueType, PixelType > PixelContainer;

  typedef typename Superclass::RegionType RegionType;

  typedef typename Superclass::SpacingType SpacingType;

  typedef typename Superclass::PointType PointType;

  typedef typename PixelContainer::Pointer      PixelContainerPointer;
  typedef typename PixelContainer::ConstPointer PixelContainerConstPointer;

  virtual void Graft(const DataObject *data) ITK_OVERRIDE;

  template< typename TCoordRep >
  bool TransformPhysicalPointToContinuousIndex(
    const Point< TCoordRep, VDimension > & point,
    ContinuousIndex< TCoordRep, VDimension > & index) const
  {
    const RegionType & region = this->GetLargestPossibleRegion();
    const double maxLateral = region.GetSize(1) - 1;

    // Convert Cartesian coordinates into angular coordinates
    const TCoordRep lateral = std::atan(point[0] / point[1]);
    const TCoordRep radius  = std::sqrt(point[0] * point[0] + point[1] * point[1] );

    // Convert the "proper" angular coordinates into index format
    index[0] = static_cast< TCoordRep >( ( ( radius - m_FirstSampleDistance )
                                           / m_RadiusSampleSize ) );
    index[1] = static_cast< TCoordRep >( ( lateral / m_LateralAngularSeparation )
                                         + ( maxLateral / 2.0 ) );

    // Now, check to see if the index is within allowed bounds
    const bool isInside = region.IsInside(index);

    return isInside;
  }

  template< typename TCoordRep >
  bool TransformPhysicalPointToIndex(
    const Point< TCoordRep, VDimension > & point,
    IndexType & index) const
  {
    const RegionType & region = this->GetLargestPossibleRegion();
    const double maxLateral = region.GetSize(1) - 1;

    // Convert Cartesian coordinates into angular coordinates
    const TCoordRep lateral = std::atan(point[0] / point[1]);
    const TCoordRep radius  = std::sqrt(point[0] * point[0] + point[1] * point[1] );

    // Convert the "proper" angular coordinates into index format
    index[0] = static_cast< IndexValueType >( ( ( radius - m_FirstSampleDistance )
                                           / m_RadiusSampleSize ) );
    index[1] = static_cast< IndexValueType >( ( lateral / m_LateralAngularSeparation )
                                         + ( maxLateral / 2.0 ) );

    // Now, check to see if the index is within allowed bounds
    const bool isInside = region.IsInside(index);

    return isInside;
  }

  template< typename TCoordRep >
  void TransformContinuousIndexToPhysicalPoint(
    const ContinuousIndex< TCoordRep, VDimension > & index,
    Point< TCoordRep, VDimension > & point) const
  {
    const RegionType & region = this->GetLargestPossibleRegion();
    const double maxLateral = region.GetSize(1) - 1;

    // Convert the index into proper angular coordinates
    const TCoordRep radius  = ( index[0] * m_RadiusSampleSize ) + m_FirstSampleDistance;
    const TCoordRep lateral = ( index[1] - ( maxLateral / 2.0 ) ) * m_LateralAngularSeparation;

    // Convert the angular coordinates into Cartesian coordinates
    point[0] = static_cast< TCoordRep >( radius * std::sin(lateral) );
    point[1] = static_cast< TCoordRep >( radius * std::cos(lateral) );
  }

  template< typename TCoordRep >
  void TransformIndexToPhysicalPoint(
    const IndexType & index,
    Point< TCoordRep, VDimension > & point) const
  {
    const RegionType & region = this->GetLargestPossibleRegion();
    const double maxLateral = region.GetSize(1) - 1;

    // Convert the index into proper angular coordinates
    const TCoordRep radius  = ( index[0] * m_RadiusSampleSize ) + m_FirstSampleDistance;
    const TCoordRep lateral = ( index[1] - ( maxLateral / 2.0 ) ) * m_LateralAngularSeparation;

    // Convert the angular coordinates into Cartesian coordinates
    point[0] = static_cast< TCoordRep >( radius * std::sin(lateral) );
    point[1] = static_cast< TCoordRep >( radius * std::cos(lateral) );
  }

  itkSetMacro(LateralAngularSeparation, double);
  itkGetConstMacro(LateralAngularSeparation, double);

  itkSetMacro(RadiusSampleSize, double);
  itkGetConstMacro(RadiusSampleSize, double);

  itkSetMacro(FirstSampleDistance, double);
  itkGetConstMacro(FirstSampleDistance, double);

  template< typename TCoordRep >
  void TransformLocalVectorToPhysicalVector(
    FixedArray< TCoordRep, VDimension > & ) const
    {}

  template< typename TCoordRep >
  void TransformPhysicalVectorToLocalVector(
    const FixedArray< TCoordRep, VDimension > & ,
    FixedArray< TCoordRep, VDimension > & ) const
    {}

protected:
  CurvilinearArraySpecialCoordinatesImage()
    {
    m_RadiusSampleSize = 1;
    m_LateralAngularSeparation   =  1 * ( 2.0 * vnl_math::pi / 360.0 ); // 1
                                                                        // degree
    m_FirstSampleDistance = 0;
    }

  virtual ~CurvilinearArraySpecialCoordinatesImage() {}
  virtual void PrintSelf(std::ostream & os, Indent indent) const ITK_OVERRIDE;

private:
  CurvilinearArraySpecialCoordinatesImage(const Self &); // purposely not implemented
  void operator=(const Self &);                          // purposely not implemented

  double m_LateralAngularSeparation;    // in radians
  double m_RadiusSampleSize;
  double m_FirstSampleDistance;
};
} // end namespace itk

#ifndef ITK_MANUAL_INSTANTIATION
#include "itkCurvilinearArraySpecialCoordinatesImage.hxx"
#endif

#endif