text: /*=========================================================================
   *
   *  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.
   *
   *=========================================================================*/
  // Disable warning for long symbol names in this file only
  
  /*
  * This is a test file for the itkImageMaskSpatialObject class.
  * The suported pixel types does not include itkRGBPixel, itkRGBAPixel, etc...
  * So far it only allows to manage images of simple types like unsigned short,
  * unsigned int, or itk::Vector<...>.
  */
  
  /*
   * This test addresses bug
   * http://public.kitware.com/Bug/view.php?id=0006340
   *
   */
  
  #include "itkImageRegionIterator.h"
  #include "itkImageMaskSpatialObject.h"
  
  
  int itkImageMaskSpatialObjectTest2(int, char* [])
  {
    const unsigned int NDimensions = 3;
    int retval=EXIT_SUCCESS;
  
    typedef itk::ImageMaskSpatialObject<NDimensions> ImageMaskSpatialObject;
    typedef ImageMaskSpatialObject::PixelType        PixelType;
    typedef itk::Image<PixelType,NDimensions>        ImageType;
    typedef itk::Image<PixelType,NDimensions>        ImageType2;
    typedef itk::ImageRegionIterator<ImageType>      Iterator;
    typedef itk::ImageRegionIterator<ImageType2>     Iterator2;
  
    // Direction was not taken into account in the image spatial object
    // explicitly test using images with directions set.
    // Also explicitly uses nonzero origin, non identity scales
    // just to see
  
    {
    ImageType::Pointer image = ImageType::New();
    ImageType::SizeType size = {{ 50, 50, 50 }};
    ImageType::PointType origin;
    origin[0]=1.5;
    origin[1]=2.1;
    origin[2]=-3;
    image->SetOrigin(origin);
  
    ImageType::SpacingType spacing;
    spacing[0]=0.5;
    spacing[1]=0.7;
    spacing[2]=1.1;
    image->SetSpacing( spacing );
    ImageType::IndexType index = {{ 0, 0, 0 }};
  
    // Does it make sense to set the direction for a non-oriented image
    // well it happens frequently enough, so better test it
    ImageType::DirectionType direction;
    direction.Fill(0.0);
    direction[0][1]=-1;
    direction[1][0]=1;
    direction[2][2]=1;
    image->SetDirection(direction);
  
    ImageType::RegionType region;
  
    region.SetSize(size);
    region.SetIndex(index);
    image->SetRegions( region );
    image->Allocate(true); // initialize buffer to zero
  
    ImageType::RegionType insideRegion;
    const unsigned int INSIDE_SIZE = 30;
    const unsigned int INSIDE_INDEX = 10;
      {
      const ImageType::SizeType insideSize   = {{ INSIDE_SIZE, INSIDE_SIZE, INSIDE_SIZE }};
      insideRegion.SetSize( insideSize );
      }
      {
      const ImageType::IndexType insideIndex = {{ INSIDE_INDEX, INSIDE_INDEX, INSIDE_INDEX }};
      insideRegion.SetIndex( insideIndex );
      }
  
  
    Iterator it( image, insideRegion );
    it.GoToBegin();
  
    while( !it.IsAtEnd() )
      {
      it.Set( itk::NumericTraits< PixelType >::max() );
      ++it;
      }
  
    ImageMaskSpatialObject::Pointer maskSO = ImageMaskSpatialObject::New();
  
    maskSO->SetImage(image);
  
    maskSO->ComputeObjectToWorldTransform();
  
    Iterator itr( image, region );
    itr.GoToBegin();
  
    while( !itr.IsAtEnd() )
      {
      const ImageType::IndexType constIndex =  itr.GetIndex();
      const bool reference = insideRegion.IsInside( constIndex );
  
      ImageType::PointType point;
      image->TransformIndexToPhysicalPoint( constIndex, point );
      const bool test      = maskSO->IsInside( point );
      if( test != reference )
        {
        std::cerr << "Error in the evaluation of IsInside() " << std::endl;
        std::cerr << "Index failed = " << constIndex << std::endl;
        std::cerr << "Point failed = " << point << std::endl;
        std::cerr << "Image is a: "<<image->GetNameOfClass()<<std::endl;
        std::cerr << "Direction is: "<<std::endl<<image->GetDirection()<<std::endl;
        retval=EXIT_FAILURE;
        break;
        }
  
      ++itr;
      }
  
    if(retval==EXIT_SUCCESS)
      {
      std::cout<<"Test with "<<image->GetNameOfClass()<<" passed."<<std::endl;
      }
  
    // Check if insideregion is properly computed at the image boundary
    {
    ImageType::IndexType startPointIndex =
      {{ INSIDE_SIZE-2,
         INSIDE_SIZE-2,
         INSIDE_SIZE-2 }};
    ImageType::IndexType endPointIndex =
      {{ INSIDE_INDEX+INSIDE_SIZE+2,
         INSIDE_INDEX+INSIDE_SIZE+2,
         INSIDE_INDEX+INSIDE_SIZE+2 }};
    ImageType::PointType startPoint;
    ImageType::PointType endPoint;
    image->TransformIndexToPhysicalPoint( startPointIndex, startPoint );
    image->TransformIndexToPhysicalPoint( endPointIndex, endPoint );
  
    // Traverse along the line that goes through mask boundaries and
    // check if the value and the mask is consistent
    const int numberOfSteps = static_cast<int> (
      std::sqrt(double(INSIDE_SIZE*INSIDE_SIZE+
                      INSIDE_SIZE*INSIDE_SIZE+
                      INSIDE_SIZE*INSIDE_SIZE))*100.0
      );
    const ImageType::SpacingType incrementVector =
      (endPoint-startPoint)/static_cast<double>(numberOfSteps);
    ImageType::PointType point=startPoint;
    for (int i=0; i<numberOfSteps; i++)
      {
      point += incrementVector;
      const bool isInside = maskSO->IsInside( point );
      double value = itk::NumericTraits< PixelType >::ZeroValue();
      maskSO->ValueAt( point, value );
      const bool isZero = (value == itk::NumericTraits< PixelType >::ZeroValue());
      if( (isInside && isZero) || (!isInside && !isZero) )
        {
        ImageType::IndexType pointIndex;
        image->TransformPhysicalPointToIndex( point, pointIndex );
        std::cerr << "Error in the evaluation ValueAt and IsInside (all the points inside the mask shall have non-zero value) " << std::endl;
        std::cerr << "isInside = " << isInside << std::endl;
        std::cerr << "value = " << value << std::endl;
        std::cerr << "Index failed = " << pointIndex << std::endl;
        std::cerr << "Point failed = " << point << std::endl;
        std::cerr << "Image is a: "<< image->GetNameOfClass()<<std::endl;
        retval= EXIT_FAILURE;
        break;
        }
      }
    }
    }
  
    {
    ImageType2::Pointer image = ImageType2::New();
    ImageType2::SizeType size = {{ 50, 50, 50 }};
    ImageType::PointType origin;
    origin[0]=1.5;
    origin[1]=2.1;
    origin[2]=-3;
    image->SetOrigin(origin);
  
    ImageType::SpacingType spacing;
    spacing[0]=0.5;
    spacing[1]=0.7;
    spacing[2]=1.1;
    image->SetSpacing( spacing );
  
    ImageType2::IndexType index = {{ 0, 0, 0 }};
    ImageType2::DirectionType direction;
    direction.Fill(0.0);
    direction[0][1]=-1;
    direction[1][0]=1;
    direction[2][2]=1;
    image->SetDirection(direction);
  
    ImageType::RegionType region;
  
    region.SetSize(size);
    region.SetIndex(index);
  
    image->SetRegions( region );
    image->Allocate(true); // initialize buffer to zero
  
    ImageType::RegionType insideRegion;
    const unsigned int INSIDE_SIZE = 30;
    const unsigned int INSIDE_INDEX = 10;
      {
      const ImageType::SizeType insideSize   = {{ INSIDE_SIZE, INSIDE_SIZE, INSIDE_SIZE }};
      insideRegion.SetSize( insideSize );
      }
      {
      const ImageType::IndexType insideIndex = {{ INSIDE_INDEX, INSIDE_INDEX, INSIDE_INDEX }};
      insideRegion.SetIndex( insideIndex );
      }
  
  
    Iterator2 it( image, insideRegion );
    it.GoToBegin();
  
    while( !it.IsAtEnd() )
      {
      it.Set( itk::NumericTraits< PixelType >::max() );
      ++it;
      }
  
    ImageMaskSpatialObject::Pointer maskSO = ImageMaskSpatialObject::New();
  
    maskSO->SetImage(image);
  
    maskSO->ComputeObjectToWorldTransform();
  
    Iterator2 itr( image, region );
    itr.GoToBegin();
  
    while( !itr.IsAtEnd() )
      {
      const ImageType::IndexType constIndex =  itr.GetIndex();
      const bool reference = insideRegion.IsInside( constIndex );
      ImageType::PointType point;
      image->TransformIndexToPhysicalPoint( constIndex, point );
      const bool test      = maskSO->IsInside( point );
      if( test != reference )
        {
        std::cerr << "Error in the evaluation of IsInside() " << std::endl;
        std::cerr << "Index failed = " << constIndex << std::endl;
        std::cerr << "Point failed = " << point << std::endl;
        std::cerr << "Image is a: "<<image->GetNameOfClass()<<std::endl;
        std::cerr << "Direction is: "<<std::endl<<image->GetDirection()<<std::endl;
        retval= EXIT_FAILURE;
        break;
        }
        ++itr;
      }
    if(retval==EXIT_SUCCESS)
      {
      std::cout<<"Test with "<<image->GetNameOfClass()<<" passed."<<std::endl;
      }
    }
  
    return retval;
  }
  