ANIMA  4.0
animaAtlasInitializer.hxx
Go to the documentation of this file.
2 
3 namespace anima
4 {
5 
6 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
8 {
9  this->SetNthInput(0, const_cast<TMaskImage*>(mask));
10 }
11 
12 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
14 {
15  this->SetNthInput(1, const_cast<TInputImage*>(image));
16  m_ImagesVector.push_back(image);
17 }
18 
19 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
21 {
22  this->SetNthInput(2, const_cast<TInputImage*>(image));
23  m_ImagesVector.push_back(image);
24 }
25 
26 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
28 {
29  this->SetNthInput(3, const_cast<TInputImage*>(image));
30  m_ImagesVector.push_back(image);
31 }
32 
33 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
35 {
36  this->SetNthInput(4, const_cast<TAtlasImage*>(image));
37  m_AtlasVector.push_back(image);
38 }
39 
40 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
42 {
43  this->SetNthInput(5, const_cast<TAtlasImage*>(image));
44  m_AtlasVector.push_back(image);
45 }
46 
47 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
49 {
50  this->SetNthInput(6, const_cast<TAtlasImage*>(image));
51  m_AtlasVector.push_back(image);
52 }
53 
54 
55 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
57 {
58  return static_cast< const TMaskImage * >
59  ( this->itk::ProcessObject::GetInput(0) );
60 }
61 
62 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
64 {
65  return static_cast< const TInputImage * >
66  ( this->itk::ProcessObject::GetInput(1) );
67 }
68 
69 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
71 {
72  return static_cast< const TInputImage * >
73  ( this->itk::ProcessObject::GetInput(2) );
74 }
75 
76 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
78 {
79  return static_cast< const TInputImage * >
80  ( this->itk::ProcessObject::GetInput(3) );
81 }
82 
83 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
85 {
86  return static_cast< const TAtlasImage * >
87  ( this->itk::ProcessObject::GetInput(4) );
88 }
89 
90 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
92 {
93  return static_cast< const TAtlasImage * >
94  ( this->itk::ProcessObject::GetInput(5) );
95 }
96 
97 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
99 {
100  return static_cast< const TAtlasImage * >
101  ( this->itk::ProcessObject::GetInput(6) );
102 }
103 
104 
105 
106 template <typename TInputImage, typename TMaskImage, typename TAtlasImage>
108 {
109  const unsigned int numberOfClasses = m_ImagesVector.size();
110 
111  /*--------------------Creates iterators --------------------*/
112  MaskConstIteratorType MaskIt(this->GetMask(),this->GetMask()->GetLargestPossibleRegion());
113 
114  std::vector<InputConstIteratorType> ImagesVectorVectorIt;
115  for ( unsigned int i = 0; i < numberOfClasses; i++ )
116  {
117  InputConstIteratorType It(m_ImagesVector[i],m_ImagesVector[i]->GetLargestPossibleRegion() );
118  ImagesVectorVectorIt.push_back(It);
119  }
120 
121  std::vector<AtlasConstIteratorType> AtlasVectorVectorIt;
122  for ( unsigned int i = 0; i < numberOfClasses; i++ )
123  {
124  AtlasConstIteratorType It(m_AtlasVector[i],m_AtlasVector[i]->GetLargestPossibleRegion() );
125  AtlasVectorVectorIt.push_back(It);
126  }
127 
128  std::vector<GaussianFunctionType::MeanVectorType> means(numberOfClasses,GaussianFunctionType::MeanVectorType( numberOfClasses ));
129  std::vector<DoubleVariableSizeMatrixType> covariances(numberOfClasses,DoubleVariableSizeMatrixType(numberOfClasses,numberOfClasses));
130 
131  m_Alphas.resize(numberOfClasses,0);
132 
133  /*--------------------Initialisation alpha, mean and covariances --------------------*/
134  for(unsigned int c = 0; c < numberOfClasses; c++)
135  {
136  for(unsigned int i = 0; i < numberOfClasses; i++)
137  {
138  means[c][i] = 0.0;
139  for(unsigned int j = 0; j < numberOfClasses; j++)
140  covariances[c](i,j) = 0.0;
141  }
142  }
143 
144  /*--------------------Computes Means --------------------*/
145  MaskIt.GoToBegin();
146  for ( unsigned int k = 0; k < numberOfClasses; k++ )
147  {
148  AtlasVectorVectorIt[k].GoToBegin();
149  ImagesVectorVectorIt[k].GoToBegin();
150  }
151 
152  while(!MaskIt.IsAtEnd())
153  {
154  if ( MaskIt.Get() != 0 )
155  {
156  for(unsigned int c = 0; c < numberOfClasses; c++)
157  {
158  //Getting factors of normalization for apriori m_ImagesVector
159  m_Alphas[c]+=static_cast<double>(AtlasVectorVectorIt[c].Get());
160 
161  //Mean values
162  for(unsigned int m = 0; m < numberOfClasses; m++)
163  {
164  (means[c])[m] += AtlasVectorVectorIt[c].Get() * ImagesVectorVectorIt[m].Get();
165  }
166  }
167  }
168  ++MaskIt;
169  for ( unsigned int k = 0; k < numberOfClasses; k++ )
170  {
171  ++AtlasVectorVectorIt[k];
172  ++ImagesVectorVectorIt[k];
173  }
174  }
175 
176  //normalization of mean
177  for(unsigned int c = 0; c < numberOfClasses; c++)
178  {
179  for(unsigned int m = 0; m < numberOfClasses; m++)
180  {
181  (means[c])[m]/=m_Alphas[c];
182  }
183  }
184 
185 
186  /*-------------------- Compute Covariances --------------------*/
187  MaskIt.GoToBegin();
188  for ( unsigned int i = 0; i < numberOfClasses; i++ )
189  {
190  ImagesVectorVectorIt[i].GoToBegin();
191  AtlasVectorVectorIt[i].GoToBegin();
192  }
193 
194  //Covariance
195  DoubleVariableSizeMatrixType x(numberOfClasses,1);
196  DoubleVariableSizeMatrixType xT (1,numberOfClasses);
197  DoubleVariableSizeMatrixType meanT(1,numberOfClasses);
198  DoubleVariableSizeMatrixType meanTTrans(numberOfClasses,1);
199 
200  while(!MaskIt.IsAtEnd())
201  {
202  if ( MaskIt.Get() != 0 )
203  {
204  for(unsigned int c = 0; c < numberOfClasses; c++)
205  {
206  for(unsigned int m = 0; m < numberOfClasses; m++)
207  {
208  x(m,0) = ImagesVectorVectorIt[m].Get();
209  xT(0,m) = x(m,0);
210  meanT(0,m) = (means[c])[m];
211  meanTTrans(m,0) = (means[c])[m];
212  }
213  covariances[c] += ( (x-meanTTrans) * (xT-meanT) ) * (static_cast<double>(AtlasVectorVectorIt[c].Get()) );
214  }
215  }
216  ++MaskIt;
217  for ( unsigned int i = 0; i < numberOfClasses; i++ )
218  {
219  ++ImagesVectorVectorIt[i];
220  ++AtlasVectorVectorIt[i];
221  }
222  }
223 
224 
225  /*-------------------- Compute Covariances Normalisation --------------------*/
226  for(unsigned int c = 0; c < numberOfClasses; c++)
227  {
228  //Covariance normalization
229  for(unsigned int m = 0; m < numberOfClasses; m++)
230  {
231  for(unsigned int n = 0; n < numberOfClasses; n++)
232  {
233  covariances[c](m,n) /= m_Alphas[c];
234  }
235  }
236  }
237 
238  double total=0.0;
239  for(unsigned int c = 0; c < numberOfClasses; c++)
240  {
241  total+=m_Alphas[c];
242  }
243 
244  for(unsigned int c = 0; c < numberOfClasses; c++)
245  {
246  GaussianFunctionType::Pointer gaussian = GaussianFunctionType::New();
247  gaussian->SetMean(means[c]);
248  gaussian->SetCovariance(covariances[c]);
249  m_GaussianModel.push_back(gaussian);
250  m_Alphas[c]/=total; // alpha normalization
251  }
252 
253 }
254 
255 }
void SetAtlasImage2(const TAtlasImage *image)
TInputImage::ConstPointer GetInputImage3()
void SetAtlasImage1(const TAtlasImage *image)
TInputImage::ConstPointer GetInputImage1()
void SetMask(const TMaskImage *mask)
void Update() ITK_OVERRIDE
calculate parameters of gaussians from images with the probability of the different classes [csf] [gm...
void SetInputImage1(const TInputImage *image)
void SetInputImage2(const TInputImage *image)
TMaskImage::ConstPointer GetMask()
itk::ImageRegionConstIterator< MaskImageType > MaskConstIteratorType
void SetInputImage3(const TInputImage *image)
TAtlasImage::ConstPointer GetAtlasImage3()
TAtlasImage::ConstPointer GetAtlasImage2()
void SetAtlasImage3(const TAtlasImage *image)
TInputImage::ConstPointer GetInputImage2()
itk::ImageRegionConstIterator< AtlasImageType > AtlasConstIteratorType
TAtlasImage::ConstPointer GetAtlasImage1()
itk::ImageRegionConstIterator< InputImageType > InputConstIteratorType
itk::VariableSizeMatrix< NumericType > DoubleVariableSizeMatrixType