27 #include "../common/MatrixUtils.h"
28 #include "../common/progress/TaskProgress.h"
29 #include "../geometry/Coord2D.h"
30 #include "../geometry/Envelope.h"
31 #include "../geometry/Point.h"
32 #include "../geometry/Polygon.h"
33 #include "../raster/Grid.h"
34 #include "../raster/PositionIterator.h"
35 #include "../raster/Utils.h"
47 #include <boost/numeric/ublas/io.hpp>
48 #include <boost/numeric/ublas/lu.hpp>
49 #include <boost/numeric/ublas/matrix.hpp>
79 m_acceptanceThreshold = 0;
101 m_covarianceMatrix(cm)
110 m_myCluster(rhs.m_myCluster),
128 m_meanVector.clear();
130 m_covarianceMatrix.clear();
132 m_covarianceInversion.clear();
139 assert(m_myCluster == 0);
142 double total_area = m_area + p->
m_area;
143 for (
unsigned i = 0; i < m_meanVector.size(); i++)
144 m_meanVector[i] = (m_meanVector[i] * m_area + p->
m_meanVector[i] * p->
m_area) / total_area;
165 unsigned int nBands = m_meanVector.size();
167 boost::numeric::ublas::matrix<double> term1(1, nBands);
168 boost::numeric::ublas::matrix<double> term2(nBands, 1);
170 for (
unsigned int i = 0; i < nBands; i++)
174 term2(i, 0) = term1(0, i);
177 term1 = prod(term1, m_covarianceInversion);
178 term1 = prod(term1, term2);
188 assert(m_myCluster == 0);
191 return m_id == rhs.
m_id;
223 const double dmax = std::numeric_limits<double>::max();
224 double chiTable[10][6] =
227 { 1.32, 2.71, 3.84, 6.64, 10.83, dmax},
228 { 2.77, 4.61, 5.99, 9.21, 13.82, dmax},
229 { 4.11, 6.25, 7.82, 11.35, 16.27, dmax},
230 { 5.39, 7.78, 9.49, 13.28, 18.47, dmax},
231 { 6.63, 9.24, 11.07, 15.09, 20.52, dmax},
232 { 7.84, 10.65, 12.59, 16.81, 22.46, dmax},
233 { 9.04, 12.02, 14.07, 18.48, 24.32, dmax},
234 {10.22, 13.36, 15.51, 20.09, 26.13, dmax},
235 {11.39, 14.68, 16.92, 21.67, 27.88, dmax},
236 {12.55, 15.99, 18.31, 23.21, 29.59, dmax}
243 if(acceptanceThreshold < 90.0)
245 else if(acceptanceThreshold < 95.0)
247 else if(acceptanceThreshold < 99.0)
249 else if(acceptanceThreshold < 99.9)
251 else if(acceptanceThreshold < 100.0)
256 return chiTable[nBands-1][it];
260 const std::vector<te::gm::Polygon*>& inputPolygons,
te::rst::Raster& outputRaster,
261 const unsigned int outputRasterBand,
const bool enableProgressInterface)
throw(te::rp::Exception)
270 for (
unsigned i = 0; i < inputPolygons.size(); i++)
274 std::vector<std::vector<double> > values_in_polygon = rattributes.getValuesFromRaster(inputRaster, *polygon, inputRasterBands);
275 std::vector<double> means;
276 for (
unsigned int b = 0; b < values_in_polygon.size(); b++)
279 means.push_back(summary.
m_mean);
282 Pattern* region =
new Pattern(i, polygon->
getArea(), means, rattributes.getCovarianceMatrix(values_in_polygon, means));
288 m_regions.insert(std::pair<double, Pattern*> (region->
m_area, region));
297 std::multimap<double, Pattern*, std::greater<double> >::iterator rit;
298 std::multimap<double, Pattern*, std::greater<double> >::iterator riti;
299 std::multimap<double, Pattern*, std::greater<double> >::iterator ritii;
304 std::set<std::pair<unsigned int, unsigned int> > compared;
318 for (riti = rit; riti !=
m_regions.end(); ++riti)
320 if (compared.find(std::pair<unsigned int, unsigned int>(rit->second->m_myCluster->m_id, riti->second->m_myCluster->m_id)) == compared.end())
322 compared.insert(std::pair<unsigned int, unsigned int>(rit->second->m_myCluster->m_id, riti->second->m_myCluster->m_id));
323 compared.insert(std::pair<unsigned int, unsigned int>(riti->second->m_myCluster->m_id, rit->second->m_myCluster->m_id));
328 if (rit->second->m_myCluster == riti->second->m_myCluster)
332 distance = rit->second->m_myCluster->getDistance(riti->second->m_myCluster);
334 distance2 = riti->second->m_myCluster->getDistance(rit->second->m_myCluster);
336 distance = (distance2 < distance? distance2: distance);
338 if (distance <= maxDistance)
343 oldid = riti->second->m_myCluster->m_id;
345 for (ritii = riti; ritii !=
m_regions.end(); ++ritii)
348 if (ritii->second->m_myCluster->m_id != oldid)
351 ritii->second->m_myCluster = rit->second->m_myCluster;
356 task_clustering.
pulse();
360 std::set<unsigned int> ulabels;
361 std::set<unsigned int>::iterator lit;
362 std::map<unsigned int, unsigned int> colormap;
365 ulabels.insert(rit->second->m_myCluster->m_id);
367 unsigned int color = 1;
368 for (lit = ulabels.begin(); lit != ulabels.end(); ++lit)
369 colormap[*lit] = color++;
372 unsigned int pattern;
380 pattern = rit->second->m_myCluster->m_id;
388 outputRaster.setValue(it.
getColumn(), it.
getRow(), colormap[pattern], outputRasterBand);
A structure to hold the set of statistics from a set of numerical values.
Raster ISOSeg Classifier strategy factory.
bool operator=(Pattern &rhs)
Return true if two clusters are equal.
AbstractParameters * clone() const
Create a clone copy of this instance.
double getDistance(Pattern *p)
Returns the Mahalanobis distance between two patterns.
~ClassifierISOSegStrategyFactory()
This class implements the strategy to iterate with spatial restriction, the iteration occurs inside a...
boost::numeric::ublas::matrix< double > m_covarianceInversion
The inversion of covariance matrix between bands.
Extraction of attributes from Raster, Bands, and Polygons.
This class can be used to inform the progress of a task.
Defines the exact sequence of characters that are acceptable.
void add(Pattern *p)
Add a region inside a cluster.
unsigned int getRow() const
Returns the current row in iterator.
An utility struct for representing 2D coordinates.
std::multimap< double, Pattern *, std::greater< double > > m_regions
A descriptive set of regions (area, features).
#define TE_TR(message)
It marks a string in order to get translated.
bool m_isInitialized
True if this instance is initialized.
bool initialize(StrategyParameters const *const strategyParams)
Initialize the classification strategy.
ClassifierISOSegStrategyFactory()
Extraction of attributes from Raster, Bands, and Polygons.
Raster Processing functions.
#define TERP_TRUE_OR_RETURN_FALSE(value, message)
Checks if value is true. For false values a warning message will be logged and a return of context wi...
Raster classifier strategy factory base class.
ISOSeg strategy for OBIA classification. The algorithm orders regions by area (larger first)...
const Parameters & operator=(const Parameters ¶ms)
static PolygonIterator end(const te::rst::Raster *r, const te::gm::Polygon *p)
Returns an iterator referring to after the end of the iterator.
An abstract class for raster data strucutures.
Raster strategy parameters base class.
ClassifierISOSegStrategy()
void pulse()
Calls setCurrentStep() function using getCurrentStep() + 1.
double getArea() const
It returns the area of this surface, as measured in the spatial reference system of this surface...
int m_id
The id of the region of the pattern.
te::rp::ClassifierStrategy * build()
Concrete factories (derived from this one) must implement this method in order to create objects...
Abstract parameters base interface.
Polygon is a subclass of CurvePolygon whose rings are defined by linear rings.
Describes a region or a cluster (group of regions with similar properties) to be used by ISOSeg metho...
~ClassifierISOSegStrategy()
static PolygonIterator begin(const te::rst::Raster *r, const te::gm::Polygon *p)
Returns an iterator referring to the first value of the band.
Raster classifier strategy base class.
bool GetInverseMatrix(const boost::numeric::ublas::matrix< T > &inputMatrix, boost::numeric::ublas::matrix< T > &outputMatrix)
Matrix inversion.
ClassifierISOSegStrategy::Parameters m_parameters
Internal execution parameters.
double m_area
The area of all regions inside a pattern.
std::vector< double > m_meanVector
The vector of mean values, 1 value per band;.
bool execute(const te::rst::Raster &inputRaster, const std::vector< unsigned int > &inputRasterBands, const std::vector< te::gm::Polygon * > &inputPolygons, te::rst::Raster &outputRaster, const unsigned int outputRasterBand, const bool enableProgressInterface)
Executes the classification strategy.
double m_acceptanceThreshold
The acceptance threshold (the closer to 100%, few clusters are created).
double getThreshold(double acceptanceThreshold, unsigned nBands)
ISOSeg strategy for segmentation-based classification.
unsigned int getColumn() const
Returns the current column in iterator.
Pattern * m_myCluster
The associated cluster of this pattern (optional).
void reset()
Clear all internal allocated resources and reset the parameters instance to its initial state...
boost::numeric::ublas::matrix< double > m_covarianceMatrix
The covariance matrix between bands.