Orfeo Toolbox  3.16
otbVegetationIndicesFunctor.h
Go to the documentation of this file.
1 /*=========================================================================
2 
3  Program: ORFEO Toolbox
4  Language: C++
5  Date: $Date$
6  Version: $Revision$
7 
8 
9  Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
10  See OTBCopyright.txt for details.
11 
12 
13  This software is distributed WITHOUT ANY WARRANTY; without even
14  the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
15  PURPOSE. See the above copyright notices for more information.
16 
17 =========================================================================*/
18 #ifndef __otbVegetationIndicesFunctor_h
19 #define __otbVegetationIndicesFunctor_h
20 
21 #include "otbMath.h"
22 #include "itkVariableLengthVector.h"
23 #include "otbBandName.h"
24 
25 namespace otb
26 {
27 
28 namespace Functor
29 {
30 
43 template<class TInput1, class TInput2, class TOutput>
44 class RAndNIRIndexBase
45 {
46 public:
49  typedef itk::VariableLengthVector<TInput1> InputVectorType;
50 
51  //operators !=
52  bool operator !=(const RAndNIRIndexBase&) const
53  {
54  return true;
55  }
56  //operator ==
57  bool operator ==(const RAndNIRIndexBase& other) const
58  {
59  return !(*this != other);
60  }
61 
62  // Operator on vector pixel type
63  inline TOutput operator ()(const InputVectorType& inputVector) const
64  {
65  return this->Evaluate(inputVector[m_RedIndex - 1], static_cast<TInput2>(inputVector[m_NIRIndex - 1]));
66  }
67 
68  // Binary operator
69  inline TOutput operator ()(const TInput1& r, const TInput2& nir) const
70  {
71  return this->Evaluate(r, nir);
72  }
74  RAndNIRIndexBase() : m_EpsilonToBeConsideredAsZero(0.0000001), m_RedIndex(3), m_NIRIndex(4) {}
76  virtual ~RAndNIRIndexBase() {}
77 
79  void SetIndex(BandName::BandName band, unsigned int channel)
80  {
81  if (band == BandName::RED)
82  {
83  m_RedIndex = channel;
84  }
85  if (band == BandName::NIR)
86  {
87  m_NIRIndex = channel;
88  }
89  }
91  unsigned int GetIndex(BandName::BandName band) const
92  {
93  if (band == BandName::RED)
94  {
95  return m_RedIndex;
96  }
97  if (band == BandName::NIR)
98  {
99  return m_NIRIndex;
100  }
101  }
102 
104  void SetRedIndex(unsigned int channel)
105  {
106  m_RedIndex = channel;
107  }
109  unsigned int GetRedIndex() const
110  {
111  return m_RedIndex;
112  }
114  void SetNIRIndex(unsigned int channel)
115  {
116  m_NIRIndex = channel;
117  }
119  unsigned int GetNIRIndex() const
120  {
121  return m_NIRIndex;
122  }
123 
125  virtual std::string GetName() const = 0;
126 
127 protected:
128  // This method must be reimplemented in subclasses to actually
129  // compute the index value
130  virtual TOutput Evaluate(const TInput1& r, const TInput2& nir) const = 0;
131  double m_EpsilonToBeConsideredAsZero;
132 
133 private:
134  unsigned int m_RedIndex;
135  unsigned int m_NIRIndex;
136 };
137 
148 template<class TInput1, class TInput2, class TInput3, class TOutput>
149 class RAndBAndNIRIndexBase
150 {
151 public:
154  typedef itk::VariableLengthVector<TInput1> InputVectorType;
155 
156  //operators !=
157  bool operator !=(const RAndBAndNIRIndexBase&) const
158  {
159  return true;
160  }
161 
162  //operator ==
163  bool operator ==(const RAndBAndNIRIndexBase& other) const
164  {
165  return !(*this != other);
166  }
167 
168  // Operator on vector pixel type
169  inline TOutput operator ()(const InputVectorType& inputVector)
170  {
171  return this->Evaluate(inputVector[m_RedIndex - 1],
172  static_cast<TInput2>(inputVector[m_BlueIndex - 1]),
173  static_cast<TInput3>(inputVector[m_NIRIndex - 1]));
174  }
175  // Binary operator
176  inline TOutput operator ()(const TInput1& r, const TInput2& b, const TInput2& nir)
177  {
178  return this->Evaluate(r, b, nir);
179  }
181  RAndBAndNIRIndexBase() : m_EpsilonToBeConsideredAsZero(0.0000001), m_RedIndex(3), m_BlueIndex(1), m_NIRIndex(4) {}
183  virtual ~RAndBAndNIRIndexBase() {}
184 
186  void SetIndex(BandName::BandName band, unsigned int channel)
187  {
188  if (band == BandName::RED)
189  {
190  m_RedIndex = channel;
191  }
192  if (band == BandName::BLUE)
193  {
194  m_BlueIndex = channel;
195  }
196  if (band == BandName::NIR)
197  {
198  m_NIRIndex = channel;
199  }
200  }
202  unsigned int GetIndex(BandName::BandName band) const
203  {
204  if (band == BandName::RED)
205  {
206  return m_RedIndex;
207  }
208  if (band == BandName::BLUE)
209  {
210  return m_BlueIndex;
211  }
212  if (band == BandName::NIR)
213  {
214  return m_NIRIndex;
215  }
216  }
217 
219  void SetRedIndex(unsigned int channel)
220  {
221  m_RedIndex = channel;
222  }
224  unsigned int GetRedIndex() const
225  {
226  return m_RedIndex;
227  }
229  void SetBlueIndex(unsigned int channel)
230  {
231  m_BlueIndex = channel;
232  }
234  unsigned int GetBlueIndex() const
235  {
236  return m_BlueIndex;
237  }
238 
240  void SetNIRIndex(unsigned int channel)
241  {
242  m_NIRIndex = channel;
243  }
245  unsigned int GetNIRIndex() const
246  {
247  return m_NIRIndex;
248  }
249 
251  virtual std::string GetName() const = 0;
252 
253 protected:
254  // This method must be reimplemented in subclasses to actually
255  // compute the index value
256  virtual TOutput Evaluate(const TInput1& r, const TInput2& b, const TInput3& nir) const = 0;
257  double m_EpsilonToBeConsideredAsZero;
258 
259 private:
260  unsigned int m_RedIndex;
261  unsigned int m_BlueIndex;
262  unsigned int m_NIRIndex;
263 };
264 
275 template<class TInput1, class TInput2, class TInput3, class TOutput>
276 class RAndGAndNIRIndexBase
277 {
278 public:
281  typedef itk::VariableLengthVector<TInput1> InputVectorType;
282 
283  //operators !=
284  bool operator !=(const RAndGAndNIRIndexBase&) const
285  {
286  return true;
287  }
288  //operator ==
289  bool operator ==(const RAndGAndNIRIndexBase& other) const
290  {
291  return !(*this != other);
292  }
293 
294  // Operator on vector pixel type
295  inline TOutput operator ()(const InputVectorType& inputVector)
296  {
297  return this->Evaluate(inputVector[m_RedIndex - 1],
298  static_cast<TInput2>(inputVector[m_GreenIndex - 1]),
299  static_cast<TInput3>(inputVector[m_NIRIndex - 1]));
300  }
301 
302  // Binary operator
303  inline TOutput operator ()(const TInput1& r, const TInput2& g, const TInput2& nir)
304  {
305  return this->Evaluate(r, g, nir);
306  }
308  RAndGAndNIRIndexBase() : m_EpsilonToBeConsideredAsZero(0.0000001), m_RedIndex(3), m_GreenIndex(2), m_NIRIndex(4) {}
310  virtual ~RAndGAndNIRIndexBase() {}
311 
313  void SetIndex(BandName::BandName band, unsigned int channel)
314  {
315  if (band == BandName::RED)
316  {
317  m_RedIndex = channel;
318  }
319  if (band == BandName::GREEN)
320  {
321  m_GreenIndex = channel;
322  }
323  if (band == BandName::NIR)
324  {
325  m_NIRIndex = channel;
326  }
327  }
329  unsigned int GetIndex(BandName::BandName band) const
330  {
331  if (band == BandName::RED)
332  {
333  return m_RedIndex;
334  }
335  if (band == BandName::GREEN)
336  {
337  return m_GreenIndex;
338  }
339  if (band == BandName::NIR)
340  {
341  return m_NIRIndex;
342  }
343  }
344 
346  void SetRedIndex(unsigned int channel)
347  {
348  m_RedIndex = channel;
349  }
351  unsigned int GetRedIndex() const
352  {
353  return m_RedIndex;
354  }
356  void SetGreenIndex(unsigned int channel)
357  {
358  m_GreenIndex = channel;
359  }
361  unsigned int GetGreenIndex() const
362  {
363  return m_GreenIndex;
364  }
365 
367  void SetNIRIndex(unsigned int channel)
368  {
369  m_NIRIndex = channel;
370  }
372  unsigned int GetNIRIndex() const
373  {
374  return m_NIRIndex;
375  }
376 
378  virtual std::string GetName() const = 0;
379 
380 protected:
381  // This method must be reimplemented in subclasses to actually
382  // compute the index value
383  virtual TOutput Evaluate(const TInput1& r, const TInput2& g, const TInput3& nir) const = 0;
384  double m_EpsilonToBeConsideredAsZero;
385 
386 private:
387  unsigned int m_RedIndex;
388  unsigned int m_GreenIndex;
389  unsigned int m_NIRIndex;
390 };
391 
400 template <class TInput1, class TInput2, class TOutput>
401 class NDVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
402 {
403 public:
404 
406  virtual std::string GetName() const
407  {
408  return "NDVI";
409  }
410 
412  NDVI() {}
414  virtual ~NDVI() {}
415  // Operator on r and nir single pixel values
416 protected:
417  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
418  {
419  double dr = static_cast<double>(r);
420  double dnir = static_cast<double>(nir);
421  if (vcl_abs(dnir + dr) < this->m_EpsilonToBeConsideredAsZero)
422  {
423  return static_cast<TOutput>(0.);
424  }
425 
426  return (static_cast<TOutput>((dnir - dr) / (dnir + dr)));
427  }
428 };
429 
438 template <class TInput1, class TInput2, class TOutput>
439 class RVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
440 {
441 public:
442 
444  virtual std::string GetName() const
445  {
446  return "RVI";
447  }
448 
449  RVI() {}
450  virtual ~RVI() {}
451 protected:
452  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
453  {
454  double dr = static_cast<double>(r);
455  double dnir = static_cast<double>(nir);
456  if (vcl_abs(dr) < this->m_EpsilonToBeConsideredAsZero)
457  {
458  return static_cast<TOutput>(0.);
459  }
460  return (static_cast<TOutput>(dnir / dr));
461  }
462 };
463 
476 template <class TInput1, class TInput2, class TOutput>
477 class PVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
478 {
479 public:
481  virtual std::string GetName() const
482  {
483  return "PVI";
484  }
485 
486  PVI() : m_A(0.90893), m_B(7.46216), m_Coeff(0.74) {}
487  virtual ~PVI() {}
489  void SetA(const double A)
490  {
491  m_A = A;
492  m_Coeff = 1. / (vcl_sqrt(m_A * m_A + 1.));
493  }
494  double GetA(void) const
495  {
496  return (m_A);
497  }
498  void SetB(const double B)
499  {
500  m_B = B;
501  }
502  double GetB(void) const
503  {
504  return (m_B);
505  }
506 protected:
507  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
508  {
509  double dnir = static_cast<double>(nir);
510  double dr = static_cast<double>(r);
511  return (static_cast<TOutput>((dnir - m_A * dr - m_B) * m_Coeff));
512  }
513 
514 private:
515 
517  double m_A;
518  double m_B;
520  double m_Coeff;
521 
522 };
523 
532 template <class TInput1, class TInput2, class TOutput>
533 class SAVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
534 {
535 public:
536 
538  virtual std::string GetName() const
539  {
540  return "SAVI";
541  }
542 
543  SAVI() : m_L(0.5) {}
544  virtual ~SAVI() {}
545 
547  void SetL(const double L)
548  {
549  m_L = L;
550  }
551  double GetL(void) const
552  {
553  return (m_L);
554  }
555 
556 protected:
557  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
558  {
559  double dnir = static_cast<double>(nir);
560  double dr = static_cast<double>(r);
561  double denominator = dnir + dr + m_L;
562  if (vcl_abs(denominator) < this->m_EpsilonToBeConsideredAsZero)
563  {
564  return static_cast<TOutput>(0.);
565  }
566  return (static_cast<TOutput>(((dnir - dr) * (1 + m_L)) / denominator));
567  }
568 
569 private:
570 
572  double m_L;
573 
574 };
575 
584 template <class TInput1, class TInput2, class TOutput>
585 class TSAVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
586 {
587 public:
588 
590  virtual std::string GetName() const
591  {
592  return "TSAVI";
593  }
594 
595  TSAVI() : m_A(0.7), m_S(0.9), m_X(0.08) {}
596  virtual ~TSAVI() {}
597 
599  void SetS(const double S)
600  {
601  m_S = S;
602  }
603  double GetS(void) const
604  {
605  return (m_S);
606  }
607  void SetA(const double A)
608  {
609  m_A = A;
610  }
611  double GetA(void) const
612  {
613  return (m_A);
614  }
616  void SetX(const double X)
617  {
618  m_X = X;
619  }
620  double GetX(void) const
621  {
622  return (m_X);
623  }
624 
625 protected:
626  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
627  {
628  double dnir = static_cast<double>(nir);
629  double dr = static_cast<double>(r);
630  double denominator = m_A * dnir + dr + m_X * (1. + m_A * m_A);
631  if (vcl_abs(denominator) < this->m_EpsilonToBeConsideredAsZero)
632  {
633  return static_cast<TOutput>(0.);
634  }
635  return (static_cast<TOutput>((m_A * (dnir - m_A * dr - m_S)) / denominator));
636  }
637 
638 private:
639 
641  double m_A;
642  double m_S;
644  double m_X;
645 
646 };
647 
656 template <class TInput1, class TInput2, class TOutput>
657 class WDVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
658 {
659 public:
661  virtual std::string GetName() const
662  {
663  return "WDVI";
664  }
665 
667  WDVI() : m_S(0.4) {}
669  virtual ~WDVI() {}
670  // Operator on r and nir single pixel values
672  void SetS(const double s)
673  {
674  m_S = s;
675  }
676  double GetS(void) const
677  {
678  return (m_S);
679  }
680 protected:
681  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
682  {
683  double dr = static_cast<double>(r);
684  double dnir = static_cast<double>(nir);
685 
686  return (dnir - m_S * dr);
687  }
688 private:
690  double m_S;
691 };
692 
702 template <class TInput1, class TInput2, class TOutput>
703 class MSAVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
704 {
705 public:
707  virtual std::string GetName() const
708  {
709  return "MSAVI";
710  }
711 
712  typedef NDVI<TInput1, TInput2, TOutput> NDVIFunctorType;
713  typedef SAVI<TInput1, TInput2, TOutput> SAVIFunctorType;
714  typedef WDVI<TInput1, TInput2, TOutput> WDVIFunctorType;
715  MSAVI() : m_S(0.4)
716  {
717  m_WDVIfunctor.SetS(m_S);
718  }
719  virtual ~MSAVI() {}
721  void SetS(const double s)
722  {
723  m_S = s;
724  m_WDVIfunctor.SetS(m_S);
725  }
726  double GetS(void) const
727  {
728  return (m_S);
729  }
730  NDVIFunctorType GetNDVI(void) const
731  {
732  return (m_NDVIfunctor);
733  }
734  WDVIFunctorType GetWDVI(void) const
735  {
736  return (m_WDVIfunctor);
737  }
738 
739 protected:
740  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
741  {
742  double dnir = static_cast<double>(nir);
743  double dr = static_cast<double>(r);
744 
745  double dNDVI = this->GetNDVI() (r, nir);
746  double dWDVI = this->GetWDVI() (r, nir);
747  double dL = 1 - 2 * m_S * dNDVI * dWDVI;
748 
749  double denominator = dnir + dr + dL;
750 
751  if (vcl_abs(denominator) < this->m_EpsilonToBeConsideredAsZero)
752  {
753  return static_cast<TOutput>(0.);
754  }
755 
756  return (static_cast<TOutput>(((dnir - dr) * (1 + dL)) / denominator));
757  }
758 
759 private:
761  double m_S;
762  NDVIFunctorType m_NDVIfunctor;
763  WDVIFunctorType m_WDVIfunctor;
764 
765 };
766 
775 template <class TInput1, class TInput2, class TOutput>
776 class MSAVI2 : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
777 {
778 public:
780  virtual std::string GetName() const
781  {
782  return "MSAVI2";
783  }
784 
785  MSAVI2() {}
786  virtual ~MSAVI2() {}
787 
788 protected:
789  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
790  {
791  double dnir = static_cast<double>(nir);
792  double dr = static_cast<double>(r);
793  double sqrt_value = (2 * dnir + 1) * (2 * dnir + 1) - 8 * (dnir - dr);
794  if (sqrt_value < 0.)
795  {
796  return static_cast<TOutput>(0.);
797  }
798  return (static_cast<TOutput>((2 * dnir + 1 - vcl_sqrt(sqrt_value)) / 2.));
799  }
800 
801 };
802 
811 template <class TInput1, class TInput2, class TOutput>
812 class GEMI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
813 {
814 public:
816  virtual std::string GetName() const
817  {
818  return "GEMI";
819  }
820 
821  GEMI() {}
822  virtual ~GEMI() {}
823 
824 protected:
825  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
826  {
827  double dnir = static_cast<double>(nir);
828  double dr = static_cast<double>(r);
829 
830  double dnu;
831  double dnumerateur_nu;
832  double ddenominateur_nu = dnir + dr + 0.5;
833  if (vcl_abs(ddenominateur_nu) < this->m_EpsilonToBeConsideredAsZero)
834  {
835  dnu = 0;
836  }
837  else
838  {
839  dnumerateur_nu = 2 * (dnir * dnir - dr * dr) + 1.5 * dnir + 0.5 * dr;
840  dnu = dnumerateur_nu / ddenominateur_nu;
841  }
842 
843  double ddenominateur_GEMI = 1 - dr;
844  if (vcl_abs(ddenominateur_GEMI) < this->m_EpsilonToBeConsideredAsZero)
845  {
846  return static_cast<TOutput>(0.);
847  }
848  return (static_cast<TOutput>((dnu * (1 - 0.25 * dnu) - (dr - 0.125)) / ddenominateur_GEMI));
849  }
850 
851 };
852 
863 template <class TInput1, class TInput2, class TInput3, class TOutput>
864 class AVI : public RAndGAndNIRIndexBase<TInput1, TInput2, TInput3, TOutput>
865 {
866 public:
868  virtual std::string GetName() const
869  {
870  return "AVI";
871  }
872 
873  AVI() : m_LambdaG(560.), m_LambdaR(660.), m_LambdaNir(830.) {}
874  virtual ~AVI() {}
876  void SetLambdaR(const double lr)
877  {
878  m_LambdaR = lr;
879  }
880  double GetLambdaR(void) const
881  {
882  return (m_LambdaR);
883  }
885  void SetLambdaG(const double lg)
886  {
887  m_LambdaG = lg;
888  }
889  double GetLambdaG(void) const
890  {
891  return (m_LambdaG);
892  }
894  void SetLambdaNir(const double lnir)
895  {
896  m_LambdaNir = lnir;
897  }
898  double GetLambdaNir(void) const
899  {
900  return (m_LambdaNir);
901  }
902 protected:
903  inline TOutput Evaluate(const TInput1& r, const TInput2& g, const TInput3& nir) const
904  {
905  double dr = static_cast<double>(r);
906  double dg = static_cast<double>(g);
907  double dnir = static_cast<double>(nir);
908 
909  double dfact1 = (m_LambdaNir - m_LambdaR) / m_LambdaR;
910  double dfact2 = (m_LambdaR - m_LambdaG) / m_LambdaR;
911  double dterm1;
912  double dterm2;
913  if (vcl_abs(dnir - dr) < this->m_EpsilonToBeConsideredAsZero)
914  {
915  dterm1 = 0;
916  }
917  else
918  {
919  dterm1 = vcl_atan(dfact1 / (dnir - dr));
920  }
921 
922  if (vcl_abs(dg - dr) < this->m_EpsilonToBeConsideredAsZero)
923  {
924  dterm2 = 0;
925  }
926  else
927  {
928  dterm2 = vcl_atan(dfact2 / (dg - dr));
929  }
930 
931  return static_cast<TOutput>(dterm1 + dterm2);
932 
933  }
934 private:
935 
937  double m_LambdaG;
938 
940  double m_LambdaR;
941 
943  double m_LambdaNir;
944 };
945 
956 template <class TInput1, class TInput2, class TInput3, class TOutput>
957 class ARVI : public RAndBAndNIRIndexBase<TInput1, TInput2, TInput3, TOutput>
958 {
959 public:
961  virtual std::string GetName() const
962  {
963  return "ARVI";
964  }
965 
966  ARVI() : m_Gamma(0.5) {}
967  virtual ~ARVI() {}
968 
970  void SetGamma(const double gamma)
971  {
972  m_Gamma = gamma;
973  }
974  double GetGamma(void) const
975  {
976  return (m_Gamma);
977  }
978 
979 protected:
980  inline TOutput Evaluate(const TInput1& r, const TInput2& b, const TInput3& nir) const
981  {
982  double dr = static_cast<double>(r);
983  double db = static_cast<double>(b);
984  double dnir = static_cast<double>(nir);
985  double RHOrb = dr - m_Gamma * (db - dr);
986  double denominator = dnir + RHOrb;
987  if (vcl_abs(denominator) < this->m_EpsilonToBeConsideredAsZero)
988  {
989  return static_cast<TOutput>(0.);
990  }
991  return (static_cast<TOutput>((dnir - RHOrb) / denominator));
992  }
993 
994 private:
995 
997  double m_Gamma;
998 };
999 
1008 template <class TInput1, class TInput2, class TInput3, class TOutput>
1009 class TSARVI : public RAndBAndNIRIndexBase<TInput1, TInput2, TInput3, TOutput>
1010 {
1011 public:
1013  virtual std::string GetName() const
1014  {
1015  return "TSARVI";
1016  }
1017 
1018  TSARVI() : m_A(0.0), m_B(0.0), m_X(0.08), m_Gamma(0.5) {}
1019  virtual ~TSARVI() {}
1020 
1022  void SetA(const double A)
1023  {
1024  m_A = A;
1025  }
1026  double GetA(void) const
1027  {
1028  return (m_A);
1029  }
1030  void SetB(const double B)
1031  {
1032  m_B = B;
1033  }
1034  double GetB(void) const
1035  {
1036  return (m_B);
1037  }
1039  void SetX(const double X)
1040  {
1041  m_X = X;
1042  }
1043  double GetX(void) const
1044  {
1045  return (m_X);
1046  }
1048  void SetGamma(const double gamma)
1049  {
1050  m_Gamma = gamma;
1051  }
1052  double GetGamma(void) const
1053  {
1054  return (m_Gamma);
1055  }
1056 
1057 protected:
1058  inline TOutput Evaluate(const TInput1& r, const TInput2& b, const TInput3& nir) const
1059  {
1060  double dr = static_cast<double>(r);
1061  double db = static_cast<double>(b);
1062  double dnir = static_cast<double>(nir);
1063  double dRB = dr - m_Gamma * (db - dr);
1064  double denominator = dRB + m_A * dnir - m_A * m_B + m_X * (1. + m_A * m_A);
1065  if (vcl_abs(denominator) < this->m_EpsilonToBeConsideredAsZero)
1066  {
1067  return static_cast<TOutput>(0.);
1068  }
1069  return (static_cast<TOutput>((m_A * (dnir - m_A * dRB - m_B)) / denominator));
1070  }
1071 
1072 private:
1073 
1075  double m_A;
1076  double m_B;
1078  double m_X;
1080  double m_Gamma;
1081 
1082 };
1083 
1094 template <class TInput1, class TInput2, class TInput3, class TOutput>
1095 class EVI : public RAndBAndNIRIndexBase<TInput1, TInput2, TInput3, TOutput>
1096 {
1097 public:
1099  virtual std::string GetName() const
1100  {
1101  return "EVI";
1102  }
1103 
1104  EVI() : m_G(2.5), m_C1(6.0), m_C2(7.5), m_L(1.0) {}
1105  virtual ~EVI() {}
1107  void SetG(const double g)
1108  {
1109  m_G = g;
1110  }
1111  double GetG(void) const
1112  {
1113  return (m_G);
1114  }
1116  void SetC1(const double c1)
1117  {
1118  m_C1 = c1;
1119  }
1120  double GetC1(void) const
1121  {
1122  return (m_C1);
1123  }
1125  void SetC2(const double c2)
1126  {
1127  m_C2 = c2;
1128  }
1129  double GetC2(void) const
1130  {
1131  return (m_C2);
1132  }
1134  void SetL(const double l)
1135  {
1136  m_L = l;
1137  }
1138  double GetL(void) const
1139  {
1140  return (m_L);
1141  }
1142 protected:
1143  inline TOutput Evaluate(const TInput1& r, const TInput2& b, const TInput3& nir) const
1144  {
1145  double dr = static_cast<double>(r);
1146  double db = static_cast<double>(b);
1147  double dnir = static_cast<double>(nir);
1148  double denominator = dnir + m_C1 * dr - m_C2 * db + m_L;
1149  if (vcl_abs(denominator) < this->m_EpsilonToBeConsideredAsZero)
1150  {
1151  return (static_cast<TOutput>(0.));
1152  }
1153  return (static_cast<TOutput>(m_G * (dnir - dr) / denominator));
1154  }
1155 
1156 private:
1157 
1159  double m_G;
1160 
1162  double m_C1;
1163 
1165  double m_C2;
1166 
1168  double m_L;
1169 };
1170 
1179 template <class TInput1, class TInput2, class TOutput>
1180 class IPVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
1181 {
1182 public:
1184  virtual std::string GetName() const
1185  {
1186  return "IPVI";
1187  }
1188 
1189  IPVI() {}
1190  virtual ~IPVI() {}
1191 
1192 protected:
1193  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
1194  {
1195  double dr = static_cast<double>(r);
1196  double dnir = static_cast<double>(nir);
1197  if (vcl_abs(dnir + dr) < this->m_EpsilonToBeConsideredAsZero)
1198  {
1199  return static_cast<TOutput>(0.);
1200  }
1201  else
1202  {
1203  return (static_cast<TOutput>(dnir / (dnir + dr)));
1204  }
1205  }
1206 };
1207 
1216 template <class TInput1, class TInput2, class TOutput>
1217 class TNDVI : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
1218 {
1219 public:
1221  virtual std::string GetName() const
1222  {
1223  return "TNDVI";
1224  }
1225 
1226  typedef NDVI<TInput1, TInput2, TOutput> NDVIFunctorType;
1227  TNDVI() {}
1228  virtual ~TNDVI() {}
1229 
1230  NDVIFunctorType GetNDVI(void) const
1231  {
1232  return (m_NDVIfunctor);
1233  }
1234 
1235 protected:
1236  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
1237  {
1238  double dval = this->GetNDVI() (r, nir) + 0.5;
1239  if (dval < 0)
1240  {
1241  return (static_cast<TOutput>(0));
1242  }
1243  else
1244  {
1245  return (static_cast<TOutput>(vcl_sqrt(dval)));
1246  }
1247  }
1248 private:
1249  NDVIFunctorType m_NDVIfunctor;
1250 };
1251 
1266 template <class TInput1, class TInput2, class TOutput>
1267 class LAIFromNDVILogarithmic : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
1268 {
1269 public:
1271  virtual std::string GetName() const
1272  {
1273  return "LAIFromNDVILogarithmic";
1274  }
1275 
1276  typedef NDVI<TInput1, TInput2, TOutput> NDVIFunctorType;
1277  LAIFromNDVILogarithmic() : m_NdviSoil(0.10), m_NdviInf(0.89), m_ExtinctionCoefficient(0.71) {}
1278  virtual ~LAIFromNDVILogarithmic() {}
1279 
1280  NDVIFunctorType GetNDVI(void) const
1281  {
1282  return (m_NDVIfunctor);
1283  }
1284 
1285  void SetNdviSoil(const double val)
1286  {
1287  m_NdviSoil = val;
1288  }
1289  double GetNdviSoil(void) const
1290  {
1291  return (m_NdviSoil);
1292  }
1293 
1294  void SetNdviInf(const double val)
1295  {
1296  m_NdviInf = val;
1297  }
1298  double GetNdviInf(void) const
1299  {
1300  return (m_NdviInf);
1301  }
1302 
1303  void SetExtinctionCoefficient(const double val)
1304  {
1305  m_ExtinctionCoefficient = val;
1306  }
1307  double GetExtinctionCoefficient(void) const
1308  {
1309  return (m_ExtinctionCoefficient);
1310  }
1311 
1312 protected:
1313  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
1314  {
1315  double dval = this->GetNDVI() (r, nir);
1316  if (dval < 0)
1317  {
1318  return (static_cast<TOutput>(0));
1319  }
1320  else
1321  {
1322  return (static_cast<TOutput>(
1323  -(1.0/m_ExtinctionCoefficient)*vcl_log((dval- m_NdviInf)/(m_NdviSoil-m_NdviInf))
1324  ));
1325  }
1326  }
1327 private:
1328  NDVIFunctorType m_NDVIfunctor;
1329  double m_NdviSoil;
1330  double m_NdviInf;
1331  double m_ExtinctionCoefficient;
1332 };
1333 
1334 
1350 template <class TInput1, class TInput2, class TOutput>
1351 class LAIFromReflectancesLinear : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
1352 {
1353 public:
1355  virtual std::string GetName() const
1356  {
1357  return "LAIFromReflectancesLinear";
1358  }
1359 
1360  typedef NDVI<TInput1, TInput2, TOutput> NDVIFunctorType;
1361  LAIFromReflectancesLinear() : m_RedCoef(-17.91), m_NirCoef(12.26) {}
1362  virtual ~LAIFromReflectancesLinear() {}
1363 
1364  NDVIFunctorType GetReflectances(void) const
1365  {
1366  return (m_NDVIfunctor);
1367  }
1368 
1369  void SetRedCoef(const double val)
1370  {
1371  m_RedCoef = val;
1372  }
1373  double GetRedCoef(void) const
1374  {
1375  return (m_RedCoef);
1376  }
1377 
1378  void SetNirCoef(const double val)
1379  {
1380  m_NirCoef = val;
1381  }
1382  double GetNirCoef(void) const
1383  {
1384  return (m_NirCoef);
1385  }
1386 
1387 protected:
1388  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
1389  {
1390  return (static_cast<TOutput>(m_RedCoef*r+m_NirCoef*nir));
1391  }
1392 private:
1393  NDVIFunctorType m_NDVIfunctor;
1394  double m_RedCoef;
1395  double m_NirCoef;
1396 };
1397 
1398 
1416  template <class TInput1, class TInput2, class TOutput>
1417  class LAIFromNDVIFormosat2Functor : public RAndNIRIndexBase<TInput1, TInput2, TOutput>
1418  {
1419  public:
1420 
1422  virtual std::string GetName() const
1423  {
1424  return "LAIFromNDVIFormosat2Functor";
1425  }
1426 
1428  LAIFromNDVIFormosat2Functor() {}
1430  virtual ~LAIFromNDVIFormosat2Functor() {}
1431  // Operator on r and nir single pixel values
1432  protected:
1433  inline TOutput Evaluate(const TInput1& r, const TInput2& nir) const
1434  {
1435  double a = 0.1519;
1436  double b = 3.9443;
1437  double c = 0.13;
1438 
1439  double dr = static_cast<double>(r);
1440  double dnir = static_cast<double>(nir);
1441  if (vcl_abs(dnir + dr) < this->m_EpsilonToBeConsideredAsZero)
1442  {
1443  return static_cast<TOutput>(0.);
1444  }
1445 
1446  return static_cast<TOutput>(a*(vcl_exp(static_cast<double>(dnir-dr)/static_cast<double>(dr+dnir)*b)-vcl_exp(c*b)));
1447  };
1448 
1449 
1450 };
1451 
1452 
1453 } // namespace Functor
1454 } // namespace otb
1455 
1456 #endif
Generated at Sun Feb 3 2013 00:54:58 for Orfeo Toolbox with doxygen 1.8.1.1