itkPSMImplicitSurfaceDomain.h

/*=========================================================================
 *
 *  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 __itkPSMImplicitSurfaceDomain_h
#define __itkPSMImplicitSurfaceDomain_h

#include "itkPSMImageDomainWithCurvature.h"
#include <ctime>

namespace itk
{
template <class T, unsigned int VDimension>
class ITK_EXPORT PSMImplicitSurfaceDomain
  : public PSMImageDomainWithCurvature<T, VDimension>
{
public:
  typedef PSMImplicitSurfaceDomain Self;
  typedef PSMImageDomainWithCurvature<T, VDimension> Superclass;
  typedef SmartPointer<Self>  Pointer;
  typedef SmartPointer<const Self> ConstPointer;
  typedef WeakPointer<const Self>  ConstWeakPointer;
  typedef typename Superclass::ImageType ImageType;
  typedef typename Superclass::PointType PointType;
  typedef T PixelType;
  
  itkNewMacro(Self);

  itkTypeMacro(PSMImplicitSurfaceDomain, ParticleClipByRegionDomain);

  itkStaticConstMacro(Dimension, unsigned int, VDimension);

  itkSetMacro(Tolerance, T);
  itkGetMacro(Tolerance, T);
  
  virtual bool ApplyConstraints(PointType &p) const;

  virtual bool ApplyVectorConstraints(vnl_vector_fixed<double, VDimension> &gradE,
                                      const PointType &pos,
                                      double maxtimestep) const;

  virtual double Distance(const PointType &, const PointType &) const;

  void SetCuttingPlane(const vnl_vector<double> &a, const vnl_vector<double> &b,
                       const vnl_vector<double> &c);

  void RemoveCuttingPlane()  { m_UseCuttingPlane = false; }

  bool IsCuttingPlaneDefined() const {return m_UseCuttingPlane;}
  
  const vnl_vector_fixed<double, VDimension> &GetCuttingPlanePoint() const
  { return m_CuttingPlanePoint; }
  const vnl_vector_fixed<double, VDimension> &GetCuttingPlaneNormal() const
  { return m_CuttingPlaneNormal; }


  void AddSphere(const vnl_vector_fixed<double,VDimension> &v, double r)
  {
    if (r > 0)
    {
      m_SphereCenterList.push_back(v);
      m_SphereRadiusList.push_back(r);
    }
  }
  
  double GetSphereRadius(unsigned int i) const
  {
    if (m_SphereRadiusList.size() > i) return m_SphereRadiusList[i];
    else return 0.0;
  }

  vnl_vector_fixed<double, VDimension> GetSphereCenter(unsigned int i) const
  {
    if (m_SphereRadiusList.size() > i) return m_SphereCenterList[i];
    else
      {  return vnl_vector_fixed<double, VDimension>(0.0,0.0,0.0);  }
  }

  unsigned int GetNumberOfSpheres() const
  {
    return m_SphereCenterList.size();
  }
    
  
protected:
  PSMImplicitSurfaceDomain() : m_Tolerance(1.0e-4), m_UseCuttingPlane(false)
    {
      srand(time(0));
    }
  void PrintSelf(std::ostream& os, Indent indent) const
  {
    Superclass::PrintSelf(os, indent);
    os << indent << "m_Tolerance = " << m_Tolerance << std::endl;
  }
  virtual ~PSMImplicitSurfaceDomain() {};

private:
  PixelType m_Tolerance;
  bool m_UseCuttingPlane;
  vnl_vector_fixed<double, VDimension> m_CuttingPlanePoint;
  vnl_vector_fixed<double, VDimension> m_CuttingPlaneNormal;
  
  std::vector< vnl_vector_fixed<double, VDimension> > m_SphereCenterList;
  std::vector< double > m_SphereRadiusList;
  
  PSMImplicitSurfaceDomain(const Self&); //purposely not implemented
  void operator=(const Self&); //purposely not implemented
};

} // end namespace itk

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

#endif