otbWaveletOperatorBase.txx 5.56 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
/*=========================================================================

  Program:   ORFEO Toolbox
  Language:  C++
  Date:      $Date$
  Version:   $Revision$


  Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
  See OTBCopyright.txt for details.

Emmanuel Christophe's avatar
Emmanuel Christophe committed
12
  Copyright (c) Institut Telecom / Telecom Bretagne. All rights reserved.
13
14
  See ITCopyright.txt for details.

Emmanuel Christophe's avatar
Emmanuel Christophe committed
15
16
     This software is distributed WITHOUT ANY WARRANTY; without even
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
17
18
19
20
21
22
     PURPOSE.  See the above copyright notices for more information.

=========================================================================*/

#ifndef __otbWaveletOperatorBase_txx
#define __otbWaveletOperatorBase_txx
Emmanuel Christophe's avatar
Emmanuel Christophe committed
23

24
25
#include "otbWaveletOperatorBase.h"

Emmanuel Christophe's avatar
Emmanuel Christophe committed
26
27
#include <cassert>

28
29
namespace otb {

OTB Bot's avatar
STYLE    
OTB Bot committed
30
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
31
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
32
void
OTB Bot's avatar
STYLE    
OTB Bot committed
33
34
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::PrintSelf(std::ostream& os, itk::Indent i) const
Emmanuel Christophe's avatar
Emmanuel Christophe committed
35
{
36
37
38
39
  os << i << "Wavelet kind : " << GetWaveletName() << "\n";
  os << i << "Up-Sampling factor " << GetUpSampleFactor() << "\n";
  Superclass::PrintSelf(os, i.GetNextIndent());
  os << i << "Wavelet coeff: [ ";
OTB Bot's avatar
STYLE    
OTB Bot committed
40
41
42
43
  for (typename Superclass::ConstIterator iter = Superclass::Begin();
       iter != Superclass::End();
       ++iter)
    {
44
    os << *iter << ' ';
OTB Bot's avatar
STYLE    
OTB Bot committed
45
    }
46
47
48
  os << "]\n";
}

OTB Bot's avatar
STYLE    
OTB Bot committed
49
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
50
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
51
void
OTB Bot's avatar
STYLE    
OTB Bot committed
52
53
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::UpSamplingCoefficients(CoefficientVector& coeff)
54
{
OTB Bot's avatar
STYLE    
OTB Bot committed
55
  if (m_UpSampleFactor <= 1) return;
56

OTB Bot's avatar
STYLE    
OTB Bot committed
57
  unsigned long radius = static_cast<unsigned long>(coeff.size()) >> 1;
58
59
  unsigned long upSampleRadius = radius * m_UpSampleFactor;

60
61
  CoefficientVector upSampledCoeff ( 2 * upSampleRadius + 1 );
  upSampledCoeff.assign(2 * upSampleRadius + 1, 0.);
62
63
  upSampledCoeff[upSampleRadius] = coeff[radius];

64
  for (unsigned int i = 1; i <= radius; ++i)
OTB Bot's avatar
STYLE    
OTB Bot committed
65
66
67
68
    {
    upSampledCoeff[upSampleRadius + m_UpSampleFactor * i] = coeff[radius + i];
    upSampledCoeff[upSampleRadius - m_UpSampleFactor * i] = coeff[radius - i];
    }
69
70
71
  coeff = upSampledCoeff;
}

OTB Bot's avatar
STYLE    
OTB Bot committed
72
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
73
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
74
void
OTB Bot's avatar
STYLE    
OTB Bot committed
75
76
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::RevertFilter(CoefficientVector& coeff)
77
78
79
80
{
  const unsigned int length = coeff.size();
  const unsigned int medianPosition = length / 2;

OTB Bot's avatar
STYLE    
OTB Bot committed
81
82
83
84
85
86
87
  CoefficientVector newCoeff(length);
  newCoeff[medianPosition] = coeff[medianPosition];
  for (unsigned int i = 1; i <= medianPosition; i++)
    {
    newCoeff[medianPosition + i] = coeff[medianPosition - i];
    newCoeff[medianPosition - i] = coeff[medianPosition + i];
    }
88
89
90
91

  coeff = newCoeff;
}

OTB Bot's avatar
STYLE    
OTB Bot committed
92
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
93
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
94
void
OTB Bot's avatar
STYLE    
OTB Bot committed
95
96
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::GenerateForwardHighPassFilterFromLowPassFilter(CoefficientVector& coeff)
97
98
99
100
{
  const unsigned int length = coeff.size();
  //const unsigned int medianPosition = length / 2;

OTB Bot's avatar
STYLE    
OTB Bot committed
101
  CoefficientVector highPassCoeff(length + 2);
102
103
104
105
106

  highPassCoeff[0] = 0.;
  highPassCoeff[1] = 0.;

  double sign = -1; //(medianPosition&1)==1 ? -1. : 1.;
OTB Bot's avatar
STYLE    
OTB Bot committed
107
108
109
  for (unsigned int i = 0; i < length; i++)
    {
    highPassCoeff[i + 2] = sign * coeff[length - 1 - i];
110
    sign *= -1.;
OTB Bot's avatar
STYLE    
OTB Bot committed
111
    }
112
113
114
115
116
  //highPassCoeff[length] = 0.;
  //highPassCoeff[length+1] = 0.;

  coeff = highPassCoeff;

OTB Bot's avatar
STYLE    
OTB Bot committed
117
118
119
120
  while (coeff[0] == coeff[coeff.size() - 1])
    {
    ReduceFilterLength(coeff);
    }
121
122
}

OTB Bot's avatar
STYLE    
OTB Bot committed
123
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
124
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
125
void
OTB Bot's avatar
STYLE    
OTB Bot committed
126
127
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::GenerateInverseHighPassFilterFromLowPassFilter(CoefficientVector& coeff)
128
129
130
{
  const unsigned int length = coeff.size();

OTB Bot's avatar
STYLE    
OTB Bot committed
131
  CoefficientVector highPassCoeff(length + 2);
132

Emmanuel Christophe's avatar
Emmanuel Christophe committed
133
  double sign = -1;
OTB Bot's avatar
STYLE    
OTB Bot committed
134
135
  for (unsigned int i = 0; i < length; i++)
    {
136
137
    highPassCoeff[i] = sign * coeff[i];
    sign *= -1.;
OTB Bot's avatar
STYLE    
OTB Bot committed
138
    }
139
  highPassCoeff[length] = 0.;
OTB Bot's avatar
STYLE    
OTB Bot committed
140
  highPassCoeff[length + 1] = 0.;
141
142
143

  coeff = highPassCoeff;

144
  while ((coeff[0] == coeff[coeff.size() - 1]) && (coeff[0] == 0.0))
OTB Bot's avatar
STYLE    
OTB Bot committed
145
146
147
    {
    ReduceFilterLength(coeff);
    }
148
149
}

OTB Bot's avatar
STYLE    
OTB Bot committed
150
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
151
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
152
void
OTB Bot's avatar
STYLE    
OTB Bot committed
153
154
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::GenerateInverseLowPassFilterFromHighPassFilter(CoefficientVector& coeff)
155
156
157
{
  const unsigned int length = coeff.size();

OTB Bot's avatar
STYLE    
OTB Bot committed
158
  CoefficientVector highPassCoeff(length + 2);
159

Emmanuel Christophe's avatar
Emmanuel Christophe committed
160
  double sign = 1;
OTB Bot's avatar
STYLE    
OTB Bot committed
161
162
  for (unsigned int i = 0; i < length; i++)
    {
163
164
    highPassCoeff[i] = sign * coeff[i];
    sign *= -1.;
OTB Bot's avatar
STYLE    
OTB Bot committed
165
    }
166
  highPassCoeff[length] = 0.;
OTB Bot's avatar
STYLE    
OTB Bot committed
167
  highPassCoeff[length + 1] = 0.;
168
169
170

  coeff = highPassCoeff;

171
  while ((coeff[0] == coeff[coeff.size() - 1]) && (coeff[0] == 0.0))
OTB Bot's avatar
STYLE    
OTB Bot committed
172
    {
Emmanuel Christophe's avatar
Emmanuel Christophe committed
173
    assert(coeff.size() > 1);
OTB Bot's avatar
STYLE    
OTB Bot committed
174
175
    ReduceFilterLength(coeff);
    }
176
177
}

OTB Bot's avatar
STYLE    
OTB Bot committed
178
template <Wavelet::Wavelet TMotherWaveletOperator,
OTB Bot's avatar
STYLE    
OTB Bot committed
179
    class TPixel, unsigned int VDimension, class TAllocator>
Emmanuel Christophe's avatar
Emmanuel Christophe committed
180
void
OTB Bot's avatar
STYLE    
OTB Bot committed
181
182
WaveletOperatorBase<TMotherWaveletOperator, TPixel, VDimension, TAllocator>
::ReduceFilterLength(CoefficientVector& coeff)
183
184
{
  const unsigned int length = coeff.size();
Emmanuel Christophe's avatar
Emmanuel Christophe committed
185
  assert(length >= 2);
OTB Bot's avatar
STYLE    
OTB Bot committed
186
187
188
189
190
  CoefficientVector newFilter(length - 2);
  for (unsigned int i = 0; i < newFilter.size(); i++)
    {
    newFilter[i] = coeff[i + 1];
    }
191
192
193
194
195
196
  coeff = newFilter;
}

} // end of namespace otb

#endif