OTB  9.0.0
Orfeo Toolbox
otbSFSTexturesFunctor.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 otbSFSTexturesFunctor_h
22 #define otbSFSTexturesFunctor_h
23 
24 #include "otbMath.h"
25 #include "itkNumericTraits.h"
26 #include <vector>
27 
28 namespace otb
29 {
56 namespace Functor
57 {
58 template <class TIter, class TOutputValue>
60 {
61 public:
63  :
64  m_SpatialThreshold(100),
67  m_Alpha(1)
68  {
69  this->SetNumberOfDirections(20); // set the step too
70  m_SelectedTextures = std::vector<bool>(6, 1);
71  }
72 
74  {
75  }
76 
77  typedef typename TIter::InternalPixelType InternalPixelType;
78  typedef typename TIter::SizeType SizeType;
79  typedef typename TIter::IndexType IndexType;
80  typedef typename TIter::OffsetType OffsetType;
81  typedef TOutputValue OutputValueType;
82  typedef std::vector<OutputValueType> OutputType;
83 
84  void SetSpatialThreshold(unsigned int thresh)
85  {
86  m_SpatialThreshold = thresh;
87  }
89  {
90  m_SpectralThreshold = thresh;
91  }
92  void SetRatioMaxConsiderationNumber(unsigned int value)
93  {
95  }
96  void SetAlpha(double alpha)
97  {
98  m_Alpha = alpha;
99  }
100  void SetNumberOfDirections(unsigned int D)
101  {
103  m_DirectionStep = CONST_PI / static_cast<double>(D);
104  }
105  void SetDirectionStep(double step)
106  {
107  m_DirectionStep = step;
108  }
109  void SetSelectedTextures(std::vector<bool>& vect)
110  {
111  m_SelectedTextures.clear();
112  m_SelectedTextures = vect;
113  }
114  void SetTextureStatus(unsigned int id, bool isSelected)
115  {
116  m_SelectedTextures[id] = isSelected;
117  }
118 
119  // Define some getters using itkGetMacro
120  itkGetMacro(SpatialThreshold, unsigned int);
121  itkGetMacro(SpectralThreshold, InternalPixelType);
122  itkGetMacro(RatioMaxConsiderationNumber, unsigned int);
123  itkGetMacro(Alpha, double);
124  itkGetMacro(NumberOfDirections, unsigned int);
125 
126  std::vector<bool> GetTexturesStatus()
127  {
128  return m_SelectedTextures;
129  }
130 
131  inline OutputType operator()(const TIter& it)
132  {
133  double length = itk::NumericTraits<double>::NonpositiveMin();
134  double width = itk::NumericTraits<double>::max();
135  double SpatialThresholdDouble = static_cast<double>(m_SpatialThreshold);
136  double NumberOfDirectionsDouble = static_cast<double>(m_NumberOfDirections);
137  double dist = 0.;
138  double angle = 0.;
139  double sdiVal = 0.;
140  double sumWMean = 0.;
141  double sum = 0.;
142  std::vector<double> di(m_NumberOfDirections, 0.);
143  std::vector<double> minSorted(m_RatioMaxConsiderationNumber, width);
144  std::vector<double> maxSorted(m_RatioMaxConsiderationNumber, length);
145  std::vector<double> sti(m_NumberOfDirections, 0.);
146  std::vector<unsigned int> lengthLine(m_NumberOfDirections, 0);
147 
148  std::vector<double>::iterator itVector;
149  OutputType out(6, 0);
150 
151  OffsetType off;
152  off.Fill(0);
153 
154  for (unsigned int d = 0; d < m_NumberOfDirections; d++)
155  {
156  // Current angle direction
157  angle = m_DirectionStep * static_cast<double>(d);
158 
159  // last offset in the direction respecting spatial threshold
160  off[0] = static_cast<int>(std::floor(SpatialThresholdDouble * std::cos(angle) + 0.5));
161  off[1] = static_cast<int>(std::floor(SpatialThresholdDouble * std::sin(angle) + 0.5));
162  // last indices in the direction respecting spectral threshold
163  OffsetType offEnd = this->FindLastOffset(it, off);
164 
165  // Check the opposite side
166  off[0] *= -1.0;
167  off[1] *= -1.0;
168  OffsetType offStart = this->FindLastOffset(it, off);
169 
170  // computes distance = dist between the 2 segment point.
171  dist = std::sqrt(std::pow(static_cast<double>(offEnd[0] - offStart[0]), 2) + std::pow(static_cast<double>(offEnd[1] - offStart[1]), 2));
172 
173  // for length computation
174  if (m_SelectedTextures[0] == true)
175  if (dist > length)
176  length = dist;
177 
178  // for width computation
179  if (m_SelectedTextures[1] == true)
180  if (dist < width)
181  width = dist;
182 
183  // for PSI computation
184  if (m_SelectedTextures[2] == true || m_SelectedTextures[5] == true)
185  sum += dist;
186 
187  // for w-mean computation
188  if (m_SelectedTextures[3] == true)
189  sdiVal = this->ComputeSDi(it, offEnd);
190 
191  // for Ratio computation
192  if (m_SelectedTextures[4] == true)
193  {
194  bool doo = false;
195  itVector = maxSorted.begin();
196  while (itVector != maxSorted.end() && doo == false)
197  {
198  if (dist > (*itVector))
199  {
200  itVector = maxSorted.insert(itVector, dist);
201  maxSorted.pop_back();
202  doo = true;
203  }
204  ++itVector;
205  }
206  doo = false;
207  itVector = minSorted.begin();
208  while (itVector != minSorted.end() && doo == false)
209  {
210  if (dist < (*itVector))
211  {
212  itVector = minSorted.insert(itVector, dist);
213  minSorted.pop_back();
214  doo = true;
215  }
216  ++itVector;
217  }
218  }
219 
220  di[d] = dist;
221  if (m_SelectedTextures[3] == true)
222  {
223  lengthLine[d] = static_cast<unsigned int>(
224  dist); // static_cast<unsigned int>( std::sqrt(std::pow(static_cast<double>(offEnd[0]), 2) + std::pow(static_cast<double>(offEnd[1]), 2)) );
225  sti[d] = sdiVal;
226  if (sdiVal != 0.)
227  sumWMean += (m_Alpha * (dist - 1) * dist /*lengthLine[n]*di[n]*/) / sdiVal;
228  }
229  }
230 
232  // length
233  if (m_SelectedTextures[0] == true)
234  out[0] = static_cast<OutputValueType>(length);
235  // width
236  if (m_SelectedTextures[1] == true)
237  out[1] = static_cast<OutputValueType>(width);
238  // PSI
239  if (m_SelectedTextures[2] == true)
240  out[2] = static_cast<OutputValueType>(sum / NumberOfDirectionsDouble);
241  // w-mean
242  if (m_SelectedTextures[3] == true)
243  out[3] = static_cast<OutputValueType>(sumWMean / NumberOfDirectionsDouble);
244  // ratio
245  if (m_SelectedTextures[4] == true)
246  {
247  double sumMin = 0;
248  double sumMax = 0;
249  for (unsigned int t = 0; t < m_RatioMaxConsiderationNumber; t++)
250  {
251  sumMin += minSorted[t];
252  sumMax += maxSorted[t];
253  }
254  if (sumMax != 0.)
255  out[4] = static_cast<OutputValueType>(std::atan(sumMin / sumMax));
256  else if (sumMax == 0. && sumMin == 0.)
257  out[4] = static_cast<OutputValueType>(1.);
258  }
259  // SD
260  if (m_SelectedTextures[5] == true)
261  {
262  double sumPSI = 0;
263  for (unsigned int n = 0; n < di.size(); ++n)
264  sumPSI += std::pow(di[n] - sumWMean / NumberOfDirectionsDouble, 2);
265  out[5] = static_cast<OutputValueType>(std::sqrt(sumPSI) / (NumberOfDirectionsDouble - 1.));
266  }
267 
268  return out;
269  }
270 
275  OffsetType FindLastOffset(const TIter& it, const OffsetType& stopOffset)
276  {
277  bool res = true;
278  int signX = this->ComputeStep(stopOffset[0]);
279  int signY = this->ComputeStep(stopOffset[1]);
281 
282  OffsetType currentOff;
283  currentOff.Fill(0);
284  currentOff[0] = signX;
285 
286  double slop = 0.;
287  if (stopOffset[0] != 0)
288  slop = static_cast<double>(stopOffset[1] / static_cast<double>(stopOffset[0]));
289 
290  bool isInside = true;
291  while (isInside == true && res == true)
292  {
293  this->ComputePointLine(currentOff, slop, signY, stopOffset[0]);
294 
295  if (std::abs(it.GetPixel(currentOff) - it.GetCenterPixel()) > m_SpectralThreshold)
296  {
297  res = false;
298  }
299  else
300  currentOff[0] += signX;
301 
302  isInside = this->CheckIsInside(signX, signY, currentOff, stopOffset);
303  }
304 
305  return currentOff;
306  }
307 
309  double ComputeSDi(const TIter& it, const OffsetType& stopOffset)
310  {
311  bool canGo = true;
312  unsigned int nbElt = 0;
313  double SDi = 0.;
314  double mean = 0.;
315  double slop = 0.;
316  if (stopOffset[0] != 0)
317  slop = static_cast<double>(stopOffset[1] / static_cast<double>(stopOffset[0]));
319 
320  int signX = this->ComputeStep(stopOffset[0]);
321  int signY = this->ComputeStep(stopOffset[1]);
322 
323  OffsetType currentOff;
324  currentOff.Fill(0);
325  currentOff[0] = signX;
326 
327  bool isInside = true;
328  // First compute mean
329  while (isInside == true && canGo == true)
330  {
331  this->ComputePointLine(currentOff, slop, signY, stopOffset[0]);
332 
333  mean += static_cast<double>(it.GetPixel(currentOff));
334  nbElt++;
335 
336  if (std::abs(it.GetPixel(currentOff) - it.GetCenterPixel()) >= m_SpectralThreshold)
337  canGo = false;
338  else
339  currentOff[0] += signX;
340 
341  isInside = this->CheckIsInside(signX, signY, currentOff, stopOffset);
342  }
343 
344  mean /= static_cast<double>(nbElt);
345  currentOff[0] = signX;
346  currentOff[1] = 0;
347  isInside = true;
348 
349  while (isInside == true && canGo == true)
350  {
351  this->ComputePointLine(currentOff, slop, signY, stopOffset[0]);
352 
353  SDi += std::pow((static_cast<double>(it.GetPixel(currentOff)) - mean), 2);
354  if (std::abs(it.GetPixel(currentOff) - it.GetCenterPixel()) >= m_SpectralThreshold)
355  canGo = false;
356  else
357  currentOff[0] += signX;
358 
359  isInside = this->CheckIsInside(signX, signY, currentOff, stopOffset);
360  }
361  return std::sqrt(SDi);
362  }
363 
365  bool CheckIsInside(const int& signX, const int& signY, const OffsetType& currentOff, const OffsetType& stopOffset)
366  {
367  bool isInside = true;
368  if (signX * currentOff[0] >= signX * stopOffset[0] && stopOffset[0] != 0)
369  isInside = false;
370  else if (signY * currentOff[1] >= signY * stopOffset[1] && stopOffset[1] != 0)
371  isInside = false;
373 
374  return isInside;
375  }
376 
381  void ComputePointLine(OffsetType& currentOff, const double& slop, const int& signY, const int& stopOffsetX)
382  {
383  if (stopOffsetX != 0)
384  currentOff[1] = static_cast<int>(std::floor(slop * static_cast<double>(currentOff[0]) + 0.5));
385  else
386  currentOff[1] += signY;
387  }
389 
393  int ComputeStep(const int& stopOffset)
394  {
395  int sign = 1;
396  if (stopOffset < 0)
397  sign = -1;
399 
400  return sign;
401  }
402 
403 protected:
405  unsigned int m_SpatialThreshold;
406 
409 
412 
414  double m_Alpha;
415 
417  unsigned int m_NumberOfDirections;
418 
421 
431  std::vector<bool> m_SelectedTextures;
432 };
433 
434 } // end namespace functor
435 } // end namespace otb
436 
437 #endif
otb::Functor::SFSTexturesFunctor::FindLastOffset
OffsetType FindLastOffset(const TIter &it, const OffsetType &stopOffset)
Definition: otbSFSTexturesFunctor.h:275
otb::Functor::SFSTexturesFunctor::InternalPixelType
TIter::InternalPixelType InternalPixelType
Definition: otbSFSTexturesFunctor.h:77
otb::Functor::SFSTexturesFunctor::SFSTexturesFunctor
SFSTexturesFunctor()
Definition: otbSFSTexturesFunctor.h:62
otb::mean
Definition: otbParserXPlugins.h:261
otb::CONST_PI
constexpr double CONST_PI
Definition: otbMath.h:49
otb::Functor::SFSTexturesFunctor::ComputeStep
int ComputeStep(const int &stopOffset)
Definition: otbSFSTexturesFunctor.h:393
otb::Functor::SFSTexturesFunctor::operator()
OutputType operator()(const TIter &it)
Definition: otbSFSTexturesFunctor.h:131
otb::Functor::SFSTexturesFunctor::m_NumberOfDirections
unsigned int m_NumberOfDirections
Definition: otbSFSTexturesFunctor.h:417
otb::Functor::SFSTexturesFunctor::m_RatioMaxConsiderationNumber
unsigned int m_RatioMaxConsiderationNumber
Definition: otbSFSTexturesFunctor.h:411
otb::Functor::SFSTexturesFunctor::IndexType
TIter::IndexType IndexType
Definition: otbSFSTexturesFunctor.h:79
otb::Functor::SFSTexturesFunctor::SetSelectedTextures
void SetSelectedTextures(std::vector< bool > &vect)
Definition: otbSFSTexturesFunctor.h:109
otbMath.h
otb
The "otb" namespace contains all Orfeo Toolbox (OTB) classes.
Definition: otbJoinContainer.h:32
otb::Functor::SFSTexturesFunctor::OffsetType
TIter::OffsetType OffsetType
Definition: otbSFSTexturesFunctor.h:80
otb::Functor::SFSTexturesFunctor::SetRatioMaxConsiderationNumber
void SetRatioMaxConsiderationNumber(unsigned int value)
Definition: otbSFSTexturesFunctor.h:92
otb::Functor::SFSTexturesFunctor::ComputePointLine
void ComputePointLine(OffsetType &currentOff, const double &slop, const int &signY, const int &stopOffsetX)
Definition: otbSFSTexturesFunctor.h:381
otb::Functor::SFSTexturesFunctor::GetTexturesStatus
std::vector< bool > GetTexturesStatus()
Definition: otbSFSTexturesFunctor.h:126
otb::Functor::SFSTexturesFunctor::m_DirectionStep
double m_DirectionStep
Definition: otbSFSTexturesFunctor.h:420
otb::Functor::SFSTexturesFunctor::m_Alpha
double m_Alpha
Definition: otbSFSTexturesFunctor.h:414
otb::Functor::SFSTexturesFunctor::SetNumberOfDirections
void SetNumberOfDirections(unsigned int D)
Definition: otbSFSTexturesFunctor.h:100
otb::Functor::SFSTexturesFunctor::~SFSTexturesFunctor
virtual ~SFSTexturesFunctor()
Definition: otbSFSTexturesFunctor.h:73
otb::Functor::SFSTexturesFunctor
Definition: otbSFSTexturesFunctor.h:59
otb::Functor::SFSTexturesFunctor::SizeType
TIter::SizeType SizeType
Definition: otbSFSTexturesFunctor.h:78
otb::Functor::SFSTexturesFunctor::SetSpatialThreshold
void SetSpatialThreshold(unsigned int thresh)
Definition: otbSFSTexturesFunctor.h:84
otb::Functor::SFSTexturesFunctor::OutputType
std::vector< OutputValueType > OutputType
Definition: otbSFSTexturesFunctor.h:82
otb::Functor::SFSTexturesFunctor::m_SelectedTextures
std::vector< bool > m_SelectedTextures
Definition: otbSFSTexturesFunctor.h:431
otb::Functor::SFSTexturesFunctor::SetTextureStatus
void SetTextureStatus(unsigned int id, bool isSelected)
Definition: otbSFSTexturesFunctor.h:114
otb::Functor::SFSTexturesFunctor::CheckIsInside
bool CheckIsInside(const int &signX, const int &signY, const OffsetType &currentOff, const OffsetType &stopOffset)
Definition: otbSFSTexturesFunctor.h:365
otb::Functor::SFSTexturesFunctor::SetSpectralThreshold
void SetSpectralThreshold(InternalPixelType thresh)
Definition: otbSFSTexturesFunctor.h:88
otb::Functor::SFSTexturesFunctor::m_SpatialThreshold
unsigned int m_SpatialThreshold
Definition: otbSFSTexturesFunctor.h:405
otb::Functor::SFSTexturesFunctor::SetDirectionStep
void SetDirectionStep(double step)
Definition: otbSFSTexturesFunctor.h:105
otb::Functor::SFSTexturesFunctor::m_SpectralThreshold
InternalPixelType m_SpectralThreshold
Definition: otbSFSTexturesFunctor.h:408
otb::Functor::SFSTexturesFunctor::OutputValueType
TOutputValue OutputValueType
Definition: otbSFSTexturesFunctor.h:81
otb::Functor::SFSTexturesFunctor::ComputeSDi
double ComputeSDi(const TIter &it, const OffsetType &stopOffset)
Definition: otbSFSTexturesFunctor.h:309
otb::Functor::SFSTexturesFunctor::SetAlpha
void SetAlpha(double alpha)
Definition: otbSFSTexturesFunctor.h:96