RTK  1.4.0
Reconstruction Toolkit
rtkConjugateGradientGetP_kPlusOneImageFilter.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright RTK 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 
19 #ifndef rtkConjugateGradientGetP_kPlusOneImageFilter_h
20 #define rtkConjugateGradientGetP_kPlusOneImageFilter_h
21 
22 #include <itkImageToImageFilter.h>
23 #include <itkAddImageFilter.h>
24 #include <itkMultiplyImageFilter.h>
25 
26 namespace rtk
27 {
34 template< typename TInputImage>
35 class ConjugateGradientGetP_kPlusOneImageFilter : public itk::ImageToImageFilter< TInputImage, TInputImage>
36 {
37 public:
38 
43  typedef typename TInputImage::RegionType OutputImageRegionType;
44  typedef itk::Image<typename TInputImage::InternalPixelType,
45  TInputImage::ImageDimension> BetaImage;
46 
48  itkNewMacro(Self)
49 
50 
52 
54  void SetR_kPlusOne(const TInputImage* R_kPlusOne);
55  void SetRk(const TInputImage* Rk);
56  void SetPk(const TInputImage* Pk);
58 
59  itkSetMacro(SquaredNormR_k, double)
60  itkSetMacro(SquaredNormR_kPlusOne, double)
61 
63  typedef itk::AddImageFilter<TInputImage> AddFilterType;
64  typedef itk::MultiplyImageFilter<TInputImage, BetaImage, TInputImage> MultiplyFilterType;
65 
66 protected:
67  ConjugateGradientGetP_kPlusOneImageFilter();
68  virtual ~ConjugateGradientGetP_kPlusOneImageFilter() ITK_OVERRIDE {}
69 
70  typename TInputImage::Pointer GetR_kPlusOne();
71  typename TInputImage::Pointer GetRk();
72  typename TInputImage::Pointer GetPk();
73 
75  void GenerateData() ITK_OVERRIDE;
76 
77  void GenerateOutputInformation() ITK_OVERRIDE;
78 
79 private:
80  ConjugateGradientGetP_kPlusOneImageFilter(const Self &); //purposely not implemented
81  void operator=(const Self &); //purposely not implemented
82 
85  double m_Betak;
86 
90 
91 };
92 } //namespace ITK
93 
94 
95 #ifndef ITK_MANUAL_INSTANTIATION
96 #include "rtkConjugateGradientGetP_kPlusOneImageFilter.hxx"
97 #endif
98 
99 #endif
itk::Image< typename TInputImage::InternalPixelType, TInputImage::ImageDimension > BetaImage
itk::ImageToImageFilter< TInputImage, TInputImage > Superclass
void SetR_kPlusOne(const TInputImage *R_kPlusOne)
#define itkSetMacro(name, type)