#pragma once#include"Constraints.h"#include"DomainType.h"#include"itkDataObject.h"#include"itkPoint.h"#include"vnl/vnl_matrix_fixed.h"#include"vnl/vnl_vector_fixed.h"namespaceshapeworks{classParticleDomain{public:typedefstd::shared_ptr<ParticleDomain>Pointer;usingPointType=itk::Point<double,3>;usingGradNType=vnl_matrix_fixed<float,3,3>;usingVectorDoubleType=vnl_vector_fixed<double,3>;usingVectorFloatType=vnl_vector_fixed<float,3>;virtualboolApplyConstraints(PointType&p,intidx,booldbg=false)const=0;// todo update should be const?virtualPointTypeUpdateParticlePosition(constPointType&point,intidx,VectorDoubleType&update)const=0;virtualvoidInvalidateParticlePosition(intidx)const{}virtualVectorDoubleTypeProjectVectorToSurfaceTangent(VectorDoubleType&gradE,constPointType&pos,intidx)const=0;virtualVectorFloatTypeSampleGradientAtPoint(constPointType&point,intidx)const=0;virtualVectorFloatTypeSampleNormalAtPoint(constPointType&point,intidx)const=0;virtualGradNTypeSampleGradNAtPoint(constPointType&p,intidx)const=0;virtualdoubleDistance(constPointType&a,intidx_a,constPointType&b,intidx_b,VectorDoubleType*out_grad=nullptr)const{if(out_grad!=nullptr){for(inti=0;i<DIMENSION;i++){(*out_grad)[i]=a[i]-b[i];}}returna.EuclideanDistanceTo(b);}virtualdoubleSquaredDistance(constPointType&a,intidx_a,constPointType&b,intidx_b)const{returna.SquaredEuclideanDistanceTo(b);}virtualboolIsWithinDistance(constPointType&a,intidx_a,constPointType&b,intidx_b,doubletest_dist,double&distance)const{distance=this->Distance(a,idx_a,b,idx_b);returndistance<test_dist;}virtualdoubleGetCurvature(constPointType&p,intidx)const=0;virtualdoubleGetSurfaceMeanCurvature()const=0;virtualdoubleGetSurfaceStdDevCurvature()const=0;virtualconstPointType&GetLowerBound()const=0;virtualconstPointType&GetUpperBound()const=0;virtualPointTypeGetZeroCrossingPoint()const=0;virtualdoubleGetSurfaceArea()const=0;virtualPointTypeGetValidLocationNear(PointTypep)const=0;virtualdoubleGetMaxDiameter()const=0;virtualvoidDeleteImages()=0;virtualvoidDeletePartialDerivativeImages()=0;virtualvoidUpdateZeroCrossingPoint()=0;boolIsDomainFixed()const{returnm_FixedDomain;}virtualshapeworks::DomainTypeGetDomainType()const=0;std::shared_ptr<shapeworks::Constraints>GetConstraints()const{returnconstraints;}// Use `random` to advance a particle and return a new positionvirtualPointTypeGetPositionAfterSplit(constPointType&pt,constVectorDoubleType&local_direction,constVectorDoubleType&global_direction,doubleepsilon)const{// todo this has been copied from itkParticleSystem::AdvancedAllParticleSplitting.// Ideally, we should compute a direction that is "consistent" depending on the domain type and use the// `UpdateParticlePosition` API to advance the particle. See ContourDomain for an example. Leaving this be for// now because we'd have to retest all MeshDomain and ImageDomain use cases if this behaviour changes.PointTypenew_pt;for(unsignedintk=0;k<3;k++){new_pt[k]=pt[k]+epsilon*local_direction[k]/5.;}returnnew_pt;}voidSetDomainID(intid){this->m_DomainID=id;}voidSetDomainName(std::stringname){this->m_DomainName=name;}protected:// is this a fixed domain or not? We start as fixed and if an image or mesh is set, we set this to falseboolm_FixedDomain{true};intm_DomainID{-1};std::stringm_DomainName;ParticleDomain(){this->constraints=std::make_shared<shapeworks::Constraints>();}virtual~ParticleDomain(){}std::shared_ptr<shapeworks::Constraints>constraints;private:ParticleDomain(constParticleDomain&);// purposely not implementedvoidoperator=(constParticleDomain&);// purposely not implemented};}// namespace shapeworks