diff --git a/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.h b/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.h
new file mode 100644
index 0000000000000000000000000000000000000000..b7564bd125ef89148d7a2dc73f4b8dade17d9950
--- /dev/null
+++ b/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.h
@@ -0,0 +1,209 @@
+/*=========================================================================
+
+ Program: ORFEO Toolbox
+ Language: C++
+ Date: $Date$
+ Version: $Revision$
+
+
+ Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
+ See OTBCopyright.txt for details.
+
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef __otbBandsStatisticsAttributesLabelMapFilter_h
+#define __otbBandsStatisticsAttributesLabelMapFilter_h
+
+#include "otbStatisticsAttributesLabelMapFilter.h"
+#include "otbMultiToMonoChannelExtractROI.h"
+
+namespace otb
+{
+namespace Functor
+{
+/** \class BandStatsAttributesLabelObjectFunctor
+* \brief Functor to compute multiple statistics attributes.
+*
+* For one label object, this functors applies the
+* StatisticsAttributesLabelObjectFunctor for each feature
+*
+* As such, it allows to compute in one pass statistics related to
+* multiple features. It is used in the
+* BandsStatisticsAttributesLabelMapFilter.
+*
+* Features can be added, removed or cleared via the appropriate
+* methods.
+*
+* \sa BandsStatisticsAttributesLabelMapFilter
+* \sa StatisticsAttributesLabelObjectFunctor
+*/
+template <class TLabelObject, class TFeatureImage>
+class BandStatsAttributesLabelObjectFunctor
+{
+public:
+ // Self typedef
+ typedef BandStatsAttributesLabelObjectFunctor Self;
+
+ /// Typedef of the feature image type
+ typedef typename TFeatureImage::PixelType FeatureType;
+
+ /// Typedef of the label object
+ typedef TLabelObject LabelObjectType;
+
+ /// Feature image const pointer
+ typedef typename TFeatureImage::ConstPointer FeatureImageConstPointer;
+
+ /// Statistics functor
+ typedef StatisticsAttributesLabelObjectFunctor
+ <TLabelObject, TFeatureImage> StatsFunctorType;
+
+ /// Map to store the functors
+ typedef std::map<std::string, StatsFunctorType> StatsFunctorsMapType;
+
+ /** Constructor */
+ BandStatsAttributesLabelObjectFunctor();
+
+ /** Destructor */
+ virtual ~BandStatsAttributesLabelObjectFunctor();
+
+ /** The comparators */
+ bool operator !=(const Self& self);
+ bool operator ==(const Self& self);
+
+ /** This is the functor implementation
+ * Calling the functor on a label object
+ * will update its statistics attributes */
+ inline void operator ()(LabelObjectType * lo) const;
+
+ /** Add a feature with the given name */
+ void AddFeature(const std::string& name, const TFeatureImage * img);
+
+ /** Remove the feature with this name if it exists */
+ bool RemoveFeature(const std::string& name);
+
+ /** Get the feature image with this name */
+ const TFeatureImage * GetFeatureImage(const std::string& name) const;
+
+ /** Clear all the features */
+ void ClearAllFeatures();
+
+ /** Get the number of features */
+ unsigned int GetNumberOfFeatures() const;
+
+ /** Set the reduced attribute set */
+ void SetReducedAttributeSet(bool flag);
+
+ /** Get the reduced attribute set */
+ bool GetReducedAttributeSet() const;
+
+private:
+ /// True to compute only a reduced attribute set
+ bool m_ReducedAttributeSet;
+
+ /// The Stat functors map
+ StatsFunctorsMapType m_StatsFunctorsMap;
+};
+} // End namespace Functor
+
+/** \class BandsStatisticsAttributesLabelMapFilter
+ * \brief This filter computes band statistics attributes for each object.
+ *
+ * Images are supposed to be compatible with otb::VectorImage
+ *
+ * This filter internally applies the
+ * StatisticsAttributesLabelMapFilter each channel independently
+ *
+ * The feature name is constructed as:
+ * 'STATS' + '::' + 'Band' + #BandIndex + '::' + StatisticName
+ *
+ * The ReducedAttributesSet flag allows to tell the internal
+ * statistics filter to compute only the main attributes (mean, variance, skewness and kurtosis).
+ *
+ * \sa MultiStatsAttributesLabelObjectFunctor AttributesMapLabelObject
+ *
+ * \ingroup ImageEnhancement MathematicalMorphologyImageFilters
+ */
+template<class TImage, class TFeatureImage>
+class ITK_EXPORT BandsStatisticsAttributesLabelMapFilter
+ : public LabelMapFeaturesFunctorImageFilter
+ <TImage,
+ typename Functor::BandStatsAttributesLabelObjectFunctor
+ <typename TImage::LabelObjectType, otb::Image<double, 2> > >
+{
+public:
+ /** Some convenient typedefs. */
+ typedef TImage ImageType;
+ typedef typename ImageType::LabelObjectType LabelObjectType;
+ typedef TFeatureImage FeatureImageType;
+ typedef typename FeatureImageType::InternalPixelType FeatureInternalPixelType;
+ typedef double InternalPrecisionType;
+ typedef Image<InternalPrecisionType, 2> InternalImageType;
+
+ /** Functor typedef */
+ typedef Functor::BandStatsAttributesLabelObjectFunctor
+ <LabelObjectType, InternalImageType> FunctorType;
+
+ /** Standard class typedefs. */
+ typedef BandsStatisticsAttributesLabelMapFilter Self;
+ typedef LabelMapFeaturesFunctorImageFilter
+ <ImageType, FunctorType> Superclass;
+ typedef itk::SmartPointer<Self> Pointer;
+ typedef itk::SmartPointer<const Self> ConstPointer;
+
+ /** ImageDimension constants */
+ itkStaticConstMacro(ImageDimension, unsigned int, TImage::ImageDimension);
+
+ /** Standard New method. */
+ itkNewMacro(Self);
+
+ /** Runtime information support. */
+ itkTypeMacro(BandsStatisticsAttributesLabelMapFilter, LabelMapFeaturesFunctorImageFilter);
+
+ // Channels
+ typedef MultiToMonoChannelExtractROI
+ <FeatureInternalPixelType, InternalPrecisionType> ChannelFilterType;
+ typedef typename ChannelFilterType::Pointer ChannelFilterPointerType;
+
+ /** Set the feature image */
+ void SetFeatureImage(const TFeatureImage *input);
+
+ /** Get the feature image */
+ const FeatureImageType * GetFeatureImage() const;
+
+ /** Set the reduced attribute set */
+ void SetReducedAttributeSet(bool flag);
+
+ /** Get the reduced attribute set */
+ bool GetReducedAttributeSet() const;
+
+ itkBooleanMacro(ReducedAttributeSet);
+
+protected:
+ /** Constructor */
+ BandsStatisticsAttributesLabelMapFilter();
+ /** Destructor */
+ ~BandsStatisticsAttributesLabelMapFilter() {}
+
+ /** Before threaded data generation */
+ virtual void BeforeThreadedGenerateData();
+
+ /** PrintSelf method */
+ void PrintSelf(std::ostream& os, itk::Indent indent) const;
+
+private:
+ BandsStatisticsAttributesLabelMapFilter(const Self &); //purposely not implemented
+ void operator =(const Self&); //purposely not implemented
+
+}; // end of class
+
+} // end namespace itk
+
+#ifndef OTB_MANUAL_INSTANTIATION
+#include "otbBandsStatisticsAttributesLabelMapFilter.txx"
+#endif
+
+#endif
diff --git a/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.txx b/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.txx
new file mode 100644
index 0000000000000000000000000000000000000000..224bbd7f904177924ddb1a121b9a885389309117
--- /dev/null
+++ b/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.txx
@@ -0,0 +1,249 @@
+/*=========================================================================
+
+ Program: ORFEO Toolbox
+ Language: C++
+ Date: $Date$
+ Version: $Revision$
+
+
+ Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
+ See OTBCopyright.txt for details.
+
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the above copyright notices for more information.
+
+=========================================================================*/
+#ifndef __otbBandsStatisticsAttributesLabelMapFilter_txx
+#define __otbBandsStatisticsAttributesLabelMapFilter_txx
+
+#include "otbBandsStatisticsAttributesLabelMapFilter.h"
+
+namespace otb
+{
+namespace Functor
+{
+/** Constructor */
+template <class TLabelObject, class TFeatureImage>
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::BandStatsAttributesLabelObjectFunctor() : m_ReducedAttributeSet(true),
+ m_StatsFunctorsMap()
+{}
+
+/** Destructor */
+template <class TLabelObject, class TFeatureImage>
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::~BandStatsAttributesLabelObjectFunctor(){}
+
+/** The comparators */
+template <class TLabelObject, class TFeatureImage>
+bool
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::operator != (const Self &self)
+ {
+ bool resp = true;
+ resp = resp && (m_ReducedAttributeSet != self.m_ReducedAttributeSet);
+ resp = resp && (m_StatsFunctorsMap != self.m_StatsFunctorsMap);
+ return resp;
+ }
+
+template <class TLabelObject, class TFeatureImage>
+bool
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::operator == (const Self &self)
+ {
+ return !(this != self);
+ }
+
+/** This is the functor implementation
+ * Calling the functor on a label object
+ * will update its statistics attributes */
+template <class TLabelObject, class TFeatureImage>
+void
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::operator() (LabelObjectType * lo) const
+{
+ // Walk every registered functors
+ for (typename StatsFunctorsMapType::const_iterator it = m_StatsFunctorsMap.begin();
+ it != m_StatsFunctorsMap.end(); ++it)
+ {
+ (it->second)(lo);
+ }
+}
+
+/** Add a feature with the given name */
+template <class TLabelObject, class TFeatureImage>
+void
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::AddFeature(const std::string& name, const TFeatureImage * img)
+{
+ // Create a new functor
+ StatsFunctorType newFunctor;
+
+ // Set the reduced attribute set option
+ newFunctor.SetReducedAttributeSet(m_ReducedAttributeSet);
+
+ // Set the feature and its name
+ newFunctor.SetFeatureName(name);
+
+ // Set the feature image
+ newFunctor.SetFeatureImage(img);
+
+ // Add it to the map
+ m_StatsFunctorsMap[name] = newFunctor;
+}
+
+/** Remove the feature with this name if it exists */
+template <class TLabelObject, class TFeatureImage>
+bool
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::RemoveFeature(const std::string& name)
+{
+ return (m_StatsFunctorsMap.erase(name) == 1);
+}
+
+/** Get the feature image with this name */
+template <class TLabelObject, class TFeatureImage>
+const TFeatureImage *
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::GetFeatureImage(const std::string& name) const
+{
+ if (m_StatsFunctorsMap.count(name) == 0)
+ {
+ itkGenericExceptionMacro(<< "No feature named " << name << " in map.");
+ }
+ return m_StatsFunctorsMap[name].GetFeatureImage();
+}
+
+/** Clear all the features */
+template <class TLabelObject, class TFeatureImage>
+void
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::ClearAllFeatures()
+{
+ m_StatsFunctorsMap.clear();
+}
+
+/** Get the number of features */
+template <class TLabelObject, class TFeatureImage>
+unsigned int
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::GetNumberOfFeatures() const
+{
+ return m_StatsFunctorsMap.size();
+}
+
+/** Set the reduced attribute set */
+template <class TLabelObject, class TFeatureImage>
+void
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::SetReducedAttributeSet(bool flag)
+{
+ // Set the flag
+ m_ReducedAttributeSet = flag;
+
+ // Set the flag to all the already existing functors
+ for (typename StatsFunctorsMapType::iterator it = m_StatsFunctorsMap.begin();
+ it != m_StatsFunctorsMap.end(); ++it)
+ {
+ it->second.SetReducedAttributeSet(m_ReducedAttributeSet);
+ }
+}
+/** Get the reduced attribute set */
+template <class TLabelObject, class TFeatureImage>
+bool
+BandStatsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
+::GetReducedAttributeSet() const
+{
+ return m_ReducedAttributeSet;
+}
+} // End namespace Functor
+
+template <class TImage, class TFeatureImage>
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::BandsStatisticsAttributesLabelMapFilter()
+{
+ this->SetNumberOfRequiredInputs(2);
+}
+
+/** Set the feature image */
+template <class TImage, class TFeatureImage>
+void
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::SetFeatureImage(const TFeatureImage *input)
+{
+ // Set the Nth input
+ this->SetNthInput(1, const_cast<TFeatureImage*>(input));
+}
+
+/** Get the feature image */
+template <class TImage, class TFeatureImage>
+const typename BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::FeatureImageType *
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::GetFeatureImage() const
+{
+ return static_cast<const TFeatureImage *>(this->itk::ProcessObject::GetInput(1));
+}
+
+/** Set the reduced attribute set */
+template <class TImage, class TFeatureImage>
+void
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::SetReducedAttributeSet(bool flag)
+{
+ if (this->GetFunctor().GetReducedAttributeSet() != flag)
+ {
+ this->GetFunctor().SetReducedAttributeSet(flag);
+ this->Modified();
+ }
+}
+
+/** Get the reduced attribute set */
+template <class TImage, class TFeatureImage>
+bool
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::GetReducedAttributeSet() const
+{
+ return this->GetFunctor().GetReducedAttributeSet();
+}
+
+template <class TImage, class TFeatureImage>
+void
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::BeforeThreadedGenerateData()
+{
+ // First call superclass implementation
+ Superclass::BeforeThreadedGenerateData();
+
+ unsigned long nbComponents = this->GetFeatureImage()->GetNumberOfComponentsPerPixel();
+
+ // Clear any previous feature
+ this->GetFunctor().ClearAllFeatures();
+
+ // Add each band of the feature image to the statistics functor
+ for (unsigned int i = 0; i < nbComponents; ++i)
+ {
+ ChannelFilterPointerType band = ChannelFilterType::New();
+ band->SetChannel(i + 1);
+ band->SetInput(this->GetFeatureImage());
+ band->GetOutput()->SetRequestedRegion(this->GetOutput()->GetRequestedRegion());
+ band->Update();
+ std::ostringstream oss;
+ oss << "Band" << i + 1; // [1..N] convention in feature naming
+ this->GetFunctor().AddFeature(oss.str(), band->GetOutput());
+ }
+
+}
+
+template <class TImage, class TFeatureImage>
+void
+BandsStatisticsAttributesLabelMapFilter<TImage, TFeatureImage>
+::PrintSelf(std::ostream& os, itk::Indent indent) const
+{
+ Superclass::PrintSelf(os, indent);
+}
+
+} // end namespace itk
+#endif
diff --git a/Code/OBIA/otbStatisticsAttributesLabelMapFilter.txx b/Code/OBIA/otbStatisticsAttributesLabelMapFilter.txx
index b3169c545a5799b0eb29f6c860eb8549f10ffce7..45d0973fe353856f0b4df6bdea3f26764834bd9c 100644
--- a/Code/OBIA/otbStatisticsAttributesLabelMapFilter.txx
+++ b/Code/OBIA/otbStatisticsAttributesLabelMapFilter.txx
@@ -78,15 +78,17 @@ StatisticsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
::operator() (LabelObjectType * lo) const
{
typename LabelObjectType::LineContainerType::const_iterator lit;
- typename LabelObjectType::LineContainerType& lineContainer = lo->GetLineContainer();
-
- FeatureType min = itk::NumericTraits<FeatureType>::max();
- FeatureType max = itk::NumericTraits<FeatureType>::NonpositiveMin();
- double sum = 0;
- double sum2 = 0;
- double sum3 = 0;
- double sum4 = 0;
- unsigned int totalFreq = 0;
+ typename LabelObjectType::LineContainerType& lineContainer = lo->GetLineContainer();
+
+ itk::OStringStream oss;
+
+ FeatureType min = itk::NumericTraits<FeatureType>::max();
+ FeatureType max = itk::NumericTraits<FeatureType>::NonpositiveMin();
+ double sum = 0;
+ double sum2 = 0;
+ double sum3 = 0;
+ double sum4 = 0;
+ unsigned int totalFreq = 0;
typename TFeatureImage::IndexType minIdx;
minIdx.Fill(0);
typename TFeatureImage::IndexType maxIdx;
@@ -104,7 +106,7 @@ StatisticsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
for (lit = lineContainer.begin(); lit != lineContainer.end(); lit++)
{
const typename TFeatureImage::IndexType& firstIdx = lit->GetIndex();
- unsigned long length = lit->GetLength();
+ unsigned long length = lit->GetLength();
long endIdx0 = firstIdx[0] + length;
for (typename TFeatureImage::IndexType idx = firstIdx; idx[0] < endIdx0; idx[0]++)
@@ -125,34 +127,38 @@ StatisticsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
}
//increase the sums
+ const double v2 = v * v;
+
sum += v;
- sum2 += vcl_pow((double) v, 2);
- sum3 += vcl_pow((double) v, 3);
- sum4 += vcl_pow((double) v, 4);
+ sum2 += v2;
+ sum3 += v2 * v;
+ sum4 += v2 * v2;
- // moments
- typename TFeatureImage::PointType physicalPosition;
- m_FeatureImage->TransformIndexToPhysicalPoint(idx, physicalPosition);
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ if (!m_ReducedAttributeSet)
{
- centerOfGravity[i] += physicalPosition[i] * v;
- centralMoments[i][i] += v * physicalPosition[i] * physicalPosition[i];
- for (unsigned int j = i + 1; j < TFeatureImage::ImageDimension; j++)
+ // moments
+ typename TFeatureImage::PointType physicalPosition;
+ m_FeatureImage->TransformIndexToPhysicalPoint(idx, physicalPosition);
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
{
- double weight = v * physicalPosition[i] * physicalPosition[j];
- centralMoments[i][j] += weight;
- centralMoments[j][i] += weight;
+ centerOfGravity[i] += physicalPosition[i] * v;
+ centralMoments[i][i] += v * physicalPosition[i] * physicalPosition[i];
+ for (unsigned int j = i + 1; j < TFeatureImage::ImageDimension; j++)
+ {
+ const double weight = v * physicalPosition[i] * physicalPosition[j];
+ centralMoments[i][j] += weight;
+ centralMoments[j][i] += weight;
+ }
}
}
-
}
}
// final computations
- double mean = sum / totalFreq;
- double variance = (sum2 - (vcl_pow(sum, 2) / totalFreq)) / (totalFreq - 1);
- double sigma = vcl_sqrt(variance);
- double mean2 = mean * mean;
+ const double mean = sum / totalFreq;
+ const double variance = (sum2 - (sum * sum / totalFreq)) / (totalFreq - 1);
+ const double sigma = vcl_sqrt(variance);
+ const double mean2 = mean * mean;
double skewness = 0;
double kurtosis = 0;
@@ -160,82 +166,10 @@ StatisticsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
if (vcl_abs(variance) > epsilon)
{
skewness = ((sum3 - 3.0 * mean * sum2) / totalFreq + 2.0 * mean * mean2) / (variance * sigma);
- kurtosis =
- ((sum4 - 4.0 * mean * sum3 + 6.0 * mean2 *
- sum2) / totalFreq - 3.0 * mean2 * mean2) / (variance * variance) - 3.0;
+ kurtosis = ((sum4 - 4.0 * mean * sum3 + 6.0 * mean2 * sum2) / totalFreq - 3.0 * mean2 * mean2) / (variance
+ * variance) - 3.0;
}
- double elongation = 0;
- if (sum != 0)
- {
- // Normalize using the total mass
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
- {
- centerOfGravity[i] /= sum;
- for (unsigned int j = 0; j < TFeatureImage::ImageDimension; j++)
- {
- centralMoments[i][j] /= sum;
- }
- }
-
- // Center the second order moments
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
- {
- for (unsigned int j = 0; j < TFeatureImage::ImageDimension; j++)
- {
- centralMoments[i][j] -= centerOfGravity[i] * centerOfGravity[j];
- }
- }
-
- // Compute principal moments and axes
- vnl_symmetric_eigensystem<double> eigen(centralMoments.GetVnlMatrix());
- vnl_diag_matrix<double> pm = eigen.D;
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
- {
- // principalMoments[i] = 4 * vcl_sqrt( pm(i,i) );
- principalMoments[i] = pm(i, i);
- }
- principalAxes = eigen.V.transpose();
-
- // Add a final reflection if needed for a proper rotation,
- // by multiplying the last row by the determinant
- vnl_real_eigensystem eigenrot(principalAxes.GetVnlMatrix());
- vnl_diag_matrix<vcl_complex<double> > eigenval = eigenrot.D;
- vcl_complex<double> det(1.0, 0.0);
-
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
- {
- det *= eigenval(i, i);
- }
-
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
- {
- principalAxes[TFeatureImage::ImageDimension - 1][i] *= std::real(det);
- }
-
- if (principalMoments[0] != 0)
- {
- // elongation = principalMoments[TFeatureImage::ImageDimension-1] / principalMoments[0];
- elongation = vcl_sqrt(principalMoments[TFeatureImage::ImageDimension - 1] / principalMoments[0]);
- }
- }
- else
- {
- // can't compute anything in that case - just set everything to a default value
- // Normalize using the total mass
- for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
- {
- centerOfGravity[i] = 0;
- principalMoments[i] = 0;
- for (unsigned int j = 0; j < TFeatureImage::ImageDimension; j++)
- {
- principalAxes[i][j] = 0;
- }
- }
- }
- itk::OStringStream oss;
-
- // finally put the values in the label object
oss.str("");
oss << "STATS::" << m_FeatureName << "::Mean";
lo->SetAttribute(oss.str().c_str(), mean);
@@ -250,55 +184,124 @@ StatisticsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
oss.str("");
oss << "STATS::" << m_FeatureName << "::Kurtosis";
- lo->SetAttribute(oss.str().c_str(), kurtosis);
+ lo->SetAttribute(oss.str().c_str(), kurtosis);
- // If we want all the features
if (!m_ReducedAttributeSet)
{
- oss.str("");
- oss << "STATS::" << m_FeatureName << "::Minimum";
- lo->SetAttribute(oss.str().c_str(), (double) min);
+ double elongation = 0;
+ if (sum != 0)
+ {
+ // Normalize using the total mass
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ {
+ centerOfGravity[i] /= sum;
+ for (unsigned int j = 0; j < TFeatureImage::ImageDimension; j++)
+ {
+ centralMoments[i][j] /= sum;
+ }
+ }
+
+ // Center the second order moments
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ {
+ for (unsigned int j = 0; j < TFeatureImage::ImageDimension; j++)
+ {
+ centralMoments[i][j] -= centerOfGravity[i] * centerOfGravity[j];
+ }
+ }
+
+ // Compute principal moments and axes
+ vnl_symmetric_eigensystem<double> eigen(centralMoments.GetVnlMatrix());
+ vnl_diag_matrix<double> pm = eigen.D;
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ {
+ // principalMoments[i] = 4 * vcl_sqrt( pm(i,i) );
+ principalMoments[i] = pm(i, i);
+ }
+ principalAxes = eigen.V.transpose();
- oss.str("");
- oss << "STATS::" << m_FeatureName << "::Maximum";
- lo->SetAttribute(oss.str().c_str(), (double) max);
+ // Add a final reflection if needed for a proper rotation,
+ // by multiplying the last row by the determinant
+ vnl_real_eigensystem eigenrot(principalAxes.GetVnlMatrix());
+ vnl_diag_matrix<vcl_complex<double> > eigenval = eigenrot.D;
+ vcl_complex<double> det(1.0, 0.0);
- oss.str("");
- oss << "STATS::" << m_FeatureName << "::Sum";
- lo->SetAttribute(oss.str().c_str(), sum);
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ {
+ det *= eigenval(i, i);
+ }
- oss.str("");
- oss << "STATS::" << m_FeatureName << "::Sigma";
- lo->SetAttribute(oss.str().c_str(), sigma);
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ {
+ principalAxes[TFeatureImage::ImageDimension - 1][i] *= std::real(det);
+ }
+
+ if (principalMoments[0] != 0)
+ {
+ // elongation = principalMoments[TFeatureImage::ImageDimension-1] / principalMoments[0];
+ elongation = vcl_sqrt(principalMoments[TFeatureImage::ImageDimension - 1] / principalMoments[0]);
+ }
+ }
+ else
+ {
+ // can't compute anything in that case - just set everything to a default value
+ // Normalize using the total mass
+ for (unsigned int i = 0; i < TFeatureImage::ImageDimension; i++)
+ {
+ centerOfGravity[i] = 0;
+ principalMoments[i] = 0;
+ for (unsigned int j = 0; j < TFeatureImage::ImageDimension; j++)
+ {
+ principalAxes[i][j] = 0;
+ }
+ }
+ }
- for (unsigned int dim = 0; dim < TFeatureImage::ImageDimension; ++dim)
- {
oss.str("");
- oss << "STATS::" << m_FeatureName << "::CenterOfGravity" << dim;
- lo->SetAttribute(oss.str().c_str(), centerOfGravity[dim]);
+ oss << "STATS::" << m_FeatureName << "::Minimum";
+ lo->SetAttribute(oss.str().c_str(), (double) min);
oss.str("");
- oss << "STATS::" << m_FeatureName << "::PrincipalMoments" << dim;
- lo->SetAttribute(oss.str().c_str(), principalMoments[dim]);
+ oss << "STATS::" << m_FeatureName << "::Maximum";
+ lo->SetAttribute(oss.str().c_str(), (double) max);
oss.str("");
- oss << "STATS::" << m_FeatureName << "::FirstMinimumIndex" << dim;
- lo->SetAttribute(oss.str().c_str(), minIdx[dim]);
+ oss << "STATS::" << m_FeatureName << "::Sum";
+ lo->SetAttribute(oss.str().c_str(), sum);
oss.str("");
- oss << "STATS::" << m_FeatureName << "::FirstMaximumIndex" << dim;
- lo->SetAttribute(oss.str().c_str(), maxIdx[dim]);
+ oss << "STATS::" << m_FeatureName << "::Sigma";
+ lo->SetAttribute(oss.str().c_str(), sigma);
- for (unsigned int dim2 = 0; dim2 < TFeatureImage::ImageDimension; ++dim2)
+ for (unsigned int dim = 0; dim < TFeatureImage::ImageDimension; ++dim)
{
oss.str("");
- oss << "STATS::" << m_FeatureName << "::PrincipalAxis" << dim << dim2;
- lo->SetAttribute(oss.str().c_str(), principalAxes(dim, dim2));
+ oss << "STATS::" << m_FeatureName << "::CenterOfGravity" << dim;
+ lo->SetAttribute(oss.str().c_str(), centerOfGravity[dim]);
+
+ oss.str("");
+ oss << "STATS::" << m_FeatureName << "::PrincipalMoments" << dim;
+ lo->SetAttribute(oss.str().c_str(), principalMoments[dim]);
+
+ oss.str("");
+ oss << "STATS::" << m_FeatureName << "::FirstMinimumIndex" << dim;
+ lo->SetAttribute(oss.str().c_str(), minIdx[dim]);
+
+ oss.str("");
+ oss << "STATS::" << m_FeatureName << "::FirstMaximumIndex" << dim;
+ lo->SetAttribute(oss.str().c_str(), maxIdx[dim]);
+
+ for (unsigned int dim2 = 0; dim2 < TFeatureImage::ImageDimension; ++dim2)
+ {
+ oss.str("");
+ oss << "STATS::" << m_FeatureName << "::PrincipalAxis" << dim << dim2;
+ lo->SetAttribute(oss.str().c_str(), principalAxes(dim, dim2));
+ }
}
}
- }
}
+
/** Set the name of the feature */
template <class TLabelObject, class TFeatureImage>
void
@@ -352,6 +355,7 @@ StatisticsAttributesLabelObjectFunctor<TLabelObject, TFeatureImage>
{
return m_ReducedAttributeSet;
}
+
} // End namespace Functor
template <class TImage, class TFeatureImage>
diff --git a/Testing/Code/OBIA/CMakeLists.txt b/Testing/Code/OBIA/CMakeLists.txt
index 46bd283112a8f1b0ce60ee9df1eaaf164bd9a1ba..8492d043e6b727a63223dc835d19d3eb9d722d6b 100644
--- a/Testing/Code/OBIA/CMakeLists.txt
+++ b/Testing/Code/OBIA/CMakeLists.txt
@@ -1,9 +1,9 @@
IF( NOT OTB_DISABLE_CXX_TESTING AND BUILD_TESTING )
-SET(BASELINE ${OTB_DATA_ROOT}/Baseline/OTB/Images)
+SET(BASELINE ${OTB_DATA_ROOT}/Baseline/OTB/Images)
SET(BASELINE_FILES ${OTB_DATA_ROOT}/Baseline/OTB/Files)
-SET(INPUTDATA ${OTB_DATA_ROOT}/Input)
-SET(TEMP ${OTBTesting_BINARY_DIR}/Temporary)
+SET(INPUTDATA ${OTB_DATA_ROOT}/Input)
+SET(TEMP ${OTBTesting_BINARY_DIR}/Temporary)
#Remote sensing images (large images )
IF(OTB_DATA_USE_LARGEINPUT)
@@ -39,21 +39,21 @@ ADD_TEST(obTuImageToLabelMapWithAttributesFilterNew ${OBIA_TESTS1}
ADD_TEST(obTvImageToLabelMapWithAttributesFilter ${OBIA_TESTS1}
otbImageToLabelMapWithAttributesFilter
- ${INPUTDATA}/calanques.tif
- ${INPUTDATA}/cala_labelled.tif)
+ ${INPUTDATA}/maur.tif
+ ${INPUTDATA}/maur_labelled.tif)
ADD_TEST(obTuLabelMapSourceNew ${OBIA_TESTS1}
- otbLabelMapSourceNew)
+ otbLabelMapSourceNew)
ADD_TEST(obTvLabelMapToVectorDataFilter ${OBIA_TESTS1}
otbLabelMapToVectorDataFilter
- ${INPUTDATA}/rcc8_mire1.png
- ${TEMP}/rcc8_mire2_vectorizer.shp)
+ ${INPUTDATA}/rcc8_mire1.png
+ ${TEMP}/rcc8_mire2_vectorizer.shp)
ADD_TEST(obTvVectorDataToLabelMapFilter ${OBIA_TESTS1}
otbVectorDataToLabelMapFilter
- ${INPUTDATA}/vectorIOexample_gis_to_vec.shp
- ${TEMP}/vectordataToLabelMap.png)
+ ${INPUTDATA}/vectorIOexample_gis_to_vec.shp
+ ${TEMP}/vectordataToLabelMap.png)
ADD_TEST(obTuLabelMapToSampleListFilterNew ${OBIA_TESTS1}
otbLabelMapToSampleListFilterNew)
@@ -66,14 +66,12 @@ ADD_TEST(obTvLabelMapToSampleListFilter ${OBIA_TESTS1}
SHAPE::Flusser09 SHAPE::Flusser10 SHAPE::Flusser11)
ADD_TEST(obTuLabelMapToVectorDataFilterNew ${OBIA_TESTS1}
-otbLabelMapToVectorDataFilterNew
-)
+ otbLabelMapToVectorDataFilterNew)
ADD_TEST(obTvLabelMapToVectorDataFilter ${OBIA_TESTS1}
-otbLabelMapToVectorDataFilter
- ${INPUTDATA}/rcc8_mire1.png
- ${TEMP}/rcc8_mire2_vectorizer.shp
-)
+ otbLabelMapToVectorDataFilter
+ ${INPUTDATA}/rcc8_mire1.png
+ ${TEMP}/rcc8_mire2_vectorizer.shp)
ADD_TEST(obTuLabelMapWithClassLabelToLabeledSampleListFilterNew ${OBIA_TESTS1}
otbLabelMapWithClassLabelToLabeledSampleListFilterNew)
@@ -86,47 +84,52 @@ ADD_TEST(obTvLabelMapWithClassLabelToLabeledSampleListFilter ${OBIA_TESTS1}
SHAPE::Flusser09 SHAPE::Flusser10 SHAPE::Flusser11)
ADD_TEST(otbLabelObjectMapVectorizer ${OBIA_TESTS1}
- otbLabelObjectMapVectorizer
- ${INPUTDATA}/rcc8_mire1.png
- rcc8_mire1_label_vectorizer.gml)
+ otbLabelObjectMapVectorizer
+ ${INPUTDATA}/rcc8_mire1.png
+ rcc8_mire1_label_vectorizer.gml)
ADD_TEST(obTuLabelObjectToPolygonFunctorNew ${OBIA_TESTS1}
-otbLabelObjectToPolygonFunctorNew
-)
+ otbLabelObjectToPolygonFunctorNew)
ADD_TEST(obTuMinMaxAttributesLabelMapFilterNew ${OBIA_TESTS1}
-otbMinMaxAttributesLabelMapFilterNew
-)
+ otbMinMaxAttributesLabelMapFilterNew)
ADD_TEST(obTvMinMaxAttributesLabelMapFilter ${OBIA_TESTS1}
- otbMinMaxAttributesLabelMapFilter
- ${INPUTDATA}/calanques.tif
- ${INPUTDATA}/cala_labelled.tif
- ${TEMP}/obTvMinMaxAttributesLabelMapFilter.txt)
+ otbMinMaxAttributesLabelMapFilter
+ ${INPUTDATA}/maur.tif
+ ${INPUTDATA}/maur_labelled.tif
+ ${TEMP}/obTvMinMaxAttributesLabelMapFilter.txt)
ADD_TEST(obTuNormalizeAttributesLabelMapFilterNew ${OBIA_TESTS1}
- otbNormalizeAttributesLabelMapFilterNew
-)
+ otbNormalizeAttributesLabelMapFilterNew)
ADD_TEST(obTvNormalizeAttributesLabelMapFilter ${OBIA_TESTS1}
- otbNormalizeAttributesLabelMapFilter
- ${INPUTDATA}/calanques.tif
- ${INPUTDATA}/cala_labelled.tif
- ${TEMP}/obTvNormalizeAttributesLabelMapFilter.txt)
+ otbNormalizeAttributesLabelMapFilter
+ ${INPUTDATA}/maur.tif
+ ${INPUTDATA}/maur_labelled.tif
+ ${TEMP}/obTvNormalizeAttributesLabelMapFilter.txt)
ADD_TEST(obTuKMeansAttributesLabelMapFilterNew ${OBIA_TESTS1}
- otbKMeansAttributesLabelMapFilterNew)
+ otbKMeansAttributesLabelMapFilterNew)
ADD_TEST(obTvKMeansAttributesLabelMapFilter ${OBIA_TESTS1}
- otbKMeansAttributesLabelMapFilter
- ${INPUTDATA}/calanques.tif
- ${INPUTDATA}/cala_labelled.tif
- ${TEMP}/obTvKMeansAttributesLabelMapFilter.txt)
+ otbKMeansAttributesLabelMapFilter
+ ${INPUTDATA}/maur.tif
+ ${INPUTDATA}/maur_labelled.tif
+ ${TEMP}/obTvKMeansAttributesLabelMapFilter.txt)
ADD_TEST(obTuRadiometricAttributesLabelMapFilterNew ${OBIA_TESTS1}
-otbRadiometricAttributesLabelMapFilterNew
-)
+ otbRadiometricAttributesLabelMapFilterNew)
+ADD_TEST(obTuBandsStatisticsAttributesLabelMapFilterNew ${OBIA_TESTS1}
+ otbBandsStatisticsAttributesLabelMapFilterNew)
+
+ADD_TEST(obTvBandsStatisticsAttributesLabelMapFilter ${OBIA_TESTS1}
+ otbBandsStatisticsAttributesLabelMapFilter
+ ${INPUTDATA}/maur.tif
+ ${INPUTDATA}/maur_labelled.tif
+ ${TEMP}/obTvBandsStatisticsAttributesLabelMapFilter.txt)
+
ADD_TEST(obTuShapeAttributesLabelMapFilterNew ${OBIA_TESTS1}
otbShapeAttributesLabelMapFilterNew)
@@ -134,26 +137,24 @@ ADD_TEST(obTuStatisticsAttributesLabelMapFilterNew ${OBIA_TESTS1}
otbStatisticsAttributesLabelMapFilterNew)
ADD_TEST(obTuVectorDataToLabelMapFilterNew ${OBIA_TESTS1}
-otbVectorDataToLabelMapFilterNew
-)
+ otbVectorDataToLabelMapFilterNew)
ADD_TEST(obTvVectorDataToLabelMapFilter ${OBIA_TESTS1}
-otbVectorDataToLabelMapFilter
- ${INPUTDATA}/vectorIOexample_gis_to_vec.shp
- ${TEMP}/vectordataToLabelMap.png
-)
+ otbVectorDataToLabelMapFilter
+ ${INPUTDATA}/vectorIOexample_gis_to_vec.shp
+ ${TEMP}/vectordataToLabelMap.png)
# OBIATests2 (need PQXX)
IF(OTB_USE_PQXX)
ADD_TEST(obTuLabelMapToGISTableFilterNew ${OBIA_TESTS2}
- otbLabelMapToGISTableFilterNew)
+ otbLabelMapToGISTableFilterNew)
ADD_TEST(obTuGISTableToLabelMapFilterNew ${OBIA_TESTS2}
- otbGISTableToLabelMapFilterNew)
+ otbGISTableToLabelMapFilterNew)
ADD_TEST(obTvLabelMapToGISTableFilter ${OBIA_TESTS2}
- otbLabelMapToGISTableFilter
+ otbLabelMapToGISTableFilter
${INPUTDATA}/rcc8_mire1.png
orfeotoolbox_test
labelmaptogis_test_table
@@ -192,6 +193,7 @@ otbMinMaxAttributesLabelMapFilter.cxx
otbNormalizeAttributesLabelMapFilter.cxx
otbKMeansAttributesLabelMapFilter.cxx
otbRadiometricAttributesLabelMapFilterNew.cxx
+otbBandsStatisticsAttributesLabelMapFilter.cxx
otbShapeAttributesLabelMapFilterNew.cxx
otbStatisticsAttributesLabelMapFilterNew.cxx
otbVectorDataToLabelMapFilter.cxx
diff --git a/Testing/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.cxx b/Testing/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..69a5037a1741adb609b24e926bdfa0408411f3e0
--- /dev/null
+++ b/Testing/Code/OBIA/otbBandsStatisticsAttributesLabelMapFilter.cxx
@@ -0,0 +1,112 @@
+/*=========================================================================
+
+ Program: ORFEO Toolbox
+ Language: C++
+ Date: $Date$
+ Version: $Revision$
+
+
+ Copyright (c) Centre National d'Etudes Spatiales. All rights reserved.
+ See OTBCopyright.txt for details.
+
+
+ This software is distributed WITHOUT ANY WARRANTY; without even
+ the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
+ PURPOSE. See the above copyright notices for more information.
+
+=========================================================================*/
+#if defined(_MSC_VER)
+#pragma warning ( disable : 4786 )
+#endif
+
+#include <fstream>
+#include <iostream>
+
+#include "otbImage.h"
+#include "otbVectorImage.h"
+#include "otbImageFileReader.h"
+#include "otbImageFileWriter.h"
+#include "otbAttributesMapLabelObject.h"
+#include "itkLabelImageToLabelMapFilter.h"
+#include "otbBandsStatisticsAttributesLabelMapFilter.h"
+
+const unsigned int Dimension = 2;
+typedef unsigned short LabelType;
+typedef double PixelType;
+
+typedef otb::AttributesMapLabelObject<LabelType, Dimension, double> LabelObjectType;
+typedef itk::LabelMap<LabelObjectType> LabelMapType;
+typedef otb::VectorImage<PixelType, Dimension> VectorImageType;
+typedef otb::Image<unsigned int,2> LabeledImageType;
+
+typedef LabelMapType::LabelObjectContainerType LabelObjectContainerType;
+typedef LabelObjectContainerType::const_iterator LabelObjectIterator;
+
+typedef otb::ImageFileReader<VectorImageType> ReaderType;
+typedef otb::ImageFileReader<LabeledImageType> LabeledReaderType;
+typedef otb::ImageFileWriter<VectorImageType> WriterType;
+typedef otb::ImageFileWriter<LabeledImageType> LabeledWriterType;
+
+typedef itk::LabelImageToLabelMapFilter<LabeledImageType,LabelMapType> LabelMapFilterType;
+typedef otb::BandsStatisticsAttributesLabelMapFilter<LabelMapType,VectorImageType> BandsStatisticsFilterType;
+
+
+int otbBandsStatisticsAttributesLabelMapFilterNew(int argc, char* argv[])
+{
+ BandsStatisticsFilterType::Pointer object = BandsStatisticsFilterType::New();
+ return EXIT_SUCCESS;
+}
+
+int otbBandsStatisticsAttributesLabelMapFilter(int argc, char* argv[])
+{
+ const char * infname = argv[1];
+ const char * lfname = argv[2];
+ const char * outfname = argv[3];
+
+ // Filters instanciation
+ ReaderType::Pointer reader = ReaderType::New();
+ LabeledReaderType::Pointer labeledReader = LabeledReaderType::New();
+ LabelMapFilterType::Pointer filter = LabelMapFilterType::New();
+ BandsStatisticsFilterType::Pointer stats = BandsStatisticsFilterType::New();
+
+ // Read inputs
+ reader->SetFileName(infname);
+ labeledReader->SetFileName(lfname);
+
+ // Make a LabelMap out of it
+ filter->SetInput(labeledReader->GetOutput());
+ filter->SetBackgroundValue(itk::NumericTraits<LabelType>::max());
+
+ //Compute band statistics attributes
+ stats->SetInput(filter->GetOutput());
+ stats->SetFeatureImage(reader->GetOutput());
+
+ stats->Update();
+
+ LabelMapType::Pointer labelMap = stats->GetOutput();
+
+ // Dump all results in the output file
+ std::ofstream outfile(outfname);
+ LabelObjectIterator it = labelMap->GetLabelObjectContainer().begin();
+ LabelObjectIterator end = labelMap->GetLabelObjectContainer().end();
+ for (; it != end; ++it)
+ {
+ LabelType label = it->first;
+ LabelObjectType::Pointer labelObject = it->second;
+
+ outfile << "Label " << label << " : " << std::endl;
+
+ std::vector<std::string> attributes = labelObject->GetAvailableAttributes();
+ std::vector<std::string>::const_iterator attrIt = attributes.begin();
+ std::vector<std::string>::const_iterator attrEnd = attributes.end();
+ for (; attrIt != attrEnd; ++attrIt)
+ {
+ LabelObjectType::AttributesValueType value = labelObject->GetAttribute(attrIt->c_str());
+ outfile << " " << *attrIt << " : "<< std::fixed << std::setprecision(6) << value << std::endl;
+ }
+ outfile << std::endl;
+ }
+
+ return EXIT_SUCCESS;
+}
+
diff --git a/Testing/Code/OBIA/otbLabelMapSVMClassifier.cxx b/Testing/Code/OBIA/otbLabelMapSVMClassifier.cxx
index 1c7ecafcec3535c200ce8f5f40b3dbefa9541557..66af9c2dc7022f008b9124f52427a8d21876ddc2 100644
--- a/Testing/Code/OBIA/otbLabelMapSVMClassifier.cxx
+++ b/Testing/Code/OBIA/otbLabelMapSVMClassifier.cxx
@@ -57,14 +57,16 @@ typedef itk::FixedArray<LabelType,1> TrainingVectorTyp
typedef itk::Statistics::ListSample<VectorType> ListSampleType;
typedef itk::Statistics::ListSample<TrainingVectorType> TrainingListSampleType;
typedef otb::LabelMapWithClassLabelToLabeledSampleListFilter<LabelMapType,ListSampleType,TrainingListSampleType>
- LabelMap2ListSampleFilterType;
+ ListSampleFilterType;
typedef otb::Functor::VariableLengthVectorToMeasurementVectorFunctor<VectorType> MeasurementVectorFunctorType;
typedef otb::SVMSampleListModelEstimator<ListSampleType,TrainingListSampleType,
MeasurementVectorFunctorType> SVMEstimatorType;
typedef otb::LabelMapSVMClassifier<LabelMapType> ClassifierType;
-typedef otb::LabelMapWithClassLabelToClassLabelImageFilter<LabelMapType,LabeledImageType> LabelMapWithClassLabelToClassLabelImageFilterType;
+typedef otb::LabelMapWithClassLabelToClassLabelImageFilter
+ <LabelMapType,LabeledImageType> ClassifImageGeneratorType;
+
LabelObjectType::Pointer makeTrainingSample(LabelMapType* labelMap, LabelType labelObjectId, LabelType classLabel)
{
@@ -86,33 +88,37 @@ int otbLabelMapSVMClassifier(int argc, char * argv[])
const char * lfname = argv[2];
const char * outfname = argv[3];
- // SmartPointer instanciation
- ReaderType::Pointer reader = ReaderType::New();
- LabeledReaderType::Pointer labeledReader = LabeledReaderType::New();
- LabelMapFilterType::Pointer filter = LabelMapFilterType::New();
- ShapeFilterType::Pointer shapeFilter = ShapeFilterType::New();
- RadiometricFilterType::Pointer radiometricFilter = RadiometricFilterType::New();
- ClassifierType::Pointer classifier = ClassifierType::New();
-
- // Inputs
+ // Filters instanciation
+ ReaderType::Pointer reader = ReaderType::New();
+ LabeledReaderType::Pointer labeledReader = LabeledReaderType::New();
+ LabelMapFilterType::Pointer filter = LabelMapFilterType::New();
+ ShapeFilterType::Pointer shapeFilter = ShapeFilterType::New();
+ RadiometricFilterType::Pointer radiometricFilter = RadiometricFilterType::New();
+ LabelMapType::Pointer trainingLabelMap = LabelMapType::New();
+ ListSampleFilterType::Pointer labelMap2SampleList = ListSampleFilterType::New();
+ SVMEstimatorType::Pointer svmEstim = SVMEstimatorType::New();
+ ClassifierType::Pointer classifier = ClassifierType::New();
+ ClassifImageGeneratorType::Pointer imGenerator = ClassifImageGeneratorType::New();
+ LabeledWriterType::Pointer writer = LabeledWriterType::New();
+
+ // Read inputs
reader->SetFileName(infname);
labeledReader->SetFileName(lfname);
- // Filter
+ // Make a LabelMap out of it
filter->SetInput(labeledReader->GetOutput());
filter->SetBackgroundValue(itk::NumericTraits<LabelType>::max());
+ //Compute shape and radimometric attributes
shapeFilter->SetInput(filter->GetOutput());
-
radiometricFilter->SetInput(shapeFilter->GetOutput());
radiometricFilter->SetFeatureImage(reader->GetOutput());
radiometricFilter->Update();
- // Build training samples
+ // Build a sub-LabelMap with class-labeled LabelObject
LabelMapType::Pointer labelMap = radiometricFilter->GetOutput();
- LabelMapType::Pointer trainingLabelMap = LabelMapType::New();
- // The following is specific to the input specified in CMakeLists
+ // The following is very specific to the input specified in CMakeLists
// water
trainingLabelMap->PushLabelObject(makeTrainingSample(labelMap, 13, 0));
// road
@@ -132,7 +138,7 @@ int otbLabelMapSVMClassifier(int argc, char * argv[])
trainingLabelMap->PushLabelObject(makeTrainingSample(labelMap, 161, 4));
trainingLabelMap->PushLabelObject(makeTrainingSample(labelMap, 46, 4));
- LabelMap2ListSampleFilterType::Pointer labelMap2SampleList = LabelMap2ListSampleFilterType::New();
+ // Make a ListSample out of trainingLabelMap
labelMap2SampleList->SetInputLabelMap(trainingLabelMap);
std::vector<std::string> attributes = labelMap->GetLabelObject(0)->GetAvailableAttributes();
@@ -144,13 +150,14 @@ int otbLabelMapSVMClassifier(int argc, char * argv[])
labelMap2SampleList->Update();
- SVMEstimatorType::Pointer svmEstim = SVMEstimatorType::New();
+ // Estimate SVM model
svmEstim->SetInputSampleList(labelMap2SampleList->GetOutputSampleList());
svmEstim->SetTrainingSampleList(labelMap2SampleList->GetOutputTrainingSampleList());
svmEstim->SetNumberOfClasses(5);
svmEstim->Modified();
svmEstim->Update();
+ // Classify using the whole LabelMap with estimated model
classifier->SetInput(labelMap);
classifier->SetModel(svmEstim->GetModel());
@@ -161,10 +168,9 @@ int otbLabelMapSVMClassifier(int argc, char * argv[])
classifier->Update();
- LabelMapWithClassLabelToClassLabelImageFilterType::Pointer imGenerator = LabelMapWithClassLabelToClassLabelImageFilterType::New();
+ // Make a labeled image with the classification result
imGenerator->SetInput(classifier->GetOutput());
- LabeledWriterType::Pointer writer = LabeledWriterType::New();
writer->SetInput(imGenerator->GetOutput());
writer->SetFileName(outfname);
writer->Update();
diff --git a/Testing/Code/OBIA/otbOBIATests1.cxx b/Testing/Code/OBIA/otbOBIATests1.cxx
index 993822ec7868a47ad06b4b2f007179e8c0ccabfd..7ecc330f96653b1f40bf5a2d5b516988b5aee1e3 100644
--- a/Testing/Code/OBIA/otbOBIATests1.cxx
+++ b/Testing/Code/OBIA/otbOBIATests1.cxx
@@ -48,6 +48,8 @@ REGISTER_TEST(otbNormalizeAttributesLabelMapFilter);
REGISTER_TEST(otbKMeansAttributesLabelMapFilterNew);
REGISTER_TEST(otbKMeansAttributesLabelMapFilter);
REGISTER_TEST(otbRadiometricAttributesLabelMapFilterNew);
+REGISTER_TEST(otbBandsStatisticsAttributesLabelMapFilterNew);
+REGISTER_TEST(otbBandsStatisticsAttributesLabelMapFilter);
REGISTER_TEST(otbShapeAttributesLabelMapFilterNew);
REGISTER_TEST(otbStatisticsAttributesLabelMapFilterNew);
REGISTER_TEST(otbVectorDataToLabelMapFilterNew);