OTB  9.0.0
Orfeo Toolbox
otbSparseUnmixingImageFilter.h
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2005-2022 Centre National d'Etudes Spatiales (CNES)
3  *
4  * This file is part of Orfeo Toolbox
5  *
6  * https://www.orfeo-toolbox.org/
7  *
8  * Licensed under the Apache License, Version 2.0 (the "License");
9  * you may not use this file except in compliance with the License.
10  * You may obtain a copy of the License at
11  *
12  * http://www.apache.org/licenses/LICENSE-2.0
13  *
14  * Unless required by applicable law or agreed to in writing, software
15  * distributed under the License is distributed on an "AS IS" BASIS,
16  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
17  * See the License for the specific language governing permissions and
18  * limitations under the License.
19  */
20 
21 #ifndef otbSparseUnmixingImageFilter_h
22 #define otbSparseUnmixingImageFilter_h
23 
24 #include "itkListSample.h"
25 #include "itkHistogram.h"
26 
27 #include "otbMacro.h"
28 #include "otbImage.h"
29 #include "otbWaveletOperator.h"
30 #include "otbWaveletFilterBank.h"
31 #include "otbWaveletTransform.h"
34 
35 namespace otb
36 {
37 
54 template <class TInputImage, class TOutputImage, unsigned int VNbInputImage, class TPrecision = double,
55  Wavelet::Wavelet TMotherWaveletOperator = Wavelet::SYMLET8>
56 class ITK_EXPORT SparseUnmixingImageFilter : public ImageToImageListFilter<TInputImage, TOutputImage>
57 {
58 public:
62  typedef itk::SmartPointer<Self> Pointer;
63  typedef itk::SmartPointer<const Self> ConstPointer;
64 
66  itkNewMacro(Self);
67 
70 
72  itkStaticConstMacro(InputImageDimension, unsigned int, TInputImage::ImageDimension);
73  itkStaticConstMacro(OutputImageDimension, unsigned int, TOutputImage::ImageDimension);
74  itkStaticConstMacro(NumberOfInputImages, unsigned int, VNbInputImage);
76 
78  typedef TInputImage InputImageType;
79  typedef TOutputImage OutputImageType;
80  typedef TPrecision PrecisionType;
81  static const Wavelet::Wavelet MotherWaveletOperatorID = TMotherWaveletOperator;
82 
85 
90  typedef typename WvltFilterType::OutputImageListType InternalImageListType;
93 
94  // typedef itk::FixedArray< PrecisionType, NumberOfInputImages-1 > AngleType;
95  typedef itk::Array<PrecisionType> AngleType;
96  typedef itk::Statistics::ListSample<AngleType> AngleListType;
97  typedef typename AngleListType::Pointer AngleListPointerType;
101 
102  typedef typename itk::Statistics::Histogram<PrecisionType> HistogramType;
103  typedef typename HistogramType::Pointer HistogramPointerType;
104  typedef typename HistogramType::SizeType HistogramSizeType;
105  typedef typename HistogramType::IndexType HistogramIndexType;
106  typedef typename HistogramType::MeasurementVectorType MeasurementVectorType;
107  typedef typename HistogramType::MeasurementType MeasurementType;
108 
114 
115  using Superclass::SetInput;
116  void SetInput(unsigned int i, const InputImageType*);
117  const InputImageType* GetInput(unsigned int i) const;
118 
119  void SetNumberOfDecomposition(unsigned int nb)
120  {
121  for (unsigned int i = 0; i < NumberOfInputImages; ++i)
122  {
123  m_WvltFilterList->GetNthElement(i)->SetNumberOfDecomposition(nb);
124  }
125  this->Modified();
126  }
127  unsigned int GetNumberOfDecomposition() const
128  {
129  return m_WvltFilterList->GetNthElement(0)->GetNumberOfDecomposition();
130  }
131 
133  {
134  m_AngleListFilter->SetThresholdValue(th);
135  this->Modified();
136  }
138  {
139  return m_AngleListFilter->GetThresholdValue();
140  }
141 
145  unsigned int GetNumberOfComponentsRequired() const
146  {
147  if (m_NumberOfComponentsRequired == 0)
148  GenerateNumberOfComponentsRequired();
149  return m_NumberOfComponentsRequired;
150  }
151  itkSetMacro(NumberOfComponentsRequired, unsigned int);
153 
154  itkSetMacro(NumberOfHistogramBins, unsigned int);
155  itkGetMacro(NumberOfHistogramBins, unsigned int);
156 
157  itkGetConstMacro(AngleList, AngleListType*);
158  itkGetConstMacro(WvltFilterList, WvltFilterListType*);
159  itkGetConstMacro(AngleListFilter, AngleListFilterType*);
160  itkGetConstMacro(Histogram, HistogramType*);
161  itkGetConstMacro(Transformer, TransformFilterType*);
162 
163 protected:
166  {
167  }
168 
169  void GenerateData() override;
170  virtual void GenerateNumberOfComponentsRequired();
171 
172 private:
173  SparseUnmixingImageFilter(const Self&) = delete;
174  void operator=(const Self&) = delete;
175 
179 
184 }; // end of class
185 
186 } // end of namespace otb
187 
188 #ifndef OTB_MANUAL_INSTANTIATION
190 #endif
191 
192 #endif
otbWaveletTransform.h
otb::SparseWvltToAngleMapperListFilter
This class select N-uple join-wvlt coeff for sparse unmixing.
Definition: otbSparseWvltToAngleMapperListFilter.h:47
otb::AngularProjectionSetImageFilter::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbAngularProjectionSetImageFilter.h:55
otb::SparseUnmixingImageFilter::m_NumberOfComponentsRequired
unsigned int m_NumberOfComponentsRequired
Definition: otbSparseUnmixingImageFilter.h:176
otb::WaveletTransform
Wavelet transformation framework.
Definition: otbWaveletTransform.h:58
otb::SparseUnmixingImageFilter::m_AngleList
AngleListPointerType m_AngleList
Definition: otbSparseUnmixingImageFilter.h:178
otb::SparseUnmixingImageFilter::AngleListFilterType
SparseWvltToAngleMapperListFilter< InternalImageListType, AngleListType, NumberOfInputImages > AngleListFilterType
Definition: otbSparseUnmixingImageFilter.h:98
otbAngularProjectionSetImageFilter.h
otb::SparseUnmixingImageFilter::FilterBankType
WaveletFilterBank< InternalImageType, InternalImageType, WaveletOperatorType, Wavelet::FORWARD > FilterBankType
Definition: otbSparseUnmixingImageFilter.h:87
otb::SparseUnmixingImageFilter::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbSparseUnmixingImageFilter.h:62
otb::SparseUnmixingImageFilter::GetNumberOfDecomposition
unsigned int GetNumberOfDecomposition() const
Definition: otbSparseUnmixingImageFilter.h:127
otb::SparseUnmixingImageFilter::WvltFilterListType
ObjectList< WvltFilterType > WvltFilterListType
Definition: otbSparseUnmixingImageFilter.h:91
otb::SparseUnmixingImageFilter::PrecisionType
TPrecision PrecisionType
Definition: otbSparseUnmixingImageFilter.h:80
otb::SparseUnmixingImageFilter::InternalImageListType
WvltFilterType::OutputImageListType InternalImageListType
Definition: otbSparseUnmixingImageFilter.h:90
otb::SparseUnmixingImageFilter::SetThresholdValue
void SetThresholdValue(PrecisionType th)
Definition: otbSparseUnmixingImageFilter.h:132
otb::SparseUnmixingImageFilter::OutputImageType
TOutputImage OutputImageType
Definition: otbSparseUnmixingImageFilter.h:79
otb::SparseUnmixingImageFilter::InputImageType
TInputImage InputImageType
Definition: otbSparseUnmixingImageFilter.h:78
otbImage.h
otb::ImageListSource
Base class for all the filters producing an otbImageList.
Definition: otbImageListSource.h:40
otb::Wavelet::Wavelet
Wavelet
Definition: otbWaveletGenerator.h:35
otb::Wavelet::SYMLET8
@ SYMLET8
Definition: otbWaveletGenerator.h:50
otb::SparseUnmixingImageFilter::Superclass
ImageToImageListFilter< TInputImage, TOutputImage > Superclass
Definition: otbSparseUnmixingImageFilter.h:61
otb::SparseUnmixingImageFilter::HistogramType
itk::Statistics::Histogram< PrecisionType > HistogramType
Definition: otbSparseUnmixingImageFilter.h:102
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::Image
Creation of an "otb" image which contains metadata.
Definition: otbImage.h:89
otb::SparseUnmixingImageFilter::WvltFilterListPointerType
WvltFilterListType::Pointer WvltFilterListPointerType
Definition: otbSparseUnmixingImageFilter.h:92
otb::WaveletOperator
A Generic NeighborhoodOperator wavelets filter set defined for templation.
Definition: otbWaveletOperator.h:56
otbMacro.h
otb::AngularProjectionSetImageFilter
Performs spherical transformation in ND space from a set of angle values.
Definition: otbAngularProjectionSetImageFilter.h:49
otb::SparseUnmixingImageFilter::m_AngleListFilter
AngleListFilterPointerType m_AngleListFilter
Definition: otbSparseUnmixingImageFilter.h:181
otb::AngularProjectionSetImageFilter::OutputImageListPointerType
OutputImageListType::Pointer OutputImageListPointerType
Definition: otbAngularProjectionSetImageFilter.h:84
otb::SparseWvltToAngleMapperListFilter::OutputSampleListType
TOutputSampleList OutputSampleListType
Definition: otbSparseWvltToAngleMapperListFilter.h:77
otb::SparseUnmixingImageFilter::ConstPointer
itk::SmartPointer< const Self > ConstPointer
Definition: otbSparseUnmixingImageFilter.h:63
otb::SparseUnmixingImageFilter::HistogramSizeType
HistogramType::SizeType HistogramSizeType
Definition: otbSparseUnmixingImageFilter.h:104
otb::ObjectList::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbObjectList.h:46
otb::SparseUnmixingImageFilter::AngleListFilterPointerType
AngleListFilterType::Pointer AngleListFilterPointerType
Definition: otbSparseUnmixingImageFilter.h:99
otb::SparseUnmixingImageFilter::m_Transformer
TransformFilterPointerType m_Transformer
Definition: otbSparseUnmixingImageFilter.h:183
otb::SparseUnmixingImageFilter::AngleType
itk::Array< PrecisionType > AngleType
Definition: otbSparseUnmixingImageFilter.h:95
otb::SparseUnmixingImageFilter::GetNumberOfComponentsRequired
unsigned int GetNumberOfComponentsRequired() const
Definition: otbSparseUnmixingImageFilter.h:145
otb::SparseUnmixingImageFilter::AngleListPointerType
AngleListType::Pointer AngleListPointerType
Definition: otbSparseUnmixingImageFilter.h:97
otb::AngularProjectionSetImageFilter::OutputImageListType
Superclass::OutputImageListType OutputImageListType
Definition: otbAngularProjectionSetImageFilter.h:83
otb::SparseUnmixingImageFilter::InternalSampleListType
AngleListFilterType::OutputSampleListType InternalSampleListType
Definition: otbSparseUnmixingImageFilter.h:100
otbWaveletOperator.h
otb::SparseUnmixingImageFilter::m_WvltFilterList
WvltFilterListPointerType m_WvltFilterList
Definition: otbSparseUnmixingImageFilter.h:180
otb::SparseUnmixingImageFilter::WvltFilterPointerType
WvltFilterType::Pointer WvltFilterPointerType
Definition: otbSparseUnmixingImageFilter.h:89
otb::SparseUnmixingImageFilter::MeasurementVectorType
HistogramType::MeasurementVectorType MeasurementVectorType
Definition: otbSparseUnmixingImageFilter.h:106
otb::SparseUnmixingImageFilter::m_NumberOfHistogramBins
unsigned int m_NumberOfHistogramBins
Definition: otbSparseUnmixingImageFilter.h:177
otb::SparseUnmixingImageFilter::AngleListType
itk::Statistics::ListSample< AngleType > AngleListType
Definition: otbSparseUnmixingImageFilter.h:96
otb::SparseUnmixingImageFilter::SetNumberOfDecomposition
void SetNumberOfDecomposition(unsigned int nb)
Definition: otbSparseUnmixingImageFilter.h:119
otbWaveletFilterBank.h
otb::SparseWvltToAngleMapperListFilter::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbSparseWvltToAngleMapperListFilter.h:53
otb::WaveletTransform::Pointer
itk::SmartPointer< Self > Pointer
Definition: otbWaveletTransform.h:64
otb::SparseUnmixingImageFilter
This class detects linear dependencies from N wavelet decompositions.
Definition: otbSparseUnmixingImageFilter.h:56
otb::SparseUnmixingImageFilter::OutputImageIterator
TransformFilterType::OutputImageIterator OutputImageIterator
Definition: otbSparseUnmixingImageFilter.h:113
otb::SparseUnmixingImageFilter::Self
SparseUnmixingImageFilter Self
Definition: otbSparseUnmixingImageFilter.h:60
otb::AngularProjectionSetImageFilter::OutputImageIterator
OutputImageListType::Iterator OutputImageIterator
Definition: otbAngularProjectionSetImageFilter.h:85
otbSparseWvltToAngleMapperListFilter.h
otb::SparseUnmixingImageFilter::m_Histogram
HistogramPointerType m_Histogram
Definition: otbSparseUnmixingImageFilter.h:182
otbSparseUnmixingImageFilter.hxx
otb::WaveletFilterBank
One level stationary wavelet transform.
Definition: otbWaveletFilterBank.h:87
otb::SparseUnmixingImageFilter::OutputImageListPointerType
TransformFilterType::OutputImageListPointerType OutputImageListPointerType
Definition: otbSparseUnmixingImageFilter.h:112
otb::SparseUnmixingImageFilter::HistogramIndexType
HistogramType::IndexType HistogramIndexType
Definition: otbSparseUnmixingImageFilter.h:105
otb::SparseUnmixingImageFilter::~SparseUnmixingImageFilter
~SparseUnmixingImageFilter() override
Definition: otbSparseUnmixingImageFilter.h:165
otb::SparseUnmixingImageFilter::HistogramPointerType
HistogramType::Pointer HistogramPointerType
Definition: otbSparseUnmixingImageFilter.h:103
otb::ObjectList
This class is a generic all-purpose wrapping around an std::vector<itk::SmartPointer<ObjectType> >.
Definition: otbObjectList.h:40
otb::SparseUnmixingImageFilter::WaveletOperatorType
WaveletOperator< MotherWaveletOperatorID, Wavelet::FORWARD, PrecisionType, InputImageDimension > WaveletOperatorType
Definition: otbSparseUnmixingImageFilter.h:86
otb::ImageToImageListFilter::InputImageType
TInputImage InputImageType
Definition: otbImageToImageListFilter.h:52
otb::SparseUnmixingImageFilter::WvltFilterType
WaveletTransform< InternalImageType, InternalImageType, FilterBankType, Wavelet::FORWARD > WvltFilterType
Definition: otbSparseUnmixingImageFilter.h:88
otb::SparseUnmixingImageFilter::OutputImageListType
TransformFilterType::OutputImageListType OutputImageListType
Definition: otbSparseUnmixingImageFilter.h:111
otb::ImageToImageListFilter
Base class for all the filters taking an image input to produce an image list.
Definition: otbImageToImageListFilter.h:39
otb::SparseUnmixingImageFilter::GetThresholdValue
PrecisionType GetThresholdValue() const
Definition: otbSparseUnmixingImageFilter.h:137
otb::SparseUnmixingImageFilter::TransformFilterPointerType
TransformFilterType::Pointer TransformFilterPointerType
Definition: otbSparseUnmixingImageFilter.h:110
otb::SparseUnmixingImageFilter::InternalImageType
Image< PrecisionType, InputImageDimension > InternalImageType
Definition: otbSparseUnmixingImageFilter.h:84
otb::SparseUnmixingImageFilter::TransformFilterType
AngularProjectionSetImageFilter< InputImageType, OutputImageType, AngleListType, PrecisionType > TransformFilterType
Definition: otbSparseUnmixingImageFilter.h:109
otb::SparseUnmixingImageFilter::MeasurementType
HistogramType::MeasurementType MeasurementType
Definition: otbSparseUnmixingImageFilter.h:107