1 /*=========================================================================
2  *
3  *  Copyright Insight Software Consortium
4  *
5  *  Licensed under the Apache License, Version 2.0 (the "License");
6  *  you may not use this file except in compliance with the License.
7  *  You may obtain a copy of the License at
8  *
9  *         http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  *  Unless required by applicable law or agreed to in writing, software
12  *  distributed under the License is distributed on an "AS IS" BASIS,
13  *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  *  See the License for the specific language governing permissions and
15  *  limitations under the License.
16  *
17  *=========================================================================*/
18 #ifndef itkConicShellInteriorExteriorSpatialFunction_hxx
19 #define itkConicShellInteriorExteriorSpatialFunction_hxx
20 
21 #include "itkConicShellInteriorExteriorSpatialFunction.h"
22 
23 namespace itk
24 {
25 template< unsigned int VDimension, typename TInput >
26 ConicShellInteriorExteriorSpatialFunction< VDimension, TInput >
ConicShellInteriorExteriorSpatialFunction()27 ::ConicShellInteriorExteriorSpatialFunction()
28 
29 {
30   m_Origin.Fill(0.0);
31   m_OriginGradient.Fill(0.0);
32 }
33 
34 template< unsigned int VDimension, typename TInput >
35 void
36 ConicShellInteriorExteriorSpatialFunction< VDimension, TInput >
SetOriginGradient(GradientType grad)37 ::SetOriginGradient(GradientType grad)
38 {
39   m_OriginGradient = grad;
40 
41   // Normalize the origin gradient
42   m_OriginGradient.GetVnlVector().normalize();
43 }
44 
45 template< unsigned int VDimension, typename TInput >
46 typename ConicShellInteriorExteriorSpatialFunction< VDimension, TInput >
47 ::OutputType
48 ConicShellInteriorExteriorSpatialFunction< VDimension, TInput >
Evaluate(const InputType & position) const49 ::Evaluate(const InputType & position) const
50 {
51   using VectorType = Vector< double, VDimension >;
52 
53   // Compute the vector from the origin to the point being tested
54   VectorType vecOriginToTest = position - m_Origin;
55 
56   // Compute the length of this vector
57   // double vecDistance = vecOriginToTest.GetVnlVector().magnitude();
58   double vecDistance = vecOriginToTest.GetNorm();
59 
60   // Check to see if this an allowed distance
61   if ( !( ( vecDistance > m_DistanceMin ) && ( vecDistance < m_DistanceMax ) ) )
62     {
63     return 0; // not inside the conic shell
64     }
65   // Normalize it
66   // vecOriginToTest.GetVnlVector().normalize();
67   vecOriginToTest.Normalize();
68 
69   // Create a temp vector to get around const problems
70   GradientType originGradient = m_OriginGradient;
71 
72   // Now compute the dot product
73   double dotprod = originGradient * vecOriginToTest;
74 
75   if ( m_Polarity == 1 )
76     {
77     dotprod = dotprod * -1;
78     }
79 
80   // Check if it meets the angle criterion
81   OutputType result;
82   if ( dotprod > ( 1 - m_Epsilon ) )
83     {
84     result = 1; // it's inside the shell
85     }
86   else
87     {
88     result = 0; // it's not inside the shell
89     }
90 
91   return result;
92 }
93 
94 template< unsigned int VDimension, typename TInput >
95 void
96 ConicShellInteriorExteriorSpatialFunction< VDimension, TInput >
PrintSelf(std::ostream & os,Indent indent) const97 ::PrintSelf(std::ostream & os, Indent indent) const
98 {
99   Superclass::PrintSelf(os, indent);
100 
101   unsigned int i;
102   os << indent << "Origin: [";
103   for ( i = 0; i < VDimension - 1; i++ )
104     {
105     os << m_Origin[i] << ", ";
106     }
107   os << "]" << std::endl;
108 
109   os << indent << "Gradient at origin: [";
110   for ( i = 0; i < VDimension - 1; i++ )
111     {
112     os << m_OriginGradient[i] << ", ";
113     }
114   os << "]" << std::endl;
115 
116   os << indent << "DistanceMin: " << m_DistanceMin << std::endl;
117   os << indent << "DistanceMax: " << m_DistanceMax << std::endl;
118   os << indent << "Epsilon: " << m_Epsilon << std::endl;
119   os << indent << "Polarity: " << m_Polarity << std::endl;
120 }
121 } // end namespace itk
122 
123 #endif
124