OTB  6.7.0
Orfeo Toolbox
otbLandsatTMIndices.h
Go to the documentation of this file.
1 /*
2  * Copyright (C) 2005-2017 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 otbLandsatTMIndices_h
22 #define otbLandsatTMIndices_h
23 
24 #include "otbMath.h"
26 #include "otbBandName.h"
27 #include "otbFuzzyVariable.h"
28 #include <vector>
29 #include <algorithm>
30 
31 namespace otb
32 {
33 
34 namespace Functor
35 {
36 
37 namespace LandsatTM
38 {
39 
41 enum SATType {L5, L7};
65 template<class TInput, class TOutput>
67 {
68 public:
69 
70  //operators !=
71  bool operator !=(const LandsatTMIndexBase&) const
72  {
73  return false;
74  }
75  //operator ==
76  bool operator ==(const LandsatTMIndexBase& other) const
77  {
78  return !(*this != other);
79  }
80 
84  virtual ~LandsatTMIndexBase() {}
85 
87  void SetIndex(BandName::LandsatTMBandNames band, unsigned int channel)
88  {
89  switch( band )
90  {
91  case BandName::TM1:
92  m_TM1 = channel;
93  break;
94  case BandName::TM2:
95  m_TM2 = channel;
96  break;
97  case BandName::TM3:
98  m_TM3 = channel;
99  break;
100  case BandName::TM4:
101  m_TM4 = channel;
102  break;
103  case BandName::TM5:
104  m_TM5 = channel;
105  break;
106  case BandName::TM60:
107  m_TM60 = channel;
108  break;
109  case BandName::TM61:
110  m_TM61 = channel;
111  break;
112  case BandName::TM62:
113  m_TM62 = channel;
114  break;
115  case BandName::TM7:
116  m_TM7 = channel;
117  break;
118  }
119  }
120 
122  unsigned int GetIndex(BandName::LandsatTMBandNames band) const
123  {
124  switch( band )
125  {
126  case BandName::TM1:
127  return m_TM1;
128  break;
129  case BandName::TM2:
130  return m_TM2;
131  break;
132  case BandName::TM3:
133  return m_TM3;
134  break;
135  case BandName::TM4:
136  return m_TM4;
137  break;
138  case BandName::TM5:
139  return m_TM5;
140  break;
141  case BandName::TM60:
142  return m_TM60;
143  break;
144  case BandName::TM61:
145  return m_TM61;
146  break;
147  case BandName::TM62:
148  return m_TM62;
149  break;
150  case BandName::TM7:
151  return m_TM7;
152  break;
153  }
154  return m_TM1;
155  }
157 
158  unsigned int GetTM1() const
159  {
160  return this->m_TM1;
161  }
162 
163  unsigned int GetTM2() const
164  {
165  return this->m_TM2;
166  }
167 
168  unsigned int GetTM3() const
169  {
170  return this->m_TM3;
171  }
172 
173  unsigned int GetTM4() const
174  {
175  return this->m_TM4;
176  }
177 
178  unsigned int GetTM5() const
179  {
180  return this->m_TM5;
181  }
182 
183  unsigned int GetTM60() const
184  {
185  return this->m_TM60;
186  }
187 
188  unsigned int GetTM61() const
189  {
190  return this->m_TM61;
191  }
192 
193  unsigned int GetTM62() const
194  {
195  return this->m_TM62;
196  }
197 
198  unsigned int GetTM7() const
199  {
200  return this->m_TM7;
201  }
202 
203  void SetSAT(SATType sat)
204  {
205 
206  this->m_SAT = sat;
207 
208  if( sat == L5 )
209  {
210  m_TM1 = 0;
211  m_TM2 = 1;
212  m_TM3 = 2;
213  m_TM4 = 3;
214  m_TM5 = 4;
215  m_TM60 = 5;
216  m_TM7 = 6;
217  m_SAT = L5;
218  }
219  else
220  {
221  m_TM1 = 0;
222  m_TM2 = 1;
223  m_TM3 = 2;
224  m_TM4 = 3;
225  m_TM5 = 4;
226  m_TM61 = 5;
227  m_TM62 = 6;
228  m_TM7 = 7;
229  m_SAT = L7;
230  }
231  }
232 
233  SATType GetSAT() const
234  {
235  return this->m_SAT;
236  }
237 
239  {
240  this->m_Degree = deg;
241  }
242 
244  {
245  return this->m_Degree;
246  }
247 
249  {
250  this->m_Reflectance = ref;
251  }
252 
254  {
255  return this->m_Reflectance;
256  }
257 
259  {
260  return this->m_EpsilonToBeConsideredAsZero;
261  }
262 
263 protected:
265 
266  TInput PrepareValues(const TInput& inputPixel)
267  {
268 
269  TInput newPix( inputPixel );
270  if( this->m_Degree == Kelvin )
271  {
272  if(this->m_SAT == L5)
273  {
274  newPix[this->m_TM60] = newPix[this->m_TM60]-273.15;
275  }
276  else
277  {
278  newPix[this->m_TM61] = newPix[this->m_TM61]-273.15;
279  newPix[this->m_TM62] = newPix[this->m_TM62]-273.15;
280  }
281  }
282  else
283  if( this->m_Degree == HundredsKelvin )
284  {
285  if(this->m_SAT == L5)
286  {
287  newPix[this->m_TM60] = newPix[this->m_TM60]/100.-273.15;
288  }
289  else
290  {
291  newPix[this->m_TM61] = newPix[this->m_TM61]/100.-273.15;
292  newPix[this->m_TM62] = newPix[this->m_TM62]/100.-273.15;
293  }
294  }
295  if( this->m_Reflectance == Thousands )
296  {
297  newPix[this->m_TM1] = newPix[this->m_TM1]/1000.;
298  newPix[this->m_TM2] = newPix[this->m_TM2]/1000.;
299  newPix[this->m_TM3] = newPix[this->m_TM3]/1000.;
300  newPix[this->m_TM4] = newPix[this->m_TM4]/1000.;
301  newPix[this->m_TM5] = newPix[this->m_TM5]/1000.;
302  newPix[this->m_TM7] = newPix[this->m_TM7]/1000.;
303  }
304 
305  return newPix;
306  }
307 
308 
310 
311  unsigned int m_TM1;
312  unsigned int m_TM2;
313  unsigned int m_TM3;
314  unsigned int m_TM4;
315  unsigned int m_TM5;
316  unsigned int m_TM60;
317  unsigned int m_TM61;
318  unsigned int m_TM62;
319  unsigned int m_TM7;
320 
324 
325 };
326 
327 
336 template <class TInput, class TOutput>
337 class LandsatTMIndex : public LandsatTMIndexBase<TInput, TOutput>
338 {
339 public:
340 
341  // Operator performing the work
342  inline TOutput operator ()(const TInput& inputPixel) const;
343 
345  virtual std::string GetName() const = 0;
346 
348  ~LandsatTMIndex() override {}
349 
350 
351 };
352 
370 template <class TInput, class TOutput>
371 class Bright : public LandsatTMIndex<TInput, TOutput>
372 {
373 public:
375  std::string GetName() const override
376  {
377  return "Bright";
378  }
379 
380  Bright() {}
381  ~Bright() override {}
382 
383  inline TOutput operator ()(const TInput& inputPixel)
384  {
385 
386  TInput newPixel(this->PrepareValues( inputPixel ));
387  double result = (newPixel[this->m_TM1]+newPixel[this->m_TM2]+2*newPixel[this->m_TM3]+2*newPixel[this->m_TM4]+newPixel[this->m_TM5]+newPixel[this->m_TM7])/8.0;
388  return static_cast<TOutput>(result);
389  }
390 
391 
392 };
393 
411 template <class TInput, class TOutput>
412 class Vis : public LandsatTMIndex<TInput, TOutput>
413 {
414 public:
416  std::string GetName() const override
417  {
418  return "Vis";
419  }
420 
421  Vis() {}
422  ~Vis() override {}
423 
424  inline TOutput operator ()(const TInput& inputPixel)
425  {
426  TInput newPixel(this->PrepareValues( inputPixel ));
427  double result = (newPixel[this->m_TM1]+newPixel[this->m_TM2]+newPixel[this->m_TM3])/3.0;
428  return static_cast<TOutput>(result);
429  }
430 
431 
432 };
433 
434 
444 template <class TInput, class TOutput>
445 class NIR : public LandsatTMIndex<TInput, TOutput>
446 {
447 public:
449  std::string GetName() const override
450  {
451  return "NIR";
452  }
453 
454  NIR() {}
455  ~NIR() override {}
456 
457  inline TOutput operator ()(const TInput& inputPixel)
458  {
459  TInput newPixel(this->PrepareValues( inputPixel ));
460  double result = newPixel[this->m_TM4];
461  return static_cast<TOutput>(result);
462  }
463 
464 
465 };
466 
476 template <class TInput, class TOutput>
477 class MIR1 : public LandsatTMIndex<TInput, TOutput>
478 {
479 public:
481  std::string GetName() const override
482  {
483  return "MIR1";
484  }
485 
486  MIR1() {}
487  ~MIR1() override {}
488 
489  inline TOutput operator ()(const TInput& inputPixel)
490  {
491  TInput newPixel(this->PrepareValues( inputPixel ));
492  double result = newPixel[this->m_TM5];
493  return static_cast<TOutput>(result);
494  }
495 
496 
497 };
498 
508 template <class TInput, class TOutput>
509 class MIR2 : public LandsatTMIndex<TInput, TOutput>
510 {
511 public:
513  std::string GetName() const override
514  {
515  return "MIR2";
516  }
517 
518  MIR2() {}
519  ~MIR2() override {}
520 
521  inline TOutput operator ()(const TInput& inputPixel)
522  {
523  TInput newPixel(this->PrepareValues( inputPixel ));
524  double result = newPixel[this->m_TM7];
525  return static_cast<TOutput>(result);
526  }
527 
528 
529 };
530 
540 template <class TInput, class TOutput>
541 class TIR : public LandsatTMIndex<TInput, TOutput>
542 {
543 public:
545  std::string GetName() const override
546  {
547  return "TIR";
548  }
549 
550  TIR() {}
551  ~TIR() override {}
552 
553  inline TOutput operator ()(const TInput& inputPixel)
554  {
555  TInput newPixel(this->PrepareValues( inputPixel ));
556  double result = newPixel[this->m_TM62];
557 
558  if( this->m_SAT == L5 )
559  result = newPixel[this->m_TM60];
560 
561  return static_cast<TOutput>(result);
562  }
563 
564 
565 };
566 
586 template <class TInput, class TOutput>
587 class MIRTIR : public LandsatTMIndex<TInput, TOutput>
588 {
589 public:
591  std::string GetName() const override
592  {
593  return "MIRTIR";
594  }
595 
596  MIRTIR() {}
597  ~MIRTIR() override {}
598 
599  inline TOutput operator ()(const TInput& inputPixel)
600  {
601  TInput newPixel(this->PrepareValues( inputPixel ));
602  double tir = newPixel[this->m_TM62];
603  double mir1 = newPixel[this->m_TM5];
604 
605  if( this->m_SAT == L5 )
606  tir = newPixel[this->m_TM60];
607 
608  double result = 255*(1 - mir1)*(tir+100)/100.;
609 
610  return static_cast<TOutput>(result);
611  }
612 
613 
614 };
615 
626 template <class TInput, class TOutput>
627 class NDVI : public LandsatTMIndex<TInput, TOutput>
628 {
629 public:
631  std::string GetName() const override
632  {
633  return "NDVI";
634  }
635 
636  NDVI() {}
637  ~NDVI() override {}
638 
639  inline TOutput operator ()(const TInput& inputPixel)
640  {
641  TInput newPixel(this->PrepareValues( inputPixel ));
642  double result = (newPixel[this->m_TM4] - newPixel[this->m_TM3])/
643  (newPixel[this->m_TM4] + newPixel[this->m_TM3] + this->m_EpsilonToBeConsideredAsZero);
644 
645  return static_cast<TOutput>(result);
646  }
647 
648 
649 };
650 
673 template <class TInput, class TOutput>
674 class NDBSI : public LandsatTMIndex<TInput, TOutput>
675 {
676 public:
678  std::string GetName() const override
679  {
680  return "NDBSI";
681  }
682 
683  NDBSI() {}
684  ~NDBSI() override {}
685 
686  inline TOutput operator ()(const TInput& inputPixel)
687  {
688  TInput newPixel(this->PrepareValues( inputPixel ));
689  double result = (newPixel[this->m_TM5] - newPixel[this->m_TM4])/
690  (newPixel[this->m_TM5] + newPixel[this->m_TM4] + this->m_EpsilonToBeConsideredAsZero);
691 
692  return static_cast<TOutput>(result);
693  }
694 
695 
696 };
697 
716 template <class TInput, class TOutput>
717 class BIO : public LandsatTMIndex<TInput, TOutput>
718 {
719 public:
721  std::string GetName() const override
722  {
723  return "BIO";
724  }
725 
726  BIO() {}
727  ~BIO() override {}
728 
729  inline TOutput operator ()(const TInput& inputPixel)
730  {
731  TInput newPixel(this->PrepareValues( inputPixel ));
732  double result = ((newPixel[this->m_TM5] + newPixel[this->m_TM3])
733  - (newPixel[this->m_TM4] + newPixel[this->m_TM1]))
734  /((newPixel[this->m_TM5] + newPixel[this->m_TM3])
735  + (newPixel[this->m_TM4] + newPixel[this->m_TM1]));
736 
737  return static_cast<TOutput>(result);
738  }
739 
740 
741 };
742 
743 
762 template <class TInput, class TOutput>
763 class NDSI : public LandsatTMIndex<TInput, TOutput>
764 {
765 public:
767  std::string GetName() const override
768  {
769  return "NDSI";
770  }
771 
772  NDSI() {}
773  ~NDSI() override {}
774 
775  inline TOutput operator ()(const TInput& inputPixel)
776  {
777  TInput newPixel(this->PrepareValues( inputPixel ));
778  double result = (newPixel[this->m_TM2] - newPixel[this->m_TM5])
779  /(newPixel[this->m_TM2] + newPixel[this->m_TM5] + this->m_EpsilonToBeConsideredAsZero);
780 
781  return static_cast<TOutput>(result);
782  }
783 
784 
785 };
786 
819 template <class TInput, class TOutput>
820 class NDSIVis : public LandsatTMIndex<TInput, TOutput>
821 {
822 public:
824  std::string GetName() const override
825  {
826  return "NDSIVis";
827  }
828 
829  NDSIVis() {}
830  ~NDSIVis() override {}
831 
832  inline TOutput operator ()(const TInput& inputPixel)
833  {
834  TInput newPixel(this->PrepareValues( inputPixel ));
835  double vis = (newPixel[this->m_TM1]+newPixel[this->m_TM2]+newPixel[this->m_TM3])/3.0;
836  double result = (vis - newPixel[this->m_TM5])/(vis + newPixel[this->m_TM5] + this->m_EpsilonToBeConsideredAsZero);
837 
838  return static_cast<TOutput>(result);
839  }
840 
841 
842 };
843 
864 template <class TInput, class TOutput>
865 class NDBBBI : public LandsatTMIndex<TInput, TOutput>
866 {
867 public:
869  std::string GetName() const override
870  {
871  return "NDBBBI";
872  }
873 
874  NDBBBI() {}
875  ~NDBBBI() override {}
876 
877  inline TOutput operator ()(const TInput& inputPixel)
878  {
879  TInput newPixel(this->PrepareValues( inputPixel ));
880  double result = (newPixel[this->m_TM1] - newPixel[this->m_TM5])
881  /(newPixel[this->m_TM1] + newPixel[this->m_TM5] + this->m_EpsilonToBeConsideredAsZero);
882  return static_cast<TOutput>(result);
883  }
884 
885 
886 };
887 
909 template <class TInput>
910 class LinguisticVariables : public LandsatTMIndexBase<TInput, itk::FixedArray<unsigned int, 11> >
911 {
912 public:
913 
914  typedef typename TInput::ValueType PrecisionType;
917 
920 
922  virtual std::string GetName() const
923  {
924  return "LandsatTM Linguistic Variables";
925  }
926 
928  {
940 
943 
944  // the thresholds are computed wrt Baraldi's paper (normalized 0-255 values)
945  m_FvBright->SetMembership(Low, minimumValue, minimumValue, 40/255., 40/255.);
946  m_FvBright->SetMembership(Medium, 40/255., 40/255., 60/255., 60/255.);
947  m_FvBright->SetMembership(High, 60/255., 60/255., maximumValue, maximumValue);
948 
949  m_FvVis->SetMembership(Low, minimumValue, minimumValue, 30/255., 30/255.);
950  m_FvVis->SetMembership(Medium, 30/255., 30/255., 50/255., 50/255.);
951  m_FvVis->SetMembership(High, 50/255., 50/255., maximumValue, maximumValue);
952 
953  m_FvNIR->SetMembership(Low, minimumValue, minimumValue, 40/255., 40/255.);
954  m_FvNIR->SetMembership(Medium, 40/255., 40/255., 60/255., 60/255.);
955  m_FvNIR->SetMembership(High, 60/255., 60/255., maximumValue, maximumValue);
956 
957  m_FvMIR1->SetMembership(Low, minimumValue, minimumValue, 40/255., 40/255.);
958  m_FvMIR1->SetMembership(Medium, 40/255., 40/255., 60/255., 60/255.);
959  m_FvMIR1->SetMembership(High, 60/255., 60/255., maximumValue, maximumValue);
960 
961  m_FvMIR2->SetMembership(Low, minimumValue, minimumValue, 30/255., 30/255.);
962  m_FvMIR2->SetMembership(Medium, 30/255., 30/255., 50/255., 50/255.);
963  m_FvMIR2->SetMembership(High, 50/255., 50/255., maximumValue, maximumValue);
964 
965  m_FvTIR->SetMembership(Low, minimumValue, minimumValue, 0, 0);
966  m_FvTIR->SetMembership(Medium, 0, 0, 28, 28);
967  m_FvTIR->SetMembership(High, 28, 28, maximumValue, maximumValue);
968 
969  m_FvMIRTIR->SetMembership(Low, minimumValue, minimumValue, 180, 180);
970  m_FvMIRTIR->SetMembership(Medium, 180, 180, 220, 220);
971  m_FvMIRTIR->SetMembership(High, 220, 220, maximumValue, maximumValue);
972 
973  m_FvNDSIVis->SetMembership(Low, minimumValue, minimumValue, 0, 0);
974  m_FvNDSIVis->SetMembership(Medium, 0, 0, 0.5, 0.5);
975  m_FvNDSIVis->SetMembership(High, 0.5, 0.5, maximumValue, maximumValue);
976 
977  m_FvNDBBBI ->SetMembership(Low, minimumValue, minimumValue, -0.20, -0.20);
978  m_FvNDBBBI->SetMembership(Medium, -0.20, -0.20, 0.10, 0.10);
979  m_FvNDBBBI->SetMembership(High, 0.10, 0.10, maximumValue, maximumValue);
980 
981  m_FvNDVI->SetMembership(Low, minimumValue, minimumValue, 0.35, 0.35);
982  m_FvNDVI->SetMembership(Medium, 0.35, 0.35, 0.6, 0.6);
983  m_FvNDVI->SetMembership(High, 0.6, 0.6, maximumValue, maximumValue);
984 
985  m_FvNDBSI->SetMembership(Low, minimumValue, minimumValue, -0.20, -0.20);
986  m_FvNDBSI->SetMembership(Medium, -0.20, -0.20, 0.0, 0.0);
987  m_FvNDBSI->SetMembership(High, 0.0, 0.0, maximumValue, maximumValue);
988  }
989  ~LinguisticVariables() override {}
990 
991  inline itk::FixedArray<unsigned int, 11> operator ()(const TInput& inputPixel)
992  {
993  TInput newPixel(this->PrepareValues( inputPixel ));
995 
996  result[ bright ] = m_FvBright->GetMaxVar(Bright<TInput, PrecisionType>()( newPixel ));
997  result[ vis ] = m_FvVis->GetMaxVar(Vis<TInput, PrecisionType>()( newPixel ));
998  result[ nir ] = m_FvNIR->GetMaxVar(NIR<TInput, PrecisionType>()( newPixel ));
999  result[ mir1 ] = m_FvMIR1->GetMaxVar(MIR1<TInput, PrecisionType>()( newPixel ));
1000 
1002  mir2F.SetSAT( this->m_SAT );
1003  result[ mir2 ] = m_FvMIR2->GetMaxVar(mir2F( newPixel ));
1004 
1006  tirF.SetSAT( this->m_SAT );
1007  result[ tir ] = m_FvTIR->GetMaxVar(tirF( newPixel ));
1008 
1010  mirtirF.SetSAT( this->m_SAT );
1011  result[ mirtir ] = m_FvMIRTIR->GetMaxVar(mirtirF( newPixel ));
1012 
1013  result[ ndsivis ] = m_FvNDSIVis->GetMaxVar(NDSIVis<TInput, PrecisionType>()( newPixel ));
1014 
1015  result[ ndbbbi ] = m_FvNDBBBI->GetMaxVar(NDBBBI<TInput, PrecisionType>()( newPixel ));
1016 
1017  result[ ndvi ] = m_FvNDVI->GetMaxVar(NDVI<TInput, PrecisionType>()( newPixel ));
1018 
1019  result[ ndbsi ] = m_FvNDBSI->GetMaxVar(NDBSI<TInput, PrecisionType>()( newPixel ));
1020 
1021  return result;
1022  }
1023 
1024 protected:
1036 
1037 
1038 };
1039 
1040 
1058 template <class TInput, class TOutput>
1059 class KernelSpectralRule : public LandsatTMIndexBase<TInput, TOutput >
1060 {
1061 public:
1062 
1063  typedef typename TInput::ValueType PrecisionType;
1064  typedef bool OutputPixelType;
1065 
1067  virtual std::string GetName() const
1068  {
1069  return "LandsatTM KernelSpectralRule";
1070  }
1071 
1072  KernelSpectralRule() : m_TV1(0.7), m_TV2(0.5) { }
1073  ~KernelSpectralRule() override {}
1074 
1076  {
1077  this->m_TV1 = tv1;
1078  }
1079 
1081  {
1082  this->m_TV2 = tv2;
1083  }
1084 
1086  {
1087  return this->m_TV1;
1088  }
1089 
1091  {
1092  return this->m_TV2;
1093  }
1094 
1095 protected:
1098 
1101 
1102  void SetMinMax(const TInput& inputPixel, PrecisionType* max13, PrecisionType* min123, PrecisionType* max123, PrecisionType* min12347, PrecisionType* max12347, PrecisionType* max234, PrecisionType* max45)
1103  {
1104  std::vector< PrecisionType > v13;
1105  v13.push_back(inputPixel[this->m_TM1]);
1106  v13.push_back(inputPixel[this->m_TM3]);
1107 
1108  *max13 = *(std::max_element ( v13.begin(), v13.end() ));
1109 
1110 
1111  std::vector< PrecisionType > v123;
1112  v123.push_back(inputPixel[this->m_TM1]);
1113  v123.push_back(inputPixel[this->m_TM2]);
1114  v123.push_back(inputPixel[this->m_TM3]);
1115 
1116  *max123 = *(std::max_element ( v123.begin(), v123.end() ));
1117  *min123 = *(std::min_element ( v123.begin(), v123.end() ));
1118 
1119 
1120  std::vector< PrecisionType > v12347;
1121  v12347.push_back(inputPixel[this->m_TM1]);
1122  v12347.push_back(inputPixel[this->m_TM2]);
1123  v12347.push_back(inputPixel[this->m_TM3]);
1124  v12347.push_back(inputPixel[this->m_TM4]);
1125  v12347.push_back(inputPixel[this->m_TM7]);
1126 
1127  *max12347 = *(std::max_element ( v12347.begin(), v12347.end() ));
1128  *min12347 = *(std::min_element ( v12347.begin(), v12347.end() ));
1129 
1130  std::vector< PrecisionType > v234;
1131  v234.push_back(inputPixel[this->m_TM2]);
1132  v234.push_back(inputPixel[this->m_TM3]);
1133  v234.push_back(inputPixel[this->m_TM4]);
1134 
1135  *max234 = *(std::max_element ( v234.begin(), v234.end() ));
1136 
1137  std::vector< PrecisionType > v45;
1138  v45.push_back(inputPixel[this->m_TM4]);
1139  v45.push_back(inputPixel[this->m_TM5]);
1140 
1141  *max45 = *(std::max_element ( v45.begin(), v45.end() ));
1142  }
1143 };
1144 
1160 template <class TInput, class TOutput>
1161 class ThickCloudsSpectralRule : public KernelSpectralRule<TInput, TOutput>
1162 {
1163 public:
1164 
1165  typedef typename TInput::ValueType PrecisionType;
1166  typedef bool OutputPixelType;
1167 
1169  std::string GetName() const override
1170  {
1171  return "LandsatTM ThickCloudsSpectralRule";
1172  }
1173 
1176 
1177  inline TOutput operator ()(const TInput& inputPixel)
1178  {
1179  TInput newPixel(this->PrepareValues( inputPixel ));
1180  PrecisionType max13;
1181  PrecisionType max123;
1182  PrecisionType min123;
1183  PrecisionType max12347;
1184  PrecisionType min12347;
1185  PrecisionType max234;
1186  PrecisionType max45;
1187  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1188 
1189  bool result = (
1190  ((min123 >= (this->m_TV1 * max123)
1191  && (max123 <= this->m_TV1 * newPixel[this->m_TM4]))
1192  || ((newPixel[this->m_TM2] >= this->m_TV1 * max13)
1193  && (max123 <= newPixel[this->m_TM4])))
1194  && (newPixel[this->m_TM5] <= this->m_TV1 * newPixel[this->m_TM4])
1195  && (newPixel[this->m_TM5] >= this->m_TV1 * max123)
1196  && (newPixel[this->m_TM7] <= this->m_TV1 * newPixel[this->m_TM4]));
1197 
1198  return static_cast<TOutput>(result);
1199  }
1200 
1201 };
1202 
1218 template <class TInput, class TOutput>
1219 class ThinCloudsSpectralRule : public KernelSpectralRule<TInput, TOutput>
1220 {
1221 public:
1222 
1223  typedef typename TInput::ValueType PrecisionType;
1224  typedef bool OutputPixelType;
1225 
1227  std::string GetName() const override
1228  {
1229  return "LandsatTM ThinCloudsSpectralRule";
1230  }
1231 
1234 
1235  inline TOutput operator ()(const TInput& inputPixel)
1236  {
1237  TInput newPixel(this->PrepareValues( inputPixel ));
1238  PrecisionType max13;
1239  PrecisionType max123;
1240  PrecisionType min123;
1241  PrecisionType max12347;
1242  PrecisionType min12347;
1243  PrecisionType max234;
1244  PrecisionType max45;
1245  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1246 
1247 
1248  bool result = (min123 >= (this->m_TV1 * max123))
1249  && (newPixel[this->m_TM4] >= max123)
1250  && !((newPixel[this->m_TM1] <= newPixel[this->m_TM2]
1251  && newPixel[this->m_TM2] <= newPixel[this->m_TM3]
1252  && newPixel[this->m_TM3] <= newPixel[this->m_TM4])
1253  && (newPixel[this->m_TM3] >= this->m_TV1 * newPixel[this->m_TM4]))
1254  && (newPixel[this->m_TM4] >= this->m_TV1 * newPixel[this->m_TM5])
1255  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM4])
1256  && (newPixel[this->m_TM5] >= this->m_TV1 * max123)
1257  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM7]);
1258 
1259  return static_cast<TOutput>(result);
1260  }
1261 
1262 };
1263 
1279 template <class TInput, class TOutput>
1280 class SnowOrIceSpectralRule : public KernelSpectralRule<TInput, TOutput>
1281 {
1282 public:
1283 
1284  typedef typename TInput::ValueType PrecisionType;
1285  typedef bool OutputPixelType;
1286 
1288  std::string GetName() const override
1289  {
1290  return "LandsatTM SnowOrIceSpectralRule";
1291  }
1292 
1295 
1296  inline TOutput operator ()(const TInput& inputPixel)
1297  {
1298  TInput newPixel(this->PrepareValues( inputPixel ));
1299  PrecisionType max13;
1300  PrecisionType max123;
1301  PrecisionType min123;
1302  PrecisionType max12347;
1303  PrecisionType min12347;
1304  PrecisionType max234;
1305  PrecisionType max45;
1306  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1307 
1308 
1309  bool result = (min123 >= (this->m_TV1 * max123))
1310  && (newPixel[this->m_TM4] >= (this->m_TV1 * max123))
1311  && (newPixel[this->m_TM5] <= this->m_TV2 * newPixel[this->m_TM4])
1312  && (newPixel[this->m_TM5] <= this->m_TV1 * min123)
1313  && (newPixel[this->m_TM5] <= this->m_TV1 * max123)
1314  && (newPixel[this->m_TM7] <= this->m_TV2 * newPixel[this->m_TM4])
1315  && (newPixel[this->m_TM7] <= this->m_TV1 * min123);
1316 
1317  return static_cast<TOutput>(result);
1318  }
1319 
1320 };
1321 
1322 
1338 template <class TInput, class TOutput>
1339 class WaterOrShadowSpectralRule : public KernelSpectralRule<TInput, TOutput>
1340 {
1341 public:
1342 
1343  typedef typename TInput::ValueType PrecisionType;
1344  typedef bool OutputPixelType;
1345 
1347  std::string GetName() const override
1348  {
1349  return "LandsatTM WaterOrShadowSpectralRule";
1350  }
1351 
1354 
1355  inline TOutput operator ()(const TInput& inputPixel)
1356  {
1357  TInput newPixel(this->PrepareValues( inputPixel ));
1358  bool result = (newPixel[this->m_TM1] >= newPixel[this->m_TM2])
1359  && (newPixel[this->m_TM2] >= newPixel[this->m_TM3])
1360  && (newPixel[this->m_TM3] >= newPixel[this->m_TM4])
1361  && (newPixel[this->m_TM4] >= newPixel[this->m_TM5])
1362  && (newPixel[this->m_TM4] >= newPixel[this->m_TM7]);
1363 
1364  return static_cast<TOutput>(result);
1365  }
1366 
1367 };
1368 
1369 
1385 template <class TInput, class TOutput>
1387 {
1388 public:
1389 
1390  typedef typename TInput::ValueType PrecisionType;
1391  typedef bool OutputPixelType;
1392 
1394  std::string GetName() const override
1395  {
1396  return "LandsatTM PitbogOrGreenhouseSpectralRule";
1397  }
1398 
1401 
1402  inline TOutput operator ()(const TInput& inputPixel)
1403  {
1404  TInput newPixel(this->PrepareValues( inputPixel ));
1405  PrecisionType max13;
1406  PrecisionType max123;
1407  PrecisionType min123;
1408  PrecisionType max12347;
1409  PrecisionType min12347;
1410  PrecisionType max234;
1411  PrecisionType max45;
1412  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1413 
1414 
1415  bool result = (newPixel[this->m_TM3] >= this->m_TV1 * newPixel[this->m_TM1])
1416  && (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM3])
1417  && (max123 <= this->m_TV1 * newPixel[this->m_TM4])
1418  && (newPixel[this->m_TM5] <= this->m_TV1 * newPixel[this->m_TM4])
1419  && (newPixel[this->m_TM3] >= this->m_TV2 * newPixel[this->m_TM5])
1420  && (min123 >= this->m_TV1 * newPixel[this->m_TM7]);
1421 
1422 
1423  return static_cast<TOutput>(result);
1424  }
1425 
1426 };
1427 
1443 template <class TInput, class TOutput>
1444 class DominantBlueSpectralRule : public KernelSpectralRule<TInput, TOutput>
1445 {
1446 public:
1447 
1448  typedef typename TInput::ValueType PrecisionType;
1449  typedef bool OutputPixelType;
1450 
1452  std::string GetName() const override
1453  {
1454  return "LandsatTM DominantBlueSpectralRule";
1455  }
1456 
1459 
1460  inline TOutput operator ()(const TInput& inputPixel)
1461  {
1462  TInput newPixel(this->PrepareValues( inputPixel ));
1463  bool result = (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM2])
1464  && (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM3])
1465  && (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM4])
1466  && (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM5])
1467  && (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM7]);
1468 
1469 
1470  return static_cast<TOutput>(result);
1471  }
1472 
1473 };
1474 
1475 
1491 template <class TInput, class TOutput>
1492 class VegetationSpectralRule : public KernelSpectralRule<TInput, TOutput>
1493 {
1494 public:
1495 
1496  typedef typename TInput::ValueType PrecisionType;
1497  typedef bool OutputPixelType;
1498 
1500  std::string GetName() const override
1501  {
1502  return "LandsatTM VegetationSpectralRule";
1503  }
1504 
1507 
1508  inline TOutput operator ()(const TInput& inputPixel)
1509  {
1510  TInput newPixel(this->PrepareValues( inputPixel ));
1511  PrecisionType max13;
1512  PrecisionType max123;
1513  PrecisionType min123;
1514  PrecisionType max12347;
1515  PrecisionType min12347;
1516  PrecisionType max234;
1517  PrecisionType max45;
1518  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1519 
1520 
1521  bool result = (newPixel[this->m_TM2] >= this->m_TV2 * newPixel[this->m_TM1])
1522  && (newPixel[this->m_TM2] >= this->m_TV1 * newPixel[this->m_TM3])
1523  && (newPixel[this->m_TM3] < this->m_TV1 * newPixel[this->m_TM4])
1524  && (newPixel[this->m_TM4] > max123)
1525  && (newPixel[this->m_TM5] < this->m_TV1 * newPixel[this->m_TM4])
1526  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM3])
1527  && (newPixel[this->m_TM7] < this->m_TV1 * newPixel[this->m_TM5]);
1528 
1529 
1530  return static_cast<TOutput>(result);
1531  }
1532 
1533 };
1534 
1535 
1551 template <class TInput, class TOutput>
1552 class RangelandSpectralRule : public KernelSpectralRule<TInput, TOutput>
1553 {
1554 public:
1555 
1556  typedef typename TInput::ValueType PrecisionType;
1557  typedef bool OutputPixelType;
1558 
1560  std::string GetName() const override
1561  {
1562  return "LandsatTM RangelandSpectralRule";
1563  }
1564 
1567 
1568  inline TOutput operator ()(const TInput& inputPixel)
1569  {
1570  TInput newPixel(this->PrepareValues( inputPixel ));
1571  PrecisionType max13;
1572  PrecisionType max123;
1573  PrecisionType min123;
1574  PrecisionType max12347;
1575  PrecisionType min12347;
1576  PrecisionType max234;
1577  PrecisionType max45;
1578  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1579 
1580 
1581  bool result = (newPixel[this->m_TM2] >= this->m_TV2 * newPixel[this->m_TM1])
1582  && (newPixel[this->m_TM2] >= this->m_TV1 * newPixel[this->m_TM3])
1583  && (newPixel[this->m_TM4] > max123)
1584  && (newPixel[this->m_TM3] < this->m_TV1 * newPixel[this->m_TM4])
1585  && (newPixel[this->m_TM4] >= this->m_TV1 * newPixel[this->m_TM5])
1586  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM4])
1587  && (newPixel[this->m_TM5] > max123)
1588  && (newPixel[this->m_TM7] < this->m_TV1 * max45)
1589  && (newPixel[this->m_TM5] >= newPixel[this->m_TM7]);
1590 
1591 
1592  return static_cast<TOutput>(result);
1593  }
1594 
1595 };
1596 
1612 template <class TInput, class TOutput>
1614 {
1615 public:
1616 
1617  typedef typename TInput::ValueType PrecisionType;
1618  typedef bool OutputPixelType;
1619 
1621  std::string GetName() const override
1622  {
1623  return "LandsatTM BarrenLandOrBuiltUpOrCloudsSpectralRule";
1624  }
1625 
1628 
1629  inline TOutput operator ()(const TInput& inputPixel)
1630  {
1631  TInput newPixel(this->PrepareValues( inputPixel ));
1632  PrecisionType max13;
1633  PrecisionType max123;
1634  PrecisionType min123;
1635  PrecisionType max12347;
1636  PrecisionType min12347;
1637  PrecisionType max234;
1638  PrecisionType max45;
1639  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1640 
1641 
1642  bool result = (newPixel[this->m_TM3] >= this->m_TV2 * newPixel[this->m_TM1])
1643  && (newPixel[this->m_TM3] >= this->m_TV1 * newPixel[this->m_TM2])
1644  && (newPixel[this->m_TM4] >= this->m_TV1 * max123)
1645  && (newPixel[this->m_TM5] >= max123)
1646  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM4])
1647  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM7])
1648  && (newPixel[this->m_TM7] >= this->m_TV2 * max45);
1649 
1650  return static_cast<TOutput>(result);
1651  }
1652 
1653 };
1654 
1670 template <class TInput, class TOutput>
1672 {
1673 public:
1674 
1675  typedef typename TInput::ValueType PrecisionType;
1676  typedef bool OutputPixelType;
1677 
1679  std::string GetName() const override
1680  {
1681  return "LandsatTM FlatResponseBarrenLandOrBuiltUpSpectralRule";
1682  }
1683 
1686 
1687  inline TOutput operator ()(const TInput& inputPixel)
1688  {
1689  TInput newPixel(this->PrepareValues( inputPixel ));
1690  PrecisionType max13;
1691  PrecisionType max123;
1692  PrecisionType min123;
1693  PrecisionType max12347;
1694  PrecisionType min12347;
1695  PrecisionType max234;
1696  PrecisionType max45;
1697  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1698 
1699 
1700  bool result = (newPixel[this->m_TM5] >= this->m_TV1 * max12347)
1701  && (min12347 >= this->m_TV2 * newPixel[this->m_TM5]);
1702 
1703  return static_cast<TOutput>(result);
1704  }
1705 
1706 };
1707 
1708 
1724 template <class TInput, class TOutput>
1726 {
1727 public:
1728 
1729  typedef typename TInput::ValueType PrecisionType;
1730  typedef bool OutputPixelType;
1731 
1733  std::string GetName() const override
1734  {
1735  return "LandsatTM ShadowWithBarrenLandSpectralRule";
1736  }
1737 
1740 
1741  inline TOutput operator ()(const TInput& inputPixel)
1742  {
1743  TInput newPixel(this->PrepareValues( inputPixel ));
1744  bool result = (newPixel[this->m_TM1] >= newPixel[this->m_TM2])
1745  && (newPixel[this->m_TM2] >= newPixel[this->m_TM3])
1746  && (newPixel[this->m_TM3] >= this->m_TV1 * newPixel[this->m_TM4])
1747  && (newPixel[this->m_TM1] >= newPixel[this->m_TM5])
1748  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM4])
1749  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM7]);
1750 
1751  return static_cast<TOutput>(result);
1752  }
1753 };
1754 
1755 
1771 template <class TInput, class TOutput>
1773 {
1774 public:
1775 
1776  typedef typename TInput::ValueType PrecisionType;
1777  typedef bool OutputPixelType;
1778 
1780  std::string GetName() const override
1781  {
1782  return "LandsatTM ShadowWithVegetationSpectralRule";
1783  }
1784 
1787 
1788  inline TOutput operator ()(const TInput& inputPixel)
1789  {
1790  TInput newPixel(this->PrepareValues( inputPixel ));
1791  bool result = (newPixel[this->m_TM1] >= newPixel[this->m_TM2])
1792  && (newPixel[this->m_TM2] >= newPixel[this->m_TM3])
1793  && (newPixel[this->m_TM1] >= this->m_TV2 * newPixel[this->m_TM4])
1794  && (newPixel[this->m_TM3] < this->m_TV1 * newPixel[this->m_TM4])
1795  && (newPixel[this->m_TM5] < this->m_TV1 * newPixel[this->m_TM4])
1796  && (newPixel[this->m_TM3] >= this->m_TV2 * newPixel[this->m_TM5])
1797  && (newPixel[this->m_TM7] < this->m_TV1 * newPixel[this->m_TM4]);
1798 
1799  return static_cast<TOutput>(result);
1800  }
1801 };
1802 
1803 
1819 template <class TInput, class TOutput>
1820 class ShadowCloudOrSnowSpectralRule : public KernelSpectralRule<TInput, TOutput>
1821 {
1822 public:
1823 
1824  typedef typename TInput::ValueType PrecisionType;
1825  typedef bool OutputPixelType;
1826 
1828  std::string GetName() const override
1829  {
1830  return "LandsatTM ShadowCloudOrSnowSpectralRule";
1831  }
1832 
1835 
1836  inline TOutput operator ()(const TInput& inputPixel)
1837  {
1838  TInput newPixel(this->PrepareValues( inputPixel ));
1839  PrecisionType max13;
1840  PrecisionType max123;
1841  PrecisionType min123;
1842  PrecisionType max12347;
1843  PrecisionType min12347;
1844  PrecisionType max234;
1845  PrecisionType max45;
1846  this->SetMinMax(newPixel, &max13, &min123, &max123, &min12347, &max12347, &max234, &max45);
1847 
1848 
1849  bool result = (newPixel[this->m_TM1] >= this->m_TV1 * max234)
1850  && (max234 >= this->m_TV1 * newPixel[this->m_TM1])
1851  && (newPixel[this->m_TM5] < newPixel[this->m_TM1])
1852  && (newPixel[this->m_TM7] < this->m_TV1 * newPixel[this->m_TM1]);
1853 
1854  return static_cast<TOutput>(result);
1855  }
1856 };
1857 
1858 
1874 template <class TInput, class TOutput>
1875 class WetlandSpectralRule : public KernelSpectralRule<TInput, TOutput>
1876 {
1877 public:
1878 
1879  typedef typename TInput::ValueType PrecisionType;
1880  typedef bool OutputPixelType;
1881 
1883  std::string GetName() const override
1884  {
1885  return "LandsatTM WetlandSpectralRule";
1886  }
1887 
1889  ~WetlandSpectralRule() override {}
1890 
1891  inline TOutput operator ()(const TInput& inputPixel)
1892  {
1893  TInput newPixel(this->PrepareValues( inputPixel ));
1894  bool result = (newPixel[this->m_TM1] >= newPixel[this->m_TM2])
1895  && (newPixel[this->m_TM2] >= newPixel[this->m_TM3])
1896  && (newPixel[this->m_TM1] >= this->m_TV1 * newPixel[this->m_TM4])
1897  && (newPixel[this->m_TM3] < newPixel[this->m_TM4])
1898  && (newPixel[this->m_TM4] >= this->m_TV1 * newPixel[this->m_TM5])
1899  && (newPixel[this->m_TM5] >= this->m_TV1 * newPixel[this->m_TM4])
1900  && (newPixel[this->m_TM3] >= this->m_TV2 * newPixel[this->m_TM5])
1901  && (newPixel[this->m_TM5] >= newPixel[this->m_TM7]);
1902 
1903  return static_cast<TOutput>(result);
1904  }
1905 };
1906 
1907 
1908 } // namespace LandsatTM
1909 } // namespace Functor
1910 } // namespace otb
1911 
1912 #endif
TOutput operator()(const TInput &inputPixel)
TOutput operator()(const TInput &inputPixel)
std::string GetName() const override
std::string GetName() const override
std::string GetName() const override
TOutput operator()(const TInput &inputPixel)
itk::FixedArray< unsigned int, 11 > operator()(const TInput &inputPixel)
void SetIndex(BandName::LandsatTMBandNames band, unsigned int channel)
TOutput operator()(const TInput &inputPixel)
TOutput operator()(const TInput &inputPixel)
std::string GetName() const override
DegreeType
Thermal band in Kelvin or Celsius.
TOutput operator()(const TInput &inputPixel)
TOutput operator()(const TInput &inputPixel)
otb::FuzzyVariable< unsigned short, PrecisionType > FuzzyVarType
std::string GetName() const override
TOutput operator()(const TInput &inputPixel)
TOutput operator()(const TInput &inputPixel)
std::string GetName() const override
TInput PrepareValues(const TInput &inputPixel)
Prepare the values so they are normalized and in C.
TOutput operator()(const TInput &inputPixel) const
bool operator==(const LandsatTMIndexBase &other) const
itk::FixedArray< unsigned int, 11 > OutputPixelType
ReflectanceType
Reflectances in thousands or in (0-1)
static ITK_CONSTEXPR_FUNC T max(const T &)
unsigned int GetIndex(BandName::LandsatTMBandNames band) const
TOutput operator()(const TInput &inputPixel)
std::string GetName() const override
TOutput operator()(const TInput &inputPixel)
bool operator!=(const LandsatTMIndexBase &) const
TOutput operator()(const TInput &inputPixel)
TOutput operator()(const TInput &inputPixel)
virtual std::string GetName() const =0
static ITK_CONSTEXPR_FUNC T NonpositiveMin()
std::string GetName() const override
TOutput operator()(const TInput &inputPixel)
std::string GetName() const override
std::string GetName() const override
static Pointer New()
Base class for Landsat-TM indices.
void SetMinMax(const TInput &inputPixel, PrecisionType *max13, PrecisionType *min123, PrecisionType *max123, PrecisionType *min12347, PrecisionType *max12347, PrecisionType *max234, PrecisionType *max45)
std::string GetName() const override
Class to represent a fuzzy N-valued variable.
std::string GetName() const override
TOutput operator()(const TInput &inputPixel)
TOutput operator()(const TInput &inputPixel)
std::string GetName() const override