Skip to content

GitLab

Commit 7ba69067 authored by remi cresson's avatar remi cresson
Browse files

ADD: new functors SlopeAndPValueLabelingFunctor and RainfallEstimatedNDVIResiduesFunctor

parent 84c4d7d2
Hide whitespace changes
Inline Side-by-side
include/otbNDVITimeSeriesFunctor.h
View file @ 7ba69067
......@@ -84,12 +84,12 @@ public:
typedef typename TInputPixel::ValueType PixelValueType;
TSReduceFunctorBase() {}
virtual ~TSReduceFunctorBase() {}
~TSReduceFunctorBase() {}
bool operator!=( const TSReduceFunctorBase & ) const { return false; }
bool operator==( const TSReduceFunctorBase & other ) const { return !(*this != other); }
virtual inline TOutputPixel operator()( const TInputPixel & A ) const = 0;
inline unsigned int GetOutputSize(){ return nbOfYears; }
inline unsigned int GetOutputSize() const { return nbOfYears; }
/* Set dates */
void SetDates(DatesType dates)
......@@ -110,11 +110,11 @@ public:
}
/* Get dates */
DatesType GetDates(){ return m_Dates; }
DatesType GetDates() const { return m_Dates; }
/* Set / Get input no-data value */
void SetInputNoDataValue(PixelValueType value) { m_InputNoDataValue = value; }
PixelValueType GetInputNoDataValue(){ return m_InputNoDataValue; }
PixelValueType GetInputNoDataValue() const { return m_InputNoDataValue; }
private:
unsigned int nbOfYears;
......@@ -132,7 +132,7 @@ private:
*
*/
template< class TInputPixel, class TOutputPixel>
class TSCumulatedRangeReduceFunctor : TSReduceFunctorBase<TInputPixel,TOutputPixel>
class TSCumulatedRangeReduceFunctor : public TSReduceFunctorBase<TInputPixel,TOutputPixel>
{
public:
......@@ -157,27 +157,23 @@ inline TOutputPixel operator ()(const TInputPixel& input) const
output.SetSize(this->GetOutputSize());
output.Fill(this->GetInputNoDataValue());
DatesType m_Dates = this->GetDates();
DatesType dates = this->GetDates();
int firstYear = m_Dates[0].year;
int firstYear = dates[0].year;
for (unsigned int i = 0; i < input.Size(); i++)
{
int index = m_Dates[i].year - firstYear;
int m = m_Dates[i].month;
int index = dates[i].year - firstYear;
PixelValueType inValue = input[i];
// Range check
if (inValue != this->GetInputNoDataValue())
{
if (m_Month1 <= m && m <= m_Month2)
if (dates[i].IsInRange(m_Month1, m_Day1, m_Month2, m_Day2))
{
int d = m_Dates[i].day;
if (m_Day1 <= d && d <= m_Day2)
{
output[index] += inValue;
}
output[index] += inValue;
}
}
}
return output;
}
......@@ -203,7 +199,7 @@ int m_Month2;
*
*/
template< class TInputPixel, class TOutputPixel>
class TSMaxReduceFunctor : TSReduceFunctorBase<TInputPixel,TOutputPixel>
class TSMaxReduceFunctor : public TSReduceFunctorBase<TInputPixel,TOutputPixel>
{
public:
......@@ -228,12 +224,12 @@ inline TOutputPixel operator ()(const TInputPixel& input) const
output.SetSize(this->GetOutputSize());
output.Fill(this->GetInputNoDataValue());
DatesType m_Dates = this->GetDates();
DatesType dates = this->GetDates();
int firstYear = m_Dates[0].year;
int firstYear = dates[0].year;
for (unsigned int i = 0; i < input.Size(); i++)
{
int index = m_Dates[i].year - firstYear;
int index = dates[i].year - firstYear;
PixelValueType inValue = input[i];
if (inValue != this->GetInputNoDataValue())
......@@ -261,7 +257,7 @@ inline TOutputPixel operator ()(const TInputPixel& input) const
*
*/
template< class TInputPixel, class TOutputPixel>
class TSAmplitudeReduceFunctor : TSReduceFunctorBase<TInputPixel,TOutputPixel>
class TSAmplitudeReduceFunctor : public TSReduceFunctorBase<TInputPixel,TOutputPixel>
{
public:
......@@ -286,17 +282,17 @@ inline TOutputPixel operator ()(const TInputPixel& input) const
output.SetSize(this->GetOutputSize());
output.Fill(this->GetInputNoDataValue());
DatesType m_Dates = this->GetDates();
DatesType dates = this->GetDates();
TOutputPixel minimums, maximums;
minimums.SetSize(this->GetOutputSize());
maximums.SetSize(this->GetOutputSize());
minimums.Fill(this->GetInputNoDataValue());
maximums.Fill(this->GetInputNoDataValue());
int firstYear = m_Dates[0].year;
int firstYear = dates[0].year;
for (unsigned int i = 0; i < input.Size(); i++)
{
int index = m_Dates[i].year - firstYear;
int index = dates[i].year - firstYear;
PixelValueType inValue = input[i];
if (inValue != this->GetInputNoDataValue())
......@@ -324,10 +320,10 @@ inline TOutputPixel operator ()(const TInputPixel& input) const
}
// Compute amplitude
int lastYear = m_Dates[m_Dates.Size()-1].year;
int lastYear = dates[dates.size()-1].year;
for (unsigned int i = 0 ; i < lastYear - firstYear +1 ; i++ )
{
if (minimums[i] != m_InputNoDataValue && maximums[i] != m_InputNoDataValue)
if (minimums[i] != this->GetInputNoDataValue() && maximums[i] != this->GetInputNoDataValue())
{
output[i] = maximums[i] - minimums[i];
}
......@@ -421,7 +417,7 @@ public:
// Variance(t)
OutputPixelValueType var_t = isqSum + nfInv * isum2;
// Variance(x)outpixel[0]
// Variance(x)
OutputPixelValueType var_x = sqSum + nfInv * sum2;
// OLS
......@@ -458,11 +454,11 @@ public:
/* Set / Get input no-data value */
void SetInputNoDataValue(InputPixelValueType value) { m_InputNoDataValue = value; }
InputPixelValueType GetInputNoDataValue(){ return m_InputNoDataValue; }
InputPixelValueType GetInputNoDataValue() const { return m_InputNoDataValue; }
/* Set / Get output no-data value */
void SetOutputNoDataValue(OutputPixelValueType value) { m_OutputNoDataValue = value; }
InputPixelValueType GetOutputNoDataValue(){ return m_OutputNoDataValue; }
OutputPixelValueType GetOutputNoDataValue() const { return m_OutputNoDataValue; }
protected:
InputPixelValueType m_InputNoDataValue;
......@@ -476,11 +472,20 @@ protected:
*
* \brief Class for computing the label corresponding to an input slope and p-value
*
*/
* Output value | Definition
* --------------------------------------
* 0 | Very Strong Neg.
* 1 | Strong Neg.
* 2 | Moderate Neg.
* 3 | Very Strong Pos.
* 4 | Strong Pos.
* 5 | Moderate Pos.
* 6 | No Significant
*
**/
template< class TInputPixel, class TLabelPixel>
class SlopeAndPValueLabelingFunctor
{
public:
enum Def
{
......@@ -573,11 +578,11 @@ public:
/* Set / Get input no-data value */
void SetInputNoDataValue(InputPixelValueType value) { m_InputNoDataValue = value; }
InputPixelValueType GetInputNoDataValue(){ return m_InputNoDataValue; }
InputPixelValueType GetInputNoDataValue() const { return m_InputNoDataValue; }
/* Set / Get output no-data value */
void SetOutputNoDataValue(LabelPixelValueType value) { m_OutputNoDataValue = value; }
InputPixelValueType GetOutputNoDataValue(){ return m_OutputNoDataValue; }
LabelPixelValueType GetOutputNoDataValue() const { return m_OutputNoDataValue; }
protected:
InputPixelValueType m_InputNoDataValue;
......@@ -587,6 +592,109 @@ protected:
}; // SlopeAndPValueLabelingFunctor
/*
* \class RainfallEstimatedNDVIResiduesFunctor
* \brief Compute the residues between the estimation of the NDVI from the rainfall,
* and the real NDVI. Returns the residues.
*/
template< class TNDVIPixel, class TRainfallPixel>
class RainfallEstimatedNDVIResiduesFunctor
{
public:
typedef typename TNDVIPixel::ValueType NDVIPixelValueType;
typedef typename TRainfallPixel::ValueType RainfallPixelValueType;
RainfallEstimatedNDVIResiduesFunctor(){}
~RainfallEstimatedNDVIResiduesFunctor(){}
// Purposely not implemented
bool operator!=( const RainfallEstimatedNDVIResiduesFunctor & ) const { return false; }
bool operator==( const RainfallEstimatedNDVIResiduesFunctor & other ) const { return !(*this != other); }
/* Set / Get NDVI input no-data value */
void SetNDVINoDataValue(NDVIPixelValueType value) { m_NDVINoDataValue = value; }
NDVIPixelValueType GetNDVINoDataValue() const { return m_NDVINoDataValue; }
/* Set / Get Rainfall input no-data value */
void SetRainfallNoDataValue(RainfallPixelValueType value) { m_RainfallNoDataValue = value; }
RainfallPixelValueType GetRainfallNoDataValue() const { return m_RainfallNoDataValue; }
/*
* Compute output pixel.
* inputNDVIPixel and inputRainfallPixel must have the same number of components.
*/
inline TNDVIPixel operator ()(const TNDVIPixel& inputNDVIPixel, const TRainfallPixel& inputRainfallPixel) const
{
unsigned int n = inputNDVIPixel.Size();
// Prepare output pixel
TNDVIPixel residues;
residues.SetSize(n);
// 1. Compute regression between Y=NDVI and X=Rainfall
NDVIPixelValueType ndviSum = 0;
NDVIPixelValueType rainSum = 0;
NDVIPixelValueType ndviSqSum = 0;
NDVIPixelValueType coSum = 0;
unsigned int count = 0;
for (unsigned int i = 0 ; i < n ; i++)
{
if (inputNDVIPixel[i] != m_NDVINoDataValue && inputRainfallPixel[i] != m_RainfallNoDataValue)
{
NDVIPixelValueType ndviVal = inputNDVIPixel[i];
NDVIPixelValueType rainVal = static_cast<NDVIPixelValueType>(inputRainfallPixel[i]);
ndviSum += ndviVal;
rainSum += rainVal;
ndviSqSum += ndviVal*ndviVal;
coSum += ndviVal*rainVal;
count++;
}
}
if (count > 1)
{
NDVIPixelValueType nf = static_cast<NDVIPixelValueType>(count);
NDVIPixelValueType nfInv = -1.0 / nf;
NDVIPixelValueType ndviSum2 = ndviSum * ndviSum;
NDVIPixelValueType prodsums = ndviSum * rainSum;
// Covariance
NDVIPixelValueType cov = coSum + nfInv * prodsums;
// Variance(ndvi)
NDVIPixelValueType var_ndvi = ndviSqSum + nfInv * ndviSum2;
// OLS Y = slope * X + zeroy
NDVIPixelValueType slope = cov / var_ndvi;
NDVIPixelValueType zeroy = nfInv *( ndviSum - slope * rainSum );
// 2. Compute residues
for (unsigned int i = 0 ; i < n ; i++)
{
if (inputNDVIPixel[i] != m_NDVINoDataValue && inputRainfallPixel[i] != m_RainfallNoDataValue)
{
residues[i] = inputNDVIPixel[i] - (zeroy + slope * inputRainfallPixel[i]);
}
else
{
residues[i] = m_NDVINoDataValue;
}
}
}
else
{
residues.Fill(m_NDVINoDataValue);
}
return residues;
}
protected:
NDVIPixelValueType m_NDVINoDataValue;
RainfallPixelValueType m_RainfallNoDataValue;
}; // RainfallEstimatedNDVIResiduesFunctor
} // namespace functor
} // namespace otb
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment