--- /dev/null
+set(the_description "Shape descriptors and matchers.")
+ocv_define_module(shape opencv_core opencv_imgproc opencv_video)
--- /dev/null
+EMD-L1
+======
+Computes the "minimal work" distance between two weighted point configurations base on the papers "EMD-L1: An efficient and Robust Algorithm
+for comparing histogram-based descriptors", by Haibin Ling and Kazunori Okuda; and "The Earth Mover's Distance is the Mallows Distance:
+Some Insights from Statistics", by Elizaveta Levina and Peter Bickel.
+
+.. ocv:function:: float EMDL1( InputArray signature1, InputArray signature2 )
+
+ :param signature1: First signature, a single column floating-point matrix. Each row is the value of the histogram in each bin.
+
+ :param signature2: Second signature of the same format and size as ``signature1``.
--- /dev/null
+Cost Matrix for Histograms Common Interface
+===========================================
+
+.. highlight:: cpp
+
+A common interface is defined to ease the implementation of some algorithms pipelines, such
+as the Shape Context Matching Algorithm. A common class is defined, so any object that implements
+a Cost Matrix builder inherits the
+:ocv:class:`HistogramCostExtractor` interface.
+
+HistogramCostExtractor
+----------------------
+.. ocv:class:: HistogramCostExtractor : public Algorithm
+
+Abstract base class for histogram cost algorithms. ::
+
+ class CV_EXPORTS_W HistogramCostExtractor : public Algorithm
+ {
+ public:
+ CV_WRAP virtual void buildCostMatrix(InputArray descriptors1, InputArray descriptors2, OutputArray costMatrix) = 0;
+
+ CV_WRAP virtual void setNDummies(int nDummies) = 0;
+ CV_WRAP virtual int getNDummies() const = 0;
+
+ CV_WRAP virtual void setDefaultCost(float defaultCost) = 0;
+ CV_WRAP virtual float getDefaultCost() const = 0;
+ };
+
+NormHistogramCostExtractor
+--------------------------
+.. ocv:class:: NormHistogramCostExtractor : public HistogramCostExtractor
+
+A norm based cost extraction. ::
+
+ class CV_EXPORTS_W NormHistogramCostExtractor : public HistogramCostExtractor
+ {
+ public:
+ CV_WRAP virtual void setNormFlag(int flag) = 0;
+ CV_WRAP virtual int getNormFlag() const = 0;
+ };
+
+ CV_EXPORTS_W Ptr<HistogramCostExtractor>
+ createNormHistogramCostExtractor(int flag=cv::DIST_L2, int nDummies=25, float defaultCost=0.2);
+
+EMDHistogramCostExtractor
+-------------------------
+.. ocv:class:: EMDHistogramCostExtractor : public HistogramCostExtractor
+
+An EMD based cost extraction. ::
+
+ class CV_EXPORTS_W EMDHistogramCostExtractor : public HistogramCostExtractor
+ {
+ public:
+ CV_WRAP virtual void setNormFlag(int flag) = 0;
+ CV_WRAP virtual int getNormFlag() const = 0;
+ };
+
+ CV_EXPORTS_W Ptr<HistogramCostExtractor>
+ createEMDHistogramCostExtractor(int flag=cv::DIST_L2, int nDummies=25, float defaultCost=0.2);
+
+ChiHistogramCostExtractor
+-------------------------
+.. ocv:class:: ChiHistogramCostExtractor : public HistogramCostExtractor
+
+An Chi based cost extraction. ::
+
+ class CV_EXPORTS_W ChiHistogramCostExtractor : public HistogramCostExtractor
+ {};
+
+ CV_EXPORTS_W Ptr<HistogramCostExtractor> createChiHistogramCostExtractor(int nDummies=25, float defaultCost=0.2);
+
+EMDL1HistogramCostExtractor
+-------------------------
+.. ocv:class:: EMDL1HistogramCostExtractor : public HistogramCostExtractor
+
+An EMD-L1 based cost extraction. ::
+
+ class CV_EXPORTS_W EMDL1HistogramCostExtractor : public HistogramCostExtractor
+ {};
+
+ CV_EXPORTS_W Ptr<HistogramCostExtractor>
+ createEMDL1HistogramCostExtractor(int nDummies=25, float defaultCost=0.2);
--- /dev/null
+**********************************
+shape. Shape Distance and Matching
+**********************************
+
+The module contains algorithms that embed a notion of shape distance.
+These algorithms may be used for shape matching and retrieval, or shape
+comparison.
+
+.. toctree::
+ :maxdepth: 2
+
+ shape_distances
+ shape_transformers
+ histogram_cost_matrix
+ emdL1
--- /dev/null
+*****
+shape
+*****
+
+The module contains algorithms that embed a notion of shape distance.
+These algorithms may be used for shape matching and retrieval, or shape
+comparison.
+
+.. toctree::
+ :maxdepth: 2
+
+ shape_distances
+ shape_transformers
+ histogram_cost_matrix
+ emdL1
--- /dev/null
+Shape Distance and Common Interfaces
+====================================
+
+.. highlight:: cpp
+
+Shape Distance algorithms in OpenCV are derivated from a common interface that allows you to
+switch between them in a practical way for solving the same problem with different methods.
+Thus, all objects that implement shape distance measures inherit the
+:ocv:class:`ShapeDistanceExtractor` interface.
+
+
+ShapeDistanceExtractor
+----------------------
+.. ocv:class:: ShapeDistanceExtractor : public Algorithm
+
+Abstract base class for shape distance algorithms. ::
+
+ class CV_EXPORTS_W ShapeDistanceExtractor : public Algorithm
+ {
+ public:
+ CV_WRAP virtual float computeDistance(InputArray contour1, InputArray contour2) = 0;
+ };
+
+ShapeDistanceExtractor::computeDistance
+---------------------------------------
+Compute the shape distance between two shapes defined by its contours.
+
+.. ocv:function:: float ShapeDistanceExtractor::computeDistance( InputArray contour1, InputArray contour2 )
+
+ :param contour1: Contour defining first shape.
+
+ :param contour2: Contour defining second shape.
+
+ShapeContextDistanceExtractor
+-----------------------------
+.. ocv:class:: ShapeContextDistanceExtractor : public ShapeDistanceExtractor
+
+Implementation of the Shape Context descriptor and matching algorithm proposed by Belongie et al. in
+"Shape Matching and Object Recognition Using Shape Contexts" (PAMI 2002).
+This implementation is packaged in a generic scheme, in order to allow you the implementation of the
+common variations of the original pipeline. ::
+
+ class CV_EXPORTS_W ShapeContextDistanceExtractor : public ShapeDistanceExtractor
+ {
+ public:
+ CV_WRAP virtual void setAngularBins(int nAngularBins) = 0;
+ CV_WRAP virtual int getAngularBins() const = 0;
+
+ CV_WRAP virtual void setRadialBins(int nRadialBins) = 0;
+ CV_WRAP virtual int getRadialBins() const = 0;
+
+ CV_WRAP virtual void setInnerRadius(float innerRadius) = 0;
+ CV_WRAP virtual float getInnerRadius() const = 0;
+
+ CV_WRAP virtual void setOuterRadius(float outerRadius) = 0;
+ CV_WRAP virtual float getOuterRadius() const = 0;
+
+ CV_WRAP virtual void setRotationInvariant(bool rotationInvariant) = 0;
+ CV_WRAP virtual bool getRotationInvariant() const = 0;
+
+ CV_WRAP virtual void setShapeContextWeight(float shapeContextWeight) = 0;
+ CV_WRAP virtual float getShapeContextWeight() const = 0;
+
+ CV_WRAP virtual void setImageAppearanceWeight(float imageAppearanceWeight) = 0;
+ CV_WRAP virtual float getImageAppearanceWeight() const = 0;
+
+ CV_WRAP virtual void setBendingEnergyWeight(float bendingEnergyWeight) = 0;
+ CV_WRAP virtual float getBendingEnergyWeight() const = 0;
+
+ CV_WRAP virtual void setImages(InputArray image1, InputArray image2) = 0;
+ CV_WRAP virtual void getImages(OutputArray image1, OutputArray image2) const = 0;
+
+ CV_WRAP virtual void setIterations(int iterations) = 0;
+ CV_WRAP virtual int getIterations() const = 0;
+
+ CV_WRAP virtual void setCostExtractor(Ptr<HistogramCostExtractor> comparer) = 0;
+ CV_WRAP virtual Ptr<HistogramCostExtractor> getCostExtractor() const = 0;
+
+ CV_WRAP virtual void setTransformAlgorithm(Ptr<ShapeTransformer> transformer) = 0;
+ CV_WRAP virtual Ptr<ShapeTransformer> getTransformAlgorithm() const = 0;
+ };
+
+ /* Complete constructor */
+ CV_EXPORTS_W Ptr<ShapeContextDistanceExtractor>
+ createShapeContextDistanceExtractor(int nAngularBins=12, int nRadialBins=4,
+ float innerRadius=0.2, float outerRadius=2, int iterations=3,
+ const Ptr<HistogramCostExtractor> &comparer = createChiHistogramCostExtractor(),
+ const Ptr<ShapeTransformer> &transformer = createThinPlateSplineShapeTransformer());
+
+ShapeContextDistanceExtractor::setAngularBins
+---------------------------------------------
+Establish the number of angular bins for the Shape Context Descriptor used in the shape matching pipeline.
+
+.. ocv:function:: void setAngularBins( int nAngularBins )
+
+ :param nAngularBins: The number of angular bins in the shape context descriptor.
+
+ShapeContextDistanceExtractor::setRadialBins
+--------------------------------------------
+Establish the number of radial bins for the Shape Context Descriptor used in the shape matching pipeline.
+
+.. ocv:function:: void setRadialBins( int nRadialBins )
+
+ :param nRadialBins: The number of radial bins in the shape context descriptor.
+
+ShapeContextDistanceExtractor::setInnerRadius
+---------------------------------------------
+Set the inner radius of the shape context descriptor.
+
+.. ocv:function:: void setInnerRadius(float innerRadius)
+
+ :param innerRadius: The value of the inner radius.
+
+ShapeContextDistanceExtractor::setOuterRadius
+---------------------------------------------
+Set the outer radius of the shape context descriptor.
+
+.. ocv:function:: void setOuterRadius(float outerRadius)
+
+ :param outerRadius: The value of the outer radius.
+
+ShapeContextDistanceExtractor::setShapeContextWeight
+----------------------------------------------------
+Set the weight of the shape context distance in the final value of the shape distance.
+The shape context distance between two shapes is defined as the symmetric sum of shape
+context matching costs over best matching points.
+The final value of the shape distance is a user-defined linear combination of the shape
+context distance, an image appearance distance, and a bending energy.
+
+.. ocv:function:: void setShapeContextWeight( float shapeContextWeight )
+
+ :param shapeContextWeight: The weight of the shape context distance in the final distance value.
+
+ShapeContextDistanceExtractor::setImageAppearanceWeight
+-------------------------------------------------------
+Set the weight of the Image Appearance cost in the final value of the shape distance.
+The image appearance cost is defined as the sum of squared brightness differences in
+Gaussian windows around corresponding image points.
+The final value of the shape distance is a user-defined linear combination of the shape
+context distance, an image appearance distance, and a bending energy.
+If this value is set to a number different from 0, is mandatory to set the images that
+correspond to each shape.
+
+.. ocv:function:: void setImageAppearanceWeight( float imageAppearanceWeight )
+
+ :param imageAppearanceWeight: The weight of the appearance cost in the final distance value.
+
+ShapeContextDistanceExtractor::setBendingEnergyWeight
+-----------------------------------------------------
+Set the weight of the Bending Energy in the final value of the shape distance.
+The bending energy definition depends on what transformation is being used to align the
+shapes.
+The final value of the shape distance is a user-defined linear combination of the shape
+context distance, an image appearance distance, and a bending energy.
+
+.. ocv:function:: void setBendingEnergyWeight( float bendingEnergyWeight )
+
+ :param bendingEnergyWeight: The weight of the Bending Energy in the final distance value.
+
+ShapeContextDistanceExtractor::setImages
+----------------------------------------
+Set the images that correspond to each shape. This images are used in the calculation of the
+Image Appearance cost.
+
+.. ocv:function:: void setImages( InputArray image1, InputArray image2 )
+
+ :param image1: Image corresponding to the shape defined by ``contours1``.
+
+ :param image2: Image corresponding to the shape defined by ``contours2``.
+
+ShapeContextDistanceExtractor::setCostExtractor
+-----------------------------------------------
+Set the algorithm used for building the shape context descriptor cost matrix.
+
+.. ocv:function:: void setCostExtractor( Ptr<HistogramCostExtractor> comparer )
+
+ :param comparer: Smart pointer to a HistogramCostExtractor, an algorithm that defines the cost matrix between descriptors.
+
+ShapeContextDistanceExtractor::setStdDev
+----------------------------------------
+Set the value of the standard deviation for the Gaussian window for the image appearance cost.
+
+.. ocv:function:: void setStdDev( float sigma )
+
+ :param sigma: Standard Deviation.
+
+ShapeContextDistanceExtractor::setTransformAlgorithm
+----------------------------------------------------
+Set the algorithm used for aligning the shapes.
+
+.. ocv:function:: void setTransformAlgorithm( Ptr<ShapeTransformer> transformer )
+
+ :param comparer: Smart pointer to a ShapeTransformer, an algorithm that defines the aligning transformation.
+
+HausdorffDistanceExtractor
+--------------------------
+.. ocv:class:: HausdorffDistanceExtractor : public ShapeDistanceExtractor
+
+A simple Hausdorff distance measure between shapes defined by contours,
+according to the paper "Comparing Images using the Hausdorff distance." by
+D.P. Huttenlocher, G.A. Klanderman, and W.J. Rucklidge. (PAMI 1993). ::
+
+ class CV_EXPORTS_W HausdorffDistanceExtractor : public ShapeDistanceExtractor
+ {
+ public:
+ CV_WRAP virtual void setDistanceFlag(int distanceFlag) = 0;
+ CV_WRAP virtual int getDistanceFlag() const = 0;
+
+ CV_WRAP virtual void setRankProportion(float rankProportion) = 0;
+ CV_WRAP virtual float getRankProportion() const = 0;
+ };
+
+ /* Constructor */
+ CV_EXPORTS_W Ptr<HausdorffDistanceExtractor> createHausdorffDistanceExtractor(int distanceFlag=cv::NORM_L2, float rankProp=0.6);
+
+HausdorffDistanceExtractor::setDistanceFlag
+-------------------------------------------
+Set the norm used to compute the Hausdorff value between two shapes. It can be L1 or L2 norm.
+
+.. ocv:function:: void setDistanceFlag( int distanceFlag )
+
+ :param distanceFlag: Flag indicating which norm is used to compute the Hausdorff distance (NORM_L1, NORM_L2).
+
+HausdorffDistanceExtractor::setRankProportion
+---------------------------------------------
+This method sets the rank proportion (or fractional value) that establish the Kth ranked value of the
+partial Hausdorff distance. Experimentally had been shown that 0.6 is a good value to compare shapes.
+
+.. ocv:function:: void setRankProportion( float rankProportion )
+
+ :param rankProportion: fractional value (between 0 and 1).
--- /dev/null
+Shape Transformers and Interfaces
+=================================
+
+.. highlight:: cpp
+
+A virtual interface that ease the use of transforming algorithms in some pipelines, such as
+the Shape Context Matching Algorithm. Thus, all objects that implement shape transformation
+techniques inherit the
+:ocv:class:`ShapeTransformer` interface.
+
+ShapeTransformer
+----------------
+.. ocv:class:: ShapeTransformer : public Algorithm
+
+Abstract base class for shape transformation algorithms. ::
+
+ class CV_EXPORTS_W ShapeTransformer : public Algorithm
+ {
+ public:
+ CV_WRAP virtual void estimateTransformation(InputArray transformingShape, InputArray targetShape,
+ std::vector<DMatch>& matches) = 0;
+
+ CV_WRAP virtual float applyTransformation(InputArray input, OutputArray output=noArray()) = 0;
+
+ CV_WRAP virtual void warpImage(InputArray transformingImage, OutputArray output,
+ int flags=INTER_LINEAR, int borderMode=BORDER_CONSTANT,
+ const Scalar& borderValue=Scalar()) const = 0;
+ };
+
+ShapeTransformer::estimateTransformation
+----------------------------------------
+Estimate the transformation parameters of the current transformer algorithm, based on point matches.
+
+.. ocv:function:: void estimateTransformation( InputArray transformingShape, InputArray targetShape, std::vector<DMatch>& matches )
+
+ :param transformingShape: Contour defining first shape.
+
+ :param targetShape: Contour defining second shape (Target).
+
+ :param matches: Standard vector of Matches between points.
+
+ShapeTransformer::applyTransformation
+-------------------------------------
+Apply a transformation, given a pre-estimated transformation parameters.
+
+.. ocv:function:: float applyTransformation( InputArray input, OutputArray output=noArray() )
+
+ :param input: Contour (set of points) to apply the transformation.
+
+ :param output: Output contour.
+
+ShapeTransformer::warpImage
+---------------------------
+Apply a transformation, given a pre-estimated transformation parameters, to an Image.
+
+.. ocv:function:: void warpImage( InputArray transformingImage, OutputArray output, int flags=INTER_LINEAR, int borderMode=BORDER_CONSTANT, const Scalar& borderValue=Scalar() )
+
+ :param transformingImage: Input image.
+
+ :param output: Output image.
+
+ :param flags: Image interpolation method.
+
+ :param borderMode: border style.
+
+ :param borderValue: border value.
+
+ThinPlateSplineShapeTransformer
+-------------------------------
+.. ocv:class:: ThinPlateSplineShapeTransformer : public Algorithm
+
+Definition of the transformation ocupied in the paper "Principal Warps: Thin-Plate Splines and Decomposition
+of Deformations", by F.L. Bookstein (PAMI 1989). ::
+
+ class CV_EXPORTS_W ThinPlateSplineShapeTransformer : public ShapeTransformer
+ {
+ public:
+ CV_WRAP virtual void setRegularizationParameter(double beta) = 0;
+ CV_WRAP virtual double getRegularizationParameter() const = 0;
+ };
+
+ /* Complete constructor */
+ CV_EXPORTS_W Ptr<ThinPlateSplineShapeTransformer>
+ createThinPlateSplineShapeTransformer(double regularizationParameter=0);
+
+ThinPlateSplineShapeTransformer::setRegularizationParameter
+-----------------------------------------------------------
+Set the regularization parameter for relaxing the exact interpolation requirements of the TPS algorithm.
+
+.. ocv:function:: void setRegularizationParameter( double beta )
+
+ :param beta: value of the regularization parameter.
+
+AffineTransformer
+-----------------
+.. ocv:class:: AffineTransformer : public Algorithm
+
+Wrapper class for the OpenCV Affine Transformation algorithm. ::
+
+ class CV_EXPORTS_W AffineTransformer : public ShapeTransformer
+ {
+ public:
+ CV_WRAP virtual void setFullAffine(bool fullAffine) = 0;
+ CV_WRAP virtual bool getFullAffine() const = 0;
+ };
+
+ /* Complete constructor */
+ CV_EXPORTS_W Ptr<AffineTransformer> createAffineTransformer(bool fullAffine);
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009-2012, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifndef __OPENCV_SHAPE_HPP__
+#define __OPENCV_SHAPE_HPP__
+
+#include "opencv2/shape/emdL1.hpp"
+#include "opencv2/shape/shape_transformer.hpp"
+#include "opencv2/shape/hist_cost.hpp"
+#include "opencv2/shape/shape_distance.hpp"
+
+namespace cv
+{
+CV_EXPORTS bool initModule_shape();
+}
+
+#endif
+
+/* End of file. */
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009-2012, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifndef __OPENCV_EMD_L1_HPP__
+#define __OPENCV_EMD_L1_HPP__
+
+#include "opencv2/core.hpp"
+
+namespace cv
+{
+/****************************************************************************************\
+* EMDL1 Function *
+\****************************************************************************************/
+
+CV_EXPORTS float EMDL1(InputArray signature1, InputArray signature2);
+
+}//namespace cv
+
+#endif
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Copyright (C) 2013, OpenCV Foundation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifndef __OPENCV_HIST_COST_HPP__
+#define __OPENCV_HIST_COST_HPP__
+
+#include "opencv2/imgproc.hpp"
+
+namespace cv
+{
+
+/*!
+ * The base class for HistogramCostExtractor.
+ */
+class CV_EXPORTS_W HistogramCostExtractor : public Algorithm
+{
+public:
+ CV_WRAP virtual void buildCostMatrix(InputArray descriptors1, InputArray descriptors2, OutputArray costMatrix) = 0;
+
+ CV_WRAP virtual void setNDummies(int nDummies) = 0;
+ CV_WRAP virtual int getNDummies() const = 0;
+
+ CV_WRAP virtual void setDefaultCost(float defaultCost) = 0;
+ CV_WRAP virtual float getDefaultCost() const = 0;
+};
+
+/*! */
+class CV_EXPORTS_W NormHistogramCostExtractor : public HistogramCostExtractor
+{
+public:
+ CV_WRAP virtual void setNormFlag(int flag) = 0;
+ CV_WRAP virtual int getNormFlag() const = 0;
+};
+
+CV_EXPORTS_W Ptr<HistogramCostExtractor>
+ createNormHistogramCostExtractor(int flag=DIST_L2, int nDummies=25, float defaultCost=0.2);
+
+/*! */
+class CV_EXPORTS_W EMDHistogramCostExtractor : public HistogramCostExtractor
+{
+public:
+ CV_WRAP virtual void setNormFlag(int flag) = 0;
+ CV_WRAP virtual int getNormFlag() const = 0;
+};
+
+CV_EXPORTS_W Ptr<HistogramCostExtractor>
+ createEMDHistogramCostExtractor(int flag=DIST_L2, int nDummies=25, float defaultCost=0.2);
+
+/*! */
+class CV_EXPORTS_W ChiHistogramCostExtractor : public HistogramCostExtractor
+{};
+
+CV_EXPORTS_W Ptr<HistogramCostExtractor> createChiHistogramCostExtractor(int nDummies=25, float defaultCost=0.2);
+
+/*! */
+class CV_EXPORTS_W EMDL1HistogramCostExtractor : public HistogramCostExtractor
+{};
+
+CV_EXPORTS_W Ptr<HistogramCostExtractor>
+ createEMDL1HistogramCostExtractor(int nDummies=25, float defaultCost=0.2);
+
+} // cv
+#endif
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Copyright (C) 2013, OpenCV Foundation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifdef __OPENCV_BUILD
+#error this is a compatibility header which should not be used inside the OpenCV library
+#endif
+
+#include "opencv2/shape.hpp"
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Copyright (C) 2013, OpenCV Foundation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifndef __OPENCV_SHAPE_SHAPE_DISTANCE_HPP__
+#define __OPENCV_SHAPE_SHAPE_DISTANCE_HPP__
+#include "opencv2/core.hpp"
+#include "opencv2/shape/hist_cost.hpp"
+#include "opencv2/shape/shape_transformer.hpp"
+
+namespace cv
+{
+
+/*!
+ * The base class for ShapeDistanceExtractor.
+ * This is just to define the common interface for
+ * shape comparisson techniques.
+ */
+class CV_EXPORTS_W ShapeDistanceExtractor : public Algorithm
+{
+public:
+ CV_WRAP virtual float computeDistance(InputArray contour1, InputArray contour2) = 0;
+};
+
+/***********************************************************************************/
+/***********************************************************************************/
+/***********************************************************************************/
+/*!
+ * Shape Context implementation.
+ * The SCD class implements SCD algorithm proposed by Belongie et al.in
+ * "Shape Matching and Object Recognition Using Shape Contexts".
+ * Implemented by Juan M. Perez for the GSOC 2013.
+ */
+class CV_EXPORTS_W ShapeContextDistanceExtractor : public ShapeDistanceExtractor
+{
+public:
+ CV_WRAP virtual void setAngularBins(int nAngularBins) = 0;
+ CV_WRAP virtual int getAngularBins() const = 0;
+
+ CV_WRAP virtual void setRadialBins(int nRadialBins) = 0;
+ CV_WRAP virtual int getRadialBins() const = 0;
+
+ CV_WRAP virtual void setInnerRadius(float innerRadius) = 0;
+ CV_WRAP virtual float getInnerRadius() const = 0;
+
+ CV_WRAP virtual void setOuterRadius(float outerRadius) = 0;
+ CV_WRAP virtual float getOuterRadius() const = 0;
+
+ CV_WRAP virtual void setRotationInvariant(bool rotationInvariant) = 0;
+ CV_WRAP virtual bool getRotationInvariant() const = 0;
+
+ CV_WRAP virtual void setShapeContextWeight(float shapeContextWeight) = 0;
+ CV_WRAP virtual float getShapeContextWeight() const = 0;
+
+ CV_WRAP virtual void setImageAppearanceWeight(float imageAppearanceWeight) = 0;
+ CV_WRAP virtual float getImageAppearanceWeight() const = 0;
+
+ CV_WRAP virtual void setBendingEnergyWeight(float bendingEnergyWeight) = 0;
+ CV_WRAP virtual float getBendingEnergyWeight() const = 0;
+
+ CV_WRAP virtual void setImages(InputArray image1, InputArray image2) = 0;
+ CV_WRAP virtual void getImages(OutputArray image1, OutputArray image2) const = 0;
+
+ CV_WRAP virtual void setIterations(int iterations) = 0;
+ CV_WRAP virtual int getIterations() const = 0;
+
+ CV_WRAP virtual void setCostExtractor(Ptr<HistogramCostExtractor> comparer) = 0;
+ CV_WRAP virtual Ptr<HistogramCostExtractor> getCostExtractor() const = 0;
+
+ CV_WRAP virtual void setStdDev(float sigma) = 0;
+ CV_WRAP virtual float getStdDev() const = 0;
+
+ CV_WRAP virtual void setTransformAlgorithm(Ptr<ShapeTransformer> transformer) = 0;
+ CV_WRAP virtual Ptr<ShapeTransformer> getTransformAlgorithm() const = 0;
+};
+
+/* Complete constructor */
+CV_EXPORTS_W Ptr<ShapeContextDistanceExtractor>
+ createShapeContextDistanceExtractor(int nAngularBins=12, int nRadialBins=4,
+ float innerRadius=0.2, float outerRadius=2, int iterations=3,
+ const Ptr<HistogramCostExtractor> &comparer = createChiHistogramCostExtractor(),
+ const Ptr<ShapeTransformer> &transformer = createThinPlateSplineShapeTransformer());
+
+/***********************************************************************************/
+/***********************************************************************************/
+/***********************************************************************************/
+/*!
+ * Hausdorff distace implementation based on
+ */
+class CV_EXPORTS_W HausdorffDistanceExtractor : public ShapeDistanceExtractor
+{
+public:
+ CV_WRAP virtual void setDistanceFlag(int distanceFlag) = 0;
+ CV_WRAP virtual int getDistanceFlag() const = 0;
+
+ CV_WRAP virtual void setRankProportion(float rankProportion) = 0;
+ CV_WRAP virtual float getRankProportion() const = 0;
+};
+
+/* Constructor */
+CV_EXPORTS_W Ptr<HausdorffDistanceExtractor> createHausdorffDistanceExtractor(int distanceFlag=cv::NORM_L2, float rankProp=0.6);
+
+} // cv
+#endif
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Copyright (C) 2013, OpenCV Foundation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifndef __OPENCV_SHAPE_SHAPE_TRANSFORM_HPP__
+#define __OPENCV_SHAPE_SHAPE_TRANSFORM_HPP__
+#include <vector>
+#include "opencv2/core.hpp"
+#include "opencv2/imgproc.hpp"
+
+namespace cv
+{
+
+/*!
+ * The base class for ShapeTransformer.
+ * This is just to define the common interface for
+ * shape transformation techniques.
+ */
+class CV_EXPORTS_W ShapeTransformer : public Algorithm
+{
+public:
+ /* Estimate, Apply Transformation and return Transforming cost*/
+ CV_WRAP virtual void estimateTransformation(InputArray transformingShape, InputArray targetShape,
+ std::vector<DMatch>& matches) = 0;
+
+ CV_WRAP virtual float applyTransformation(InputArray input, OutputArray output=noArray()) = 0;
+
+ CV_WRAP virtual void warpImage(InputArray transformingImage, OutputArray output,
+ int flags=INTER_LINEAR, int borderMode=BORDER_CONSTANT,
+ const Scalar& borderValue=Scalar()) const = 0;
+};
+
+/***********************************************************************************/
+/***********************************************************************************/
+/*!
+ * Thin Plate Spline Transformation
+ * Implementation of the TPS transformation
+ * according to "Principal Warps: Thin-Plate Splines and the
+ * Decomposition of Deformations" by Juan Manuel Perez for the GSOC 2013
+ */
+
+class CV_EXPORTS_W ThinPlateSplineShapeTransformer : public ShapeTransformer
+{
+public:
+ CV_WRAP virtual void setRegularizationParameter(double beta) = 0;
+ CV_WRAP virtual double getRegularizationParameter() const = 0;
+};
+
+/* Complete constructor */
+CV_EXPORTS_W Ptr<ThinPlateSplineShapeTransformer>
+ createThinPlateSplineShapeTransformer(double regularizationParameter=0);
+
+/***********************************************************************************/
+/***********************************************************************************/
+/*!
+ * Affine Transformation as a derivated from ShapeTransformer
+ */
+
+class CV_EXPORTS_W AffineTransformer : public ShapeTransformer
+{
+public:
+ CV_WRAP virtual void setFullAffine(bool fullAffine) = 0;
+ CV_WRAP virtual bool getFullAffine() const = 0;
+};
+
+/* Complete constructor */
+CV_EXPORTS_W Ptr<AffineTransformer> createAffineTransformer(bool fullAffine);
+
+} // cv
+#endif
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "precomp.hpp"
+
+namespace cv
+{
+
+class AffineTransformerImpl : public AffineTransformer
+{
+public:
+ /* Constructors */
+ AffineTransformerImpl()
+ {
+ fullAffine = true;
+ name_ = "ShapeTransformer.AFF";
+ }
+
+ AffineTransformerImpl(bool _fullAffine)
+ {
+ fullAffine = _fullAffine;
+ name_ = "ShapeTransformer.AFF";
+ }
+
+ /* Destructor */
+ ~AffineTransformerImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual void estimateTransformation(InputArray transformingShape, InputArray targetShape, std::vector<DMatch> &matches);
+ virtual float applyTransformation(InputArray input, OutputArray output=noArray());
+ virtual void warpImage(InputArray transformingImage, OutputArray output,
+ int flags, int borderMode, const Scalar& borderValue) const;
+
+ //! Setters/Getters
+ virtual void setFullAffine(bool _fullAffine) {fullAffine=_fullAffine;}
+ virtual bool getFullAffine() const {return fullAffine;}
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "affine_type" << fullAffine;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ fullAffine = (int)fn["affine_type"];
+ }
+
+private:
+ bool fullAffine;
+ Mat affineMat;
+ float transformCost;
+
+protected:
+ String name_;
+};
+
+void AffineTransformerImpl::warpImage(InputArray transformingImage, OutputArray output,
+ int flags, int borderMode, const Scalar& borderValue) const
+{
+ CV_Assert(!affineMat.empty());
+ warpAffine(transformingImage, output, affineMat, transformingImage.getMat().size(), flags, borderMode, borderValue);
+}
+
+
+static Mat _localAffineEstimate(const std::vector<Point2f>& shape1, const std::vector<Point2f>& shape2,
+ bool fullAfine)
+{
+ Mat out(2,3,CV_32F);
+ int siz=2*shape1.size();
+
+ if (fullAfine)
+ {
+ Mat matM(siz, 6, CV_32F);
+ Mat matP(siz,1,CV_32F);
+ int contPt=0;
+ for (int ii=0; ii<siz; ii++)
+ {
+ Mat therow = Mat::zeros(1,6,CV_32F);
+ if (ii%2==0)
+ {
+ therow.at<float>(0,0)=shape1[contPt].x;
+ therow.at<float>(0,1)=shape1[contPt].y;
+ therow.at<float>(0,2)=1;
+ therow.row(0).copyTo(matM.row(ii));
+ matP.at<float>(ii,0) = shape2[contPt].x;
+ }
+ else
+ {
+ therow.at<float>(0,3)=shape1[contPt].x;
+ therow.at<float>(0,4)=shape1[contPt].y;
+ therow.at<float>(0,5)=1;
+ therow.row(0).copyTo(matM.row(ii));
+ matP.at<float>(ii,0) = shape2[contPt].y;
+ contPt++;
+ }
+ }
+ Mat sol;
+ solve(matM, matP, sol, DECOMP_SVD);
+ out = sol.reshape(0,2);
+ }
+ else
+ {
+ Mat matM(siz, 4, CV_32F);
+ Mat matP(siz,1,CV_32F);
+ int contPt=0;
+ for (int ii=0; ii<siz; ii++)
+ {
+ Mat therow = Mat::zeros(1,4,CV_32F);
+ if (ii%2==0)
+ {
+ therow.at<float>(0,0)=shape1[contPt].x;
+ therow.at<float>(0,1)=shape1[contPt].y;
+ therow.at<float>(0,2)=1;
+ therow.row(0).copyTo(matM.row(ii));
+ matP.at<float>(ii,0) = shape2[contPt].x;
+ }
+ else
+ {
+ therow.at<float>(0,0)=-shape1[contPt].y;
+ therow.at<float>(0,1)=shape1[contPt].x;
+ therow.at<float>(0,3)=1;
+ therow.row(0).copyTo(matM.row(ii));
+ matP.at<float>(ii,0) = shape2[contPt].y;
+ contPt++;
+ }
+ }
+ Mat sol;
+ solve(matM, matP, sol, DECOMP_SVD);
+ out.at<float>(0,0)=sol.at<float>(0,0);
+ out.at<float>(0,1)=sol.at<float>(1,0);
+ out.at<float>(0,2)=sol.at<float>(2,0);
+ out.at<float>(1,0)=-sol.at<float>(1,0);
+ out.at<float>(1,1)=sol.at<float>(0,0);
+ out.at<float>(1,2)=sol.at<float>(3,0);
+ }
+ return out;
+}
+
+void AffineTransformerImpl::estimateTransformation(InputArray _pts1, InputArray _pts2, std::vector<DMatch>& _matches)
+{
+ Mat pts1 = _pts1.getMat();
+ Mat pts2 = _pts2.getMat();
+ CV_Assert((pts1.channels()==2) & (pts1.cols>0) & (pts2.channels()==2) & (pts2.cols>0));
+ CV_Assert(_matches.size()>1);
+
+ if (pts1.type() != CV_32F)
+ pts1.convertTo(pts1, CV_32F);
+ if (pts2.type() != CV_32F)
+ pts2.convertTo(pts2, CV_32F);
+
+ // Use only valid matchings //
+ std::vector<DMatch> matches;
+ for (size_t i=0; i<_matches.size(); i++)
+ {
+ if (_matches[i].queryIdx<pts1.cols &&
+ _matches[i].trainIdx<pts2.cols)
+ {
+ matches.push_back(_matches[i]);
+ }
+ }
+
+ // Organizing the correspondent points in vector style //
+ std::vector<Point2f> shape1; // transforming shape
+ std::vector<Point2f> shape2; // target shape
+ for (size_t i=0; i<matches.size(); i++)
+ {
+ Point2f pt1=pts1.at<Point2f>(0,matches[i].queryIdx);
+ shape1.push_back(pt1);
+
+ Point2f pt2=pts2.at<Point2f>(0,matches[i].trainIdx);
+ shape2.push_back(pt2);
+ }
+
+ // estimateRigidTransform //
+ Mat affine;
+ estimateRigidTransform(shape1, shape2, fullAffine).convertTo(affine, CV_32F);
+
+ if (affine.empty())
+ affine=_localAffineEstimate(shape1, shape2, fullAffine); //In case there is not good solution, just give a LLS based one
+
+ affineMat = affine;
+}
+
+float AffineTransformerImpl::applyTransformation(InputArray inPts, OutputArray outPts)
+{
+ Mat pts1 = inPts.getMat();
+ CV_Assert((pts1.channels()==2) & (pts1.cols>0));
+
+ //Apply transformation in the complete set of points
+ Mat fAffine;
+ transform(pts1, fAffine, affineMat);
+
+ // Ensambling output //
+ if (outPts.needed())
+ {
+ outPts.create(1,fAffine.cols, CV_32FC2);
+ Mat outMat = outPts.getMat();
+ for (int i=0; i<fAffine.cols; i++)
+ outMat.at<Point2f>(0,i)=fAffine.at<Point2f>(0,i);
+ }
+
+ // Updating Transform Cost //
+ Mat Af(2, 2, CV_32F);
+ Af.at<float>(0,0)=affineMat.at<float>(0,0);
+ Af.at<float>(0,1)=affineMat.at<float>(1,0);
+ Af.at<float>(1,0)=affineMat.at<float>(0,1);
+ Af.at<float>(1,1)=affineMat.at<float>(1,1);
+ SVD mysvd(Af, SVD::NO_UV);
+ Mat singVals=mysvd.w;
+ transformCost=std::log((singVals.at<float>(0,0)+FLT_MIN)/(singVals.at<float>(1,0)+FLT_MIN));
+
+ return transformCost;
+}
+
+Ptr <AffineTransformer> createAffineTransformer(bool fullAffine)
+{
+ return Ptr<AffineTransformer>( new AffineTransformerImpl(fullAffine) );
+}
+
+} // cv
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+/*
+ * Implementation of an optimized EMD for histograms based in
+ * the papers "EMD-L1: An efficient and Robust Algorithm
+ * for comparing histogram-based descriptors", by Haibin Ling and
+ * Kazunori Okuda; and "The Earth Mover's Distance is the Mallows
+ * Distance: Some Insights from Statistics", by Elizaveta Levina and
+ * Peter Bickel, based on HAIBIN LING AND KAZUNORI OKADA implementation.
+ */
+
+#include "precomp.hpp"
+#include "emdL1_def.hpp"
+
+
+/****************************************************************************************\
+* EMDL1 Class *
+\****************************************************************************************/
+
+float EmdL1::getEMDL1(cv::Mat &sig1, cv::Mat &sig2)
+{
+ // Initialization
+ CV_Assert((sig1.rows==sig2.rows) & (sig1.cols==sig2.cols) & (!sig1.empty()) & (!sig2.empty()));
+ if(!initBaseTrees(sig1.rows, 1))
+ return -1;
+
+ float *H1=new float[sig1.rows], *H2 = new float[sig2.rows];
+ for (int ii=0; ii<sig1.rows; ii++)
+ {
+ H1[ii]=sig1.at<float>(ii,0);
+ H2[ii]=sig2.at<float>(ii,0);
+ }
+
+ fillBaseTrees(H1,H2); // Initialize histograms
+ greedySolution(); // Construct an initial Basic Feasible solution
+ initBVTree(); // Initialize BVTree
+
+ // Iteration
+ bool bOptimal = false;
+ m_nItr = 0;
+ while(!bOptimal && m_nItr<nMaxIt)
+ {
+ // Derive U=(u_ij) for row i and column j
+ if(m_nItr==0) updateSubtree(m_pRoot);
+ else updateSubtree(m_pEnter->pChild);
+
+ // Optimality test
+ bOptimal = isOptimal();
+
+ // Find new solution
+ if(!bOptimal)
+ findNewSolution();
+ ++m_nItr;
+ }
+ delete [] H1;
+ delete [] H2;
+ // Output the total flow
+ return compuTotalFlow();
+}
+
+void EmdL1::setMaxIteration(int _nMaxIt)
+{
+ nMaxIt=_nMaxIt;
+}
+
+//-- SubFunctions called in the EMD algorithm
+bool EmdL1::initBaseTrees(int n1, int n2, int n3)
+{
+ if(binsDim1==n1 && binsDim2==n2 && binsDim3==n3)
+ return true;
+ binsDim1 = n1;
+ binsDim2 = n2;
+ binsDim3 = n3;
+ if(binsDim1==0 || binsDim2==0) dimension = 0;
+ else dimension = (binsDim3==0)?2:3;
+
+ if(dimension==2)
+ {
+ m_Nodes.resize(binsDim1);
+ m_EdgesUp.resize(binsDim1);
+ m_EdgesRight.resize(binsDim1);
+ for(int i1=0; i1<binsDim1; i1++)
+ {
+ m_Nodes[i1].resize(binsDim2);
+ m_EdgesUp[i1].resize(binsDim2);
+ m_EdgesRight[i1].resize(binsDim2);
+ }
+ m_NBVEdges.resize(binsDim1*binsDim2*4+2);
+ m_auxQueue.resize(binsDim1*binsDim2+2);
+ m_fromLoop.resize(binsDim1*binsDim2+2);
+ m_toLoop.resize(binsDim1*binsDim2+2);
+ }
+ else if(dimension==3)
+ {
+ m_3dNodes.resize(binsDim1);
+ m_3dEdgesUp.resize(binsDim1);
+ m_3dEdgesRight.resize(binsDim1);
+ m_3dEdgesDeep.resize(binsDim1);
+ for(int i1=0; i1<binsDim1; i1++)
+ {
+ m_3dNodes[i1].resize(binsDim2);
+ m_3dEdgesUp[i1].resize(binsDim2);
+ m_3dEdgesRight[i1].resize(binsDim2);
+ m_3dEdgesDeep[i1].resize(binsDim2);
+ for(int i2=0; i2<binsDim2; i2++)
+ {
+ m_3dNodes[i1][i2].resize(binsDim3);
+ m_3dEdgesUp[i1][i2].resize(binsDim3);
+ m_3dEdgesRight[i1][i2].resize(binsDim3);
+ m_3dEdgesDeep[i1][i2].resize(binsDim3);
+ }
+ }
+ m_NBVEdges.resize(binsDim1*binsDim2*binsDim3*6+4);
+ m_auxQueue.resize(binsDim1*binsDim2*binsDim3+4);
+ m_fromLoop.resize(binsDim1*binsDim2*binsDim3+4);
+ m_toLoop.resize(binsDim1*binsDim2*binsDim3+2);
+ }
+ else
+ return false;
+
+ return true;
+}
+
+bool EmdL1::fillBaseTrees(float *H1, float *H2)
+{
+ //- Set global counters
+ m_pRoot = NULL;
+ // Graph initialization
+ float *p1 = H1;
+ float *p2 = H2;
+ if(dimension==2)
+ {
+ for(int c=0; c<binsDim2; c++)
+ {
+ for(int r=0; r<binsDim1; r++)
+ {
+ //- initialize nodes and links
+ m_Nodes[r][c].pos[0] = r;
+ m_Nodes[r][c].pos[1] = c;
+ m_Nodes[r][c].d = *(p1++)-*(p2++);
+ m_Nodes[r][c].pParent = NULL;
+ m_Nodes[r][c].pChild = NULL;
+ m_Nodes[r][c].iLevel = -1;
+
+ //- initialize edges
+ // to the right
+ m_EdgesRight[r][c].pParent = &(m_Nodes[r][c]);
+ m_EdgesRight[r][c].pChild = &(m_Nodes[r][(c+1)%binsDim2]);
+ m_EdgesRight[r][c].flow = 0;
+ m_EdgesRight[r][c].iDir = 1;
+ m_EdgesRight[r][c].pNxt = NULL;
+
+ // to the upward
+ m_EdgesUp[r][c].pParent = &(m_Nodes[r][c]);
+ m_EdgesUp[r][c].pChild = &(m_Nodes[(r+1)%binsDim1][c]);
+ m_EdgesUp[r][c].flow = 0;
+ m_EdgesUp[r][c].iDir = 1;
+ m_EdgesUp[r][c].pNxt = NULL;
+ }
+ }
+ }
+ else if(dimension==3)
+ {
+ for(int z=0; z<binsDim3; z++)
+ {
+ for(int c=0; c<binsDim2; c++)
+ {
+ for(int r=0; r<binsDim1; r++)
+ {
+ //- initialize nodes and edges
+ m_3dNodes[r][c][z].pos[0] = r;
+ m_3dNodes[r][c][z].pos[1] = c;
+ m_3dNodes[r][c][z].pos[2] = z;
+ m_3dNodes[r][c][z].d = *(p1++)-*(p2++);
+ m_3dNodes[r][c][z].pParent = NULL;
+ m_3dNodes[r][c][z].pChild = NULL;
+ m_3dNodes[r][c][z].iLevel = -1;
+
+ //- initialize edges
+ // to the upward
+ m_3dEdgesUp[r][c][z].pParent= &(m_3dNodes[r][c][z]);
+ m_3dEdgesUp[r][c][z].pChild = &(m_3dNodes[(r+1)%binsDim1][c][z]);
+ m_3dEdgesUp[r][c][z].flow = 0;
+ m_3dEdgesUp[r][c][z].iDir = 1;
+ m_3dEdgesUp[r][c][z].pNxt = NULL;
+
+ // to the right
+ m_3dEdgesRight[r][c][z].pParent = &(m_3dNodes[r][c][z]);
+ m_3dEdgesRight[r][c][z].pChild = &(m_3dNodes[r][(c+1)%binsDim2][z]);
+ m_3dEdgesRight[r][c][z].flow = 0;
+ m_3dEdgesRight[r][c][z].iDir = 1;
+ m_3dEdgesRight[r][c][z].pNxt = NULL;
+
+ // to the deep
+ m_3dEdgesDeep[r][c][z].pParent = &(m_3dNodes[r][c][z]);
+ m_3dEdgesDeep[r][c][z].pChild = &(m_3dNodes[r][c])[(z+1)%binsDim3];
+ m_3dEdgesDeep[r][c][z].flow = 0;
+ m_3dEdgesDeep[r][c][z].iDir = 1;
+ m_3dEdgesDeep[r][c][z].pNxt = NULL;
+ }
+ }
+ }
+ }
+ return true;
+}
+
+bool EmdL1::greedySolution()
+{
+ return dimension==2?greedySolution2():greedySolution3();
+}
+
+bool EmdL1::greedySolution2()
+{
+ //- Prepare auxiliary array, D=H1-H2
+ int c,r;
+ floatArray2D D(binsDim1);
+ for(r=0; r<binsDim1; r++)
+ {
+ D[r].resize(binsDim2);
+ for(c=0; c<binsDim2; c++) D[r][c] = m_Nodes[r][c].d;
+ }
+ // compute integrated values along each dimension
+ std::vector<float> d2s(binsDim2);
+ d2s[0] = 0;
+ for(c=0; c<binsDim2-1; c++)
+ {
+ d2s[c+1] = d2s[c];
+ for(r=0; r<binsDim1; r++) d2s[c+1]-= D[r][c];
+ }
+
+ std::vector<float> d1s(binsDim1);
+ d1s[0] = 0;
+ for(r=0; r<binsDim1-1; r++)
+ {
+ d1s[r+1] = d1s[r];
+ for(c=0; c<binsDim2; c++) d1s[r+1]-= D[r][c];
+ }
+
+ //- Greedy algorithm for initial solution
+ cvPEmdEdge pBV;
+ float dFlow;
+ bool bUpward = false;
+ nNBV = 0; // number of NON-BV edges
+
+ for(c=0; c<binsDim2-1; c++)
+ for(r=0; r<binsDim1; r++)
+ {
+ dFlow = D[r][c];
+ bUpward = (r<binsDim1-1) && (fabs(dFlow+d2s[c+1]) > fabs(dFlow+d1s[r+1])); // Move upward or right
+
+ // modify basic variables, record BV and related values
+ if(bUpward)
+ {
+ // move to up
+ pBV = &(m_EdgesUp[r][c]);
+ m_NBVEdges[nNBV++] = &(m_EdgesRight[r][c]);
+ D[r+1][c] += dFlow; // auxilary matrix maintanence
+ d1s[r+1] += dFlow; // auxilary matrix maintanence
+ }
+ else
+ {
+ // move to right, no other choice
+ pBV = &(m_EdgesRight[r][c]);
+ if(r<binsDim1-1)
+ m_NBVEdges[nNBV++] = &(m_EdgesUp[r][c]);
+
+ D[r][c+1] += dFlow; // auxilary matrix maintanence
+ d2s[c+1] += dFlow; // auxilary matrix maintanence
+ }
+ pBV->pParent->pChild = pBV;
+ pBV->flow = fabs(dFlow);
+ pBV->iDir = dFlow>0; // 1:outward, 0:inward
+ }
+
+ //- rightmost column, no choice but move upward
+ c = binsDim2-1;
+ for(r=0; r<binsDim1-1; r++)
+ {
+ dFlow = D[r][c];
+ pBV = &(m_EdgesUp[r][c]);
+ D[r+1][c] += dFlow; // auxilary matrix maintanence
+ pBV->pParent->pChild= pBV;
+ pBV->flow = fabs(dFlow);
+ pBV->iDir = dFlow>0; // 1:outward, 0:inward
+ }
+ return true;
+}
+
+bool EmdL1::greedySolution3()
+{
+ //- Prepare auxiliary array, D=H1-H2
+ int i1,i2,i3;
+ std::vector<floatArray2D> D(binsDim1);
+ for(i1=0; i1<binsDim1; i1++)
+ {
+ D[i1].resize(binsDim2);
+ for(i2=0; i2<binsDim2; i2++)
+ {
+ D[i1][i2].resize(binsDim3);
+ for(i3=0; i3<binsDim3; i3++)
+ D[i1][i2][i3] = m_3dNodes[i1][i2][i3].d;
+ }
+ }
+
+ // compute integrated values along each dimension
+ std::vector<float> d1s(binsDim1);
+ d1s[0] = 0;
+ for(i1=0; i1<binsDim1-1; i1++)
+ {
+ d1s[i1+1] = d1s[i1];
+ for(i2=0; i2<binsDim2; i2++)
+ {
+ for(i3=0; i3<binsDim3; i3++)
+ d1s[i1+1] -= D[i1][i2][i3];
+ }
+ }
+
+ std::vector<float> d2s(binsDim2);
+ d2s[0] = 0;
+ for(i2=0; i2<binsDim2-1; i2++)
+ {
+ d2s[i2+1] = d2s[i2];
+ for(i1=0; i1<binsDim1; i1++)
+ {
+ for(i3=0; i3<binsDim3; i3++)
+ d2s[i2+1] -= D[i1][i2][i3];
+ }
+ }
+
+ std::vector<float> d3s(binsDim3);
+ d3s[0] = 0;
+ for(i3=0; i3<binsDim3-1; i3++)
+ {
+ d3s[i3+1] = d3s[i3];
+ for(i1=0; i1<binsDim1; i1++)
+ {
+ for(i2=0; i2<binsDim2; i2++)
+ d3s[i3+1] -= D[i1][i2][i3];
+ }
+ }
+
+ //- Greedy algorithm for initial solution
+ cvPEmdEdge pBV;
+ float dFlow, f1,f2,f3;
+ nNBV = 0; // number of NON-BV edges
+ for(i3=0; i3<binsDim3; i3++)
+ {
+ for(i2=0; i2<binsDim2; i2++)
+ {
+ for(i1=0; i1<binsDim1; i1++)
+ {
+ if(i3==binsDim3-1 && i2==binsDim2-1 && i1==binsDim1-1) break;
+
+ //- determine which direction to move, either right or upward
+ dFlow = D[i1][i2][i3];
+ f1 = i1<(binsDim1-1)?fabs(dFlow+d1s[i1+1]):VHIGH;
+ f2 = i2<(binsDim2-1)?fabs(dFlow+d2s[i2+1]):VHIGH;
+ f3 = i3<(binsDim3-1)?fabs(dFlow+d3s[i3+1]):VHIGH;
+
+ if(f1<f2 && f1<f3)
+ {
+ pBV = &(m_3dEdgesUp[i1][i2][i3]); // up
+ if(i2<binsDim2-1) m_NBVEdges[nNBV++] = &(m_3dEdgesRight[i1][i2][i3]); // right
+ if(i3<binsDim3-1) m_NBVEdges[nNBV++] = &(m_3dEdgesDeep[i1][i2][i3]); // deep
+ D[i1+1][i2][i3] += dFlow; // maintain auxilary matrix
+ d1s[i1+1] += dFlow;
+ }
+ else if(f2<f3)
+ {
+ pBV = &(m_3dEdgesRight[i1][i2][i3]); // right
+ if(i1<binsDim1-1) m_NBVEdges[nNBV++] = &(m_3dEdgesUp[i1][i2][i3]); // up
+ if(i3<binsDim3-1) m_NBVEdges[nNBV++] = &(m_3dEdgesDeep[i1][i2][i3]); // deep
+ D[i1][i2+1][i3] += dFlow; // maintain auxilary matrix
+ d2s[i2+1] += dFlow;
+ }
+ else
+ {
+ pBV = &(m_3dEdgesDeep[i1][i2][i3]); // deep
+ if(i2<binsDim2-1) m_NBVEdges[nNBV++] = &(m_3dEdgesRight[i1][i2][i3]); // right
+ if(i1<binsDim1-1) m_NBVEdges[nNBV++] = &(m_3dEdgesUp[i1][i2][i3]); // up
+ D[i1][i2][i3+1] += dFlow; // maintain auxilary matrix
+ d3s[i3+1] += dFlow;
+ }
+
+ pBV->flow = fabs(dFlow);
+ pBV->iDir = dFlow>0; // 1:outward, 0:inward
+ pBV->pParent->pChild= pBV;
+ }
+ }
+ }
+ return true;
+}
+
+void EmdL1::initBVTree()
+{
+ // initialize BVTree from the initial BF solution
+ //- Using the center of the graph as the root
+ int r = (int)(0.5*binsDim1-.5);
+ int c = (int)(0.5*binsDim2-.5);
+ int z = (int)(0.5*binsDim3-.5);
+ m_pRoot = dimension==2 ? &(m_Nodes[r][c]) : &(m_3dNodes[r][c][z]);
+ m_pRoot->u = 0;
+ m_pRoot->iLevel = 0;
+ m_pRoot->pParent= NULL;
+ m_pRoot->pPEdge = NULL;
+
+ //- Prepare a queue
+ m_auxQueue[0] = m_pRoot;
+ int nQueue = 1; // length of queue
+ int iQHead = 0; // head of queue
+
+ //- Recursively build subtrees
+ cvPEmdEdge pCurE=NULL, pNxtE=NULL;
+ cvPEmdNode pCurN=NULL, pNxtN=NULL;
+ int nBin = binsDim1*binsDim2*std::max(binsDim3,1);
+ while(iQHead<nQueue && nQueue<nBin)
+ {
+ pCurN = m_auxQueue[iQHead++]; // pop out from queue
+ r = pCurN->pos[0];
+ c = pCurN->pos[1];
+ z = pCurN->pos[2];
+
+ // check connection from itself
+ pCurE = pCurN->pChild; // the initial child from initial solution
+ if(pCurE)
+ {
+ pNxtN = pCurE->pChild;
+ pNxtN->pParent = pCurN;
+ pNxtN->pPEdge = pCurE;
+ m_auxQueue[nQueue++] = pNxtN;
+ }
+
+ // check four neighbor nodes
+ int nNB = dimension==2?4:6;
+ for(int k=0;k<nNB;k++)
+ {
+ if(dimension==2)
+ {
+ if(k==0 && c>0) pNxtN = &(m_Nodes[r][c-1]); // left
+ else if(k==1 && r>0) pNxtN = &(m_Nodes[r-1][c]); // down
+ else if(k==2 && c<binsDim2-1) pNxtN = &(m_Nodes[r][c+1]); // right
+ else if(k==3 && r<binsDim1-1) pNxtN = &(m_Nodes[r+1][c]); // up
+ else continue;
+ }
+ else if(dimension==3)
+ {
+ if(k==0 && c>0) pNxtN = &(m_3dNodes[r][c-1][z]); // left
+ else if(k==1 && c<binsDim2-1) pNxtN = &(m_3dNodes[r][c+1][z]); // right
+ else if(k==2 && r>0) pNxtN = &(m_3dNodes[r-1][c][z]); // down
+ else if(k==3 && r<binsDim1-1) pNxtN = &(m_3dNodes[r+1][c][z]); // up
+ else if(k==4 && z>0) pNxtN = &(m_3dNodes[r][c][z-1]); // shallow
+ else if(k==5 && z<binsDim3-1) pNxtN = &(m_3dNodes[r][c][z+1]); // deep
+ else continue;
+ }
+ if(pNxtN != pCurN->pParent)
+ {
+ pNxtE = pNxtN->pChild;
+ if(pNxtE && pNxtE->pChild==pCurN) // has connection
+ {
+ pNxtN->pParent = pCurN;
+ pNxtN->pPEdge = pNxtE;
+ pNxtN->pChild = NULL;
+ m_auxQueue[nQueue++] = pNxtN;
+
+ pNxtE->pParent = pCurN; // reverse direction
+ pNxtE->pChild = pNxtN;
+ pNxtE->iDir = !pNxtE->iDir;
+
+ if(pCurE) pCurE->pNxt = pNxtE; // add to edge list
+ else pCurN->pChild = pNxtE;
+ pCurE = pNxtE;
+ }
+ }
+ }
+ }
+}
+
+void EmdL1::updateSubtree(cvPEmdNode pRoot)
+{
+ // Initialize auxiliary queue
+ m_auxQueue[0] = pRoot;
+ int nQueue = 1; // queue length
+ int iQHead = 0; // head of queue
+
+ // BFS browing
+ cvPEmdNode pCurN=NULL,pNxtN=NULL;
+ cvPEmdEdge pCurE=NULL;
+ while(iQHead<nQueue)
+ {
+ pCurN = m_auxQueue[iQHead++]; // pop out from queue
+ pCurE = pCurN->pChild;
+
+ // browsing all children
+ while(pCurE)
+ {
+ pNxtN = pCurE->pChild;
+ pNxtN->iLevel = pCurN->iLevel+1;
+ pNxtN->u = pCurE->iDir ? (pCurN->u - 1) : (pCurN->u + 1);
+ pCurE = pCurE->pNxt;
+ m_auxQueue[nQueue++] = pNxtN;
+ }
+ }
+}
+
+bool EmdL1::isOptimal()
+{
+ int iC, iMinC = 0;
+ cvPEmdEdge pE;
+ m_pEnter = NULL;
+ m_iEnter = -1;
+
+ // test each NON-BV edges
+ for(int k=0; k<nNBV; ++k)
+ {
+ pE = m_NBVEdges[k];
+ iC = 1 - pE->pParent->u + pE->pChild->u;
+ if(iC<iMinC)
+ {
+ iMinC = iC;
+ m_iEnter= k;
+ }
+ else
+ {
+ // Try reversing the direction
+ iC = 1 + pE->pParent->u - pE->pChild->u;
+ if(iC<iMinC)
+ {
+ iMinC = iC;
+ m_iEnter= k;
+ }
+ }
+ }
+
+ if(m_iEnter>=0)
+ {
+ m_pEnter = m_NBVEdges[m_iEnter];
+ if(iMinC == (1 - m_pEnter->pChild->u + m_pEnter->pParent->u)) {
+ // reverse direction
+ cvPEmdNode pN = m_pEnter->pParent;
+ m_pEnter->pParent = m_pEnter->pChild;
+ m_pEnter->pChild = pN;
+ }
+
+ m_pEnter->iDir = 1;
+ }
+ return m_iEnter==-1;
+}
+
+void EmdL1::findNewSolution()
+{
+ // Find loop formed by adding the Enter BV edge.
+ findLoopFromEnterBV();
+ // Modify flow values along the loop
+ cvPEmdEdge pE = NULL;
+ float minFlow = m_pLeave->flow;
+ int k;
+ for(k=0; k<m_iFrom; k++)
+ {
+ pE = m_fromLoop[k];
+ if(pE->iDir) pE->flow += minFlow; // outward
+ else pE->flow -= minFlow; // inward
+ }
+ for(k=0; k<m_iTo; k++)
+ {
+ pE = m_toLoop[k];
+ if(pE->iDir) pE->flow -= minFlow; // outward
+ else pE->flow += minFlow; // inward
+ }
+
+ // Update BV Tree, removing the Leaving-BV edge
+ cvPEmdNode pLParentN = m_pLeave->pParent;
+ cvPEmdNode pLChildN = m_pLeave->pChild;
+ cvPEmdEdge pPreE = pLParentN->pChild;
+ if(pPreE==m_pLeave)
+ {
+ pLParentN->pChild = m_pLeave->pNxt; // Leaving-BV is the first child
+ }
+ else
+ {
+ while(pPreE->pNxt != m_pLeave)
+ pPreE = pPreE->pNxt;
+ pPreE->pNxt = m_pLeave->pNxt; // remove Leaving-BV from child list
+ }
+ pLChildN->pParent = NULL;
+ pLChildN->pPEdge = NULL;
+
+ m_NBVEdges[m_iEnter]= m_pLeave; // put the leaving-BV into the NBV array
+
+ // Add the Enter BV edge
+ cvPEmdNode pEParentN = m_pEnter->pParent;
+ cvPEmdNode pEChildN = m_pEnter->pChild;
+ m_pEnter->flow = minFlow;
+ m_pEnter->pNxt = pEParentN->pChild; // insert the Enter BV as the first child
+ pEParentN->pChild = m_pEnter; // of its parent
+
+ // Recursively update the tree start from pEChildN
+ cvPEmdNode pPreN = pEParentN;
+ cvPEmdNode pCurN = pEChildN;
+ cvPEmdNode pNxtN;
+ cvPEmdEdge pNxtE, pPreE0;
+ pPreE = m_pEnter;
+ while(pCurN)
+ {
+ pNxtN = pCurN->pParent;
+ pNxtE = pCurN->pPEdge;
+ pCurN->pParent = pPreN;
+ pCurN->pPEdge = pPreE;
+ if(pNxtN)
+ {
+ // remove the edge from pNxtN's child list
+ if(pNxtN->pChild==pNxtE)
+ {
+ pNxtN->pChild = pNxtE->pNxt; // first child
+ }
+ else
+ {
+ pPreE0 = pNxtN->pChild;
+ while(pPreE0->pNxt != pNxtE)
+ pPreE0 = pPreE0->pNxt;
+ pPreE0->pNxt = pNxtE->pNxt; // remove Leaving-BV from child list
+ }
+ // reverse the parent-child direction
+ pNxtE->pParent = pCurN;
+ pNxtE->pChild = pNxtN;
+ pNxtE->iDir = !pNxtE->iDir;
+ pNxtE->pNxt = pCurN->pChild;
+ pCurN->pChild = pNxtE;
+ pPreE = pNxtE;
+ pPreN = pCurN;
+ }
+ pCurN = pNxtN;
+ }
+
+ // Update U at the child of the Enter BV
+ pEChildN->u = m_pEnter->iDir?(pEParentN->u-1):(pEParentN->u + 1);
+ pEChildN->iLevel = pEParentN->iLevel+1;
+}
+
+void EmdL1::findLoopFromEnterBV()
+{
+ // Initialize Leaving-BV edge
+ float minFlow = VHIGH;
+ cvPEmdEdge pE = NULL;
+ int iLFlag = 0; // 0: in the FROM list, 1: in the TO list
+
+ // Using two loop list to store the loop nodes
+ cvPEmdNode pFrom = m_pEnter->pParent;
+ cvPEmdNode pTo = m_pEnter->pChild;
+ m_iFrom = 0;
+ m_iTo = 0;
+ m_pLeave = NULL;
+
+ // Trace back to make pFrom and pTo at the same level
+ while(pFrom->iLevel > pTo->iLevel)
+ {
+ pE = pFrom->pPEdge;
+ m_fromLoop[m_iFrom++] = pE;
+ if(!pE->iDir && pE->flow<minFlow)
+ {
+ minFlow = pE->flow;
+ m_pLeave = pE;
+ iLFlag = 0; // 0: in the FROM list
+ }
+ pFrom = pFrom->pParent;
+ }
+
+ while(pTo->iLevel > pFrom->iLevel)
+ {
+ pE = pTo->pPEdge;
+ m_toLoop[m_iTo++] = pE;
+ if(pE->iDir && pE->flow<minFlow)
+ {
+ minFlow = pE->flow;
+ m_pLeave = pE;
+ iLFlag = 1; // 1: in the TO list
+ }
+ pTo = pTo->pParent;
+ }
+
+ // Trace pTo and pFrom simultaneously till find their common ancester
+ while(pTo!=pFrom)
+ {
+ pE = pFrom->pPEdge;
+ m_fromLoop[m_iFrom++] = pE;
+ if(!pE->iDir && pE->flow<minFlow)
+ {
+ minFlow = pE->flow;
+ m_pLeave = pE;
+ iLFlag = 0; // 0: in the FROM list, 1: in the TO list
+ }
+ pFrom = pFrom->pParent;
+
+ pE = pTo->pPEdge;
+ m_toLoop[m_iTo++] = pE;
+ if(pE->iDir && pE->flow<minFlow)
+ {
+ minFlow = pE->flow;
+ m_pLeave = pE;
+ iLFlag = 1; // 0: in the FROM list, 1: in the TO list
+ }
+ pTo = pTo->pParent;
+ }
+
+ // Reverse the direction of the Enter BV edge if necessary
+ if(iLFlag==0)
+ {
+ cvPEmdNode pN = m_pEnter->pParent;
+ m_pEnter->pParent = m_pEnter->pChild;
+ m_pEnter->pChild = pN;
+ m_pEnter->iDir = !m_pEnter->iDir;
+ }
+}
+
+float EmdL1::compuTotalFlow()
+{
+ // Computing the total flow as the final distance
+ float f = 0;
+
+ // Initialize auxiliary queue
+ m_auxQueue[0] = m_pRoot;
+ int nQueue = 1; // length of queue
+ int iQHead = 0; // head of queue
+
+ // BFS browing the tree
+ cvPEmdNode pCurN=NULL,pNxtN=NULL;
+ cvPEmdEdge pCurE=NULL;
+ while(iQHead<nQueue)
+ {
+ pCurN = m_auxQueue[iQHead++]; // pop out from queue
+ pCurE = pCurN->pChild;
+
+ // browsing all children
+ while(pCurE)
+ {
+ f += pCurE->flow;
+ pNxtN = pCurE->pChild;
+ pCurE = pCurE->pNxt;
+ m_auxQueue[nQueue++] = pNxtN;
+ }
+ }
+ return f;
+}
+
+/****************************************************************************************\
+* EMDL1 Function *
+\****************************************************************************************/
+
+float cv::EMDL1(InputArray _signature1, InputArray _signature2)
+{
+ Mat signature1 = _signature1.getMat(), signature2 = _signature2.getMat();
+ EmdL1 emdl1;
+ return emdl1.getEMDL1(signature1, signature2);
+}
+
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include <stdlib.h>
+#include <math.h>
+#include <vector>
+
+#define VHIGH 1e10;
+/****************************************************************************************\
+* For EMDL1 Framework *
+\****************************************************************************************/
+typedef struct cvEMDEdge* cvPEmdEdge;
+typedef struct cvEMDNode* cvPEmdNode;
+struct cvEMDNode
+{
+ int pos[3]; // grid position
+ float d; // initial value
+ int u;
+ // tree maintainance
+ int iLevel; // level in the tree, 0 means root
+ cvPEmdNode pParent; // pointer to its parent
+ cvPEmdEdge pChild;
+ cvPEmdEdge pPEdge; // point to the edge coming out from its parent
+};
+struct cvEMDEdge
+{
+ float flow; // initial value
+ int iDir; // 1:outward, 0:inward
+ // tree maintainance
+ cvPEmdNode pParent; // point to its parent
+ cvPEmdNode pChild; // the child node
+ cvPEmdEdge pNxt; // next child/edge
+};
+typedef std::vector<cvEMDNode> cvEMDNodeArray;
+typedef std::vector<cvEMDEdge> cvEMDEdgeArray;
+typedef std::vector<cvEMDNodeArray> cvEMDNodeArray2D;
+typedef std::vector<cvEMDEdgeArray> cvEMDEdgeArray2D;
+typedef std::vector<float> floatArray;
+typedef std::vector<floatArray> floatArray2D;
+
+/****************************************************************************************\
+* EMDL1 Class *
+\****************************************************************************************/
+class EmdL1
+{
+public:
+ EmdL1()
+ {
+ m_pRoot = NULL;
+ binsDim1 = 0;
+ binsDim2 = 0;
+ binsDim3 = 0;
+ dimension = 0;
+ nMaxIt = 500;
+ }
+
+ ~EmdL1()
+ {
+ }
+
+ float getEMDL1(cv::Mat &sig1, cv::Mat &sig2);
+ void setMaxIteration(int _nMaxIt);
+
+private:
+ //-- SubFunctions called in the EMD algorithm
+ bool initBaseTrees(int n1=0, int n2=0, int n3=0);
+ bool fillBaseTrees(float *H1, float *H2);
+ bool greedySolution();
+ bool greedySolution2();
+ bool greedySolution3();
+ void initBVTree();
+ void updateSubtree(cvPEmdNode pRoot);
+ bool isOptimal();
+ void findNewSolution();
+ void findLoopFromEnterBV();
+ float compuTotalFlow();
+
+private:
+ int dimension;
+ int binsDim1, binsDim2, binsDim3; // the hitogram contains m_n1 rows and m_n2 columns
+ int nNBV; // number of Non-Basic Variables (NBV)
+ int nMaxIt;
+ cvEMDNodeArray2D m_Nodes; // all nodes
+ cvEMDEdgeArray2D m_EdgesRight; // all edges to right
+ cvEMDEdgeArray2D m_EdgesUp; // all edges to upward
+ std::vector<cvEMDNodeArray2D> m_3dNodes; // all nodes for 3D
+ std::vector<cvEMDEdgeArray2D> m_3dEdgesRight; // all edges to right, 3D
+ std::vector<cvEMDEdgeArray2D> m_3dEdgesUp; // all edges to upward, 3D
+ std::vector<cvEMDEdgeArray2D> m_3dEdgesDeep; // all edges to deep, 3D
+ std::vector<cvPEmdEdge> m_NBVEdges; // pointers to all NON-BV edges
+ std::vector<cvPEmdNode> m_auxQueue; // auxiliary node queue
+ cvPEmdNode m_pRoot; // root of the BV Tree
+ cvPEmdEdge m_pEnter; // Enter BV edge
+ int m_iEnter; // Enter BV edge, index in m_NBVEdges
+ cvPEmdEdge m_pLeave; // Leave BV edge
+ int m_nItr; // number of iteration
+ // auxiliary variables for searching a new loop
+ std::vector<cvPEmdEdge> m_fromLoop;
+ std::vector<cvPEmdEdge> m_toLoop;
+ int m_iFrom;
+ int m_iTo;
+};
+
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// Intel License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000, Intel Corporation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of Intel Corporation may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "precomp.hpp"
+
+namespace cv
+{
+
+class HausdorffDistanceExtractorImpl : public HausdorffDistanceExtractor
+{
+public:
+ /* Constructor */
+ HausdorffDistanceExtractorImpl(int _distanceFlag = NORM_L1, float _rankProportion=0.6)
+ {
+ distanceFlag = _distanceFlag;
+ rankProportion = _rankProportion;
+ name_ = "ShapeDistanceExtractor.HAU";
+ }
+
+ /* Destructor */
+ ~HausdorffDistanceExtractorImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual float computeDistance(InputArray contour1, InputArray contour2);
+
+ //! Setters/Getters
+ virtual void setDistanceFlag(int _distanceFlag) {distanceFlag=_distanceFlag;}
+ virtual int getDistanceFlag() const {return distanceFlag;}
+
+ virtual void setRankProportion(float _rankProportion)
+ {
+ CV_Assert((_rankProportion>0) & (_rankProportion<=1));
+ rankProportion=_rankProportion;
+ }
+ virtual float getRankProportion() const {return rankProportion;}
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "distance" << distanceFlag
+ << "rank" << rankProportion;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ distanceFlag = (int)fn["distance"];
+ rankProportion = (int)fn["rank"];
+ }
+
+private:
+ int distanceFlag;
+ float rankProportion;
+
+protected:
+ String name_;
+};
+
+//! Hausdorff distance for a pair of set of points
+static float _apply(const Mat &set1, const Mat &set2, int distType, double propRank)
+{
+ // Building distance matrix //
+ Mat disMat(set1.cols, set2.cols, CV_32F);
+ int K = int(propRank*(disMat.rows-1));
+
+ for (int r=0; r<disMat.rows; r++)
+ {
+ for (int c=0; c<disMat.cols; c++)
+ {
+ Point2f diff = set1.at<Point2f>(0,r)-set2.at<Point2f>(0,c);
+ disMat.at<float>(r,c) = norm(Mat(diff), distType);
+ }
+ }
+
+ Mat shortest(disMat.rows,1,CV_32F);
+ for (int ii=0; ii<disMat.rows; ii++)
+ {
+ Mat therow = disMat.row(ii);
+ double mindis;
+ minMaxIdx(therow, &mindis);
+ shortest.at<float>(ii,0) = float(mindis);
+ }
+ Mat sorted;
+ cv::sort(shortest, sorted, SORT_EVERY_ROW | SORT_DESCENDING);
+ return sorted.at<float>(K,0);
+}
+
+float HausdorffDistanceExtractorImpl::computeDistance(InputArray contour1, InputArray contour2)
+{
+ Mat set1=contour1.getMat(), set2=contour2.getMat();
+ if (set1.type() != CV_32F)
+ set1.convertTo(set1, CV_32F);
+ if (set2.type() != CV_32F)
+ set2.convertTo(set2, CV_32F);
+ CV_Assert((set1.channels()==2) & (set1.cols>0));
+ CV_Assert((set2.channels()==2) & (set2.cols>0));
+ return std::max( _apply(set1, set2, distanceFlag, rankProportion),
+ _apply(set2, set1, distanceFlag, rankProportion) );
+}
+
+Ptr <HausdorffDistanceExtractor> createHausdorffDistanceExtractor(int distanceFlag, float rankProp)
+{
+ return Ptr<HausdorffDistanceExtractor>(new HausdorffDistanceExtractorImpl(distanceFlag, rankProp));
+}
+
+} // cv
+
+
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// Intel License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000, Intel Corporation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of Intel Corporation may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "precomp.hpp"
+
+namespace cv
+{
+
+/*! */
+class NormHistogramCostExtractorImpl : public NormHistogramCostExtractor
+{
+public:
+ /* Constructors */
+ NormHistogramCostExtractorImpl(int _flag, int _nDummies, float _defaultCost)
+ {
+ flag=_flag;
+ nDummies=_nDummies;
+ defaultCost=_defaultCost;
+ name_ = "HistogramCostExtractor.NOR";
+ }
+
+ /* Destructor */
+ ~NormHistogramCostExtractorImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual void buildCostMatrix(InputArray descriptors1, InputArray descriptors2, OutputArray costMatrix);
+
+ //! Setters/Getters
+ void setNDummies(int _nDummies)
+ {
+ nDummies=_nDummies;
+ }
+
+ int getNDummies() const
+ {
+ return nDummies;
+ }
+
+ void setDefaultCost(float _defaultCost)
+ {
+ defaultCost=_defaultCost;
+ }
+
+ float getDefaultCost() const
+ {
+ return defaultCost;
+ }
+
+ virtual void setNormFlag(int _flag)
+ {
+ flag=_flag;
+ }
+
+ virtual int getNormFlag() const
+ {
+ return flag;
+ }
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "flag" << flag
+ << "dummies" << nDummies
+ << "default" << defaultCost;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ flag = (int)fn["flag"];
+ nDummies = (int)fn["dummies"];
+ defaultCost = (float)fn["default"];
+ }
+
+private:
+ int flag;
+ int nDummies;
+ float defaultCost;
+
+protected:
+ String name_;
+};
+
+void NormHistogramCostExtractorImpl::buildCostMatrix(InputArray _descriptors1, InputArray _descriptors2, OutputArray _costMatrix)
+{
+ // size of the costMatrix with dummies //
+ Mat descriptors1=_descriptors1.getMat();
+ Mat descriptors2=_descriptors2.getMat();
+ int costrows = std::max(descriptors1.rows, descriptors2.rows)+nDummies;
+ _costMatrix.create(costrows, costrows, CV_32F);
+ Mat costMatrix=_costMatrix.getMat();
+
+
+ // Obtain copies of the descriptors //
+ cv::Mat scd1 = descriptors1.clone();
+ cv::Mat scd2 = descriptors2.clone();
+
+ // row normalization //
+ for(int i=0; i<scd1.rows; i++)
+ {
+ scd1.row(i)/=(sum(scd1.row(i))[0]+FLT_EPSILON);
+ }
+ for(int i=0; i<scd2.rows; i++)
+ {
+ scd2.row(i)/=(sum(scd2.row(i))[0]+FLT_EPSILON);
+ }
+
+ // Compute the Cost Matrix //
+ for(int i=0; i<costrows; i++)
+ {
+ for(int j=0; j<costrows; j++)
+ {
+ if (i<scd1.rows && j<scd2.rows)
+ {
+ Mat columnDiff = scd1.row(i)-scd2.row(j);
+ costMatrix.at<float>(i,j)=norm(columnDiff, flag);
+ }
+ else
+ {
+ costMatrix.at<float>(i,j)=defaultCost;
+ }
+ }
+ }
+}
+
+Ptr <HistogramCostExtractor> createNormHistogramCostExtractor(int flag, int nDummies, float defaultCost)
+{
+ return Ptr <HistogramCostExtractor>( new NormHistogramCostExtractorImpl(flag, nDummies, defaultCost) );
+}
+
+/*! */
+class EMDHistogramCostExtractorImpl : public EMDHistogramCostExtractor
+{
+public:
+ /* Constructors */
+ EMDHistogramCostExtractorImpl(int _flag, int _nDummies, float _defaultCost)
+ {
+ flag=_flag;
+ nDummies=_nDummies;
+ defaultCost=_defaultCost;
+ name_ = "HistogramCostExtractor.EMD";
+ }
+
+ /* Destructor */
+ ~EMDHistogramCostExtractorImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual void buildCostMatrix(InputArray descriptors1, InputArray descriptors2, OutputArray costMatrix);
+
+ //! Setters/Getters
+ void setNDummies(int _nDummies)
+ {
+ nDummies=_nDummies;
+ }
+
+ int getNDummies() const
+ {
+ return nDummies;
+ }
+
+ void setDefaultCost(float _defaultCost)
+ {
+ defaultCost=_defaultCost;
+ }
+
+ float getDefaultCost() const
+ {
+ return defaultCost;
+ }
+
+ virtual void setNormFlag(int _flag)
+ {
+ flag=_flag;
+ }
+
+ virtual int getNormFlag() const
+ {
+ return flag;
+ }
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "flag" << flag
+ << "dummies" << nDummies
+ << "default" << defaultCost;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ flag = (int)fn["flag"];
+ nDummies = (int)fn["dummies"];
+ defaultCost = (float)fn["default"];
+ }
+
+private:
+ int flag;
+ int nDummies;
+ float defaultCost;
+
+protected:
+ String name_;
+};
+
+void EMDHistogramCostExtractorImpl::buildCostMatrix(InputArray _descriptors1, InputArray _descriptors2, OutputArray _costMatrix)
+{
+ // size of the costMatrix with dummies //
+ Mat descriptors1=_descriptors1.getMat();
+ Mat descriptors2=_descriptors2.getMat();
+ int costrows = std::max(descriptors1.rows, descriptors2.rows)+nDummies;
+ _costMatrix.create(costrows, costrows, CV_32F);
+ Mat costMatrix=_costMatrix.getMat();
+
+ // Obtain copies of the descriptors //
+ cv::Mat scd1=descriptors1.clone();
+ cv::Mat scd2=descriptors2.clone();
+
+ // row normalization //
+ for(int i=0; i<scd1.rows; i++)
+ {
+ cv::Mat row = scd1.row(i);
+ scd1.row(i)/=(sum(row)[0]+FLT_EPSILON);
+ }
+ for(int i=0; i<scd2.rows; i++)
+ {
+ cv::Mat row = scd2.row(i);
+ scd2.row(i)/=(sum(row)[0]+FLT_EPSILON);
+ }
+
+ // Compute the Cost Matrix //
+ for(int i=0; i<costrows; i++)
+ {
+ for(int j=0; j<costrows; j++)
+ {
+ if (i<scd1.rows && j<scd2.rows)
+ {
+ cv::Mat sig1(scd1.cols,2,CV_32F), sig2(scd2.cols,2,CV_32F);
+ sig1.col(0)=scd1.row(i).t();
+ sig2.col(0)=scd2.row(j).t();
+ for (int k=0; k<sig1.rows; k++)
+ {
+ sig1.at<float>(k,1)=k;
+ }
+ for (int k=0; k<sig2.rows; k++)
+ {
+ sig2.at<float>(k,1)=k;
+ }
+
+ costMatrix.at<float>(i,j) = cv::EMD(sig1, sig2, flag);
+ }
+ else
+ {
+ costMatrix.at<float>(i,j) = defaultCost;
+ }
+ }
+ }
+}
+
+Ptr <HistogramCostExtractor> createEMDHistogramCostExtractor(int flag, int nDummies, float defaultCost)
+{
+ return Ptr <HistogramCostExtractor>( new EMDHistogramCostExtractorImpl(flag, nDummies, defaultCost) );
+}
+
+/*! */
+class ChiHistogramCostExtractorImpl : public ChiHistogramCostExtractor
+{
+public:
+ /* Constructors */
+ ChiHistogramCostExtractorImpl(int _nDummies, float _defaultCost)
+ {
+ name_ = "HistogramCostExtractor.CHI";
+ nDummies=_nDummies;
+ defaultCost=_defaultCost;
+ }
+
+ /* Destructor */
+ ~ChiHistogramCostExtractorImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual void buildCostMatrix(InputArray descriptors1, InputArray descriptors2, OutputArray costMatrix);
+
+ //! setters / getters
+ void setNDummies(int _nDummies)
+ {
+ nDummies=_nDummies;
+ }
+
+ int getNDummies() const
+ {
+ return nDummies;
+ }
+
+ void setDefaultCost(float _defaultCost)
+ {
+ defaultCost=_defaultCost;
+ }
+
+ float getDefaultCost() const
+ {
+ return defaultCost;
+ }
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "dummies" << nDummies
+ << "default" << defaultCost;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ nDummies = (int)fn["dummies"];
+ defaultCost = (float)fn["default"];
+ }
+
+protected:
+ String name_;
+ int nDummies;
+ float defaultCost;
+};
+
+void ChiHistogramCostExtractorImpl::buildCostMatrix(InputArray _descriptors1, InputArray _descriptors2, OutputArray _costMatrix)
+{
+ // size of the costMatrix with dummies //
+ Mat descriptors1=_descriptors1.getMat();
+ Mat descriptors2=_descriptors2.getMat();
+ int costrows = std::max(descriptors1.rows, descriptors2.rows)+nDummies;
+ _costMatrix.create(costrows, costrows, CV_32FC1);
+ Mat costMatrix=_costMatrix.getMat();
+
+ // Obtain copies of the descriptors //
+ cv::Mat scd1=descriptors1.clone();
+ cv::Mat scd2=descriptors2.clone();
+
+ // row normalization //
+ for(int i=0; i<scd1.rows; i++)
+ {
+ cv::Mat row = scd1.row(i);
+ scd1.row(i)/=(sum(row)[0]+FLT_EPSILON);
+ }
+ for(int i=0; i<scd2.rows; i++)
+ {
+ cv::Mat row = scd2.row(i);
+ scd2.row(i)/=(sum(row)[0]+FLT_EPSILON);
+ }
+
+ // Compute the Cost Matrix //
+ for(int i=0; i<costrows; i++)
+ {
+ for(int j=0; j<costrows; j++)
+ {
+ if (i<scd1.rows && j<scd2.rows)
+ {
+ float csum = 0;
+ for(int k=0; k<scd2.cols; k++)
+ {
+ float resta=scd1.at<float>(i,k)-scd2.at<float>(j,k);
+ float suma=scd1.at<float>(i,k)+scd2.at<float>(j,k);
+ csum += resta*resta/(FLT_EPSILON+suma);
+ }
+ costMatrix.at<float>(i,j)=csum/2;
+ }
+ else
+ {
+ costMatrix.at<float>(i,j)=defaultCost;
+ }
+ }
+ }
+}
+
+Ptr <HistogramCostExtractor> createChiHistogramCostExtractor(int nDummies, float defaultCost)
+{
+ return Ptr <HistogramCostExtractor>( new ChiHistogramCostExtractorImpl(nDummies, defaultCost) );
+}
+
+/*! */
+class EMDL1HistogramCostExtractorImpl : public EMDL1HistogramCostExtractor
+{
+public:
+ /* Constructors */
+ EMDL1HistogramCostExtractorImpl(int _nDummies, float _defaultCost)
+ {
+ name_ = "HistogramCostExtractor.CHI";
+ nDummies=_nDummies;
+ defaultCost=_defaultCost;
+ }
+
+ /* Destructor */
+ ~EMDL1HistogramCostExtractorImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual void buildCostMatrix(InputArray descriptors1, InputArray descriptors2, OutputArray costMatrix);
+
+ //! setters / getters
+ void setNDummies(int _nDummies)
+ {
+ nDummies=_nDummies;
+ }
+
+ int getNDummies() const
+ {
+ return nDummies;
+ }
+
+ void setDefaultCost(float _defaultCost)
+ {
+ defaultCost=_defaultCost;
+ }
+
+ float getDefaultCost() const
+ {
+ return defaultCost;
+ }
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "dummies" << nDummies
+ << "default" << defaultCost;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ nDummies = (int)fn["dummies"];
+ defaultCost = (float)fn["default"];
+ }
+
+protected:
+ String name_;
+ int nDummies;
+ float defaultCost;
+};
+
+void EMDL1HistogramCostExtractorImpl::buildCostMatrix(InputArray _descriptors1, InputArray _descriptors2, OutputArray _costMatrix)
+{
+ // size of the costMatrix with dummies //
+ Mat descriptors1=_descriptors1.getMat();
+ Mat descriptors2=_descriptors2.getMat();
+ int costrows = std::max(descriptors1.rows, descriptors2.rows)+nDummies;
+ _costMatrix.create(costrows, costrows, CV_32F);
+ Mat costMatrix=_costMatrix.getMat();
+
+ // Obtain copies of the descriptors //
+ cv::Mat scd1=descriptors1.clone();
+ cv::Mat scd2=descriptors2.clone();
+
+ // row normalization //
+ for(int i=0; i<scd1.rows; i++)
+ {
+ cv::Mat row = scd1.row(i);
+ scd1.row(i)/=(sum(row)[0]+FLT_EPSILON);
+ }
+ for(int i=0; i<scd2.rows; i++)
+ {
+ cv::Mat row = scd2.row(i);
+ scd2.row(i)/=(sum(row)[0]+FLT_EPSILON);
+ }
+
+ // Compute the Cost Matrix //
+ for(int i=0; i<costrows; i++)
+ {
+ for(int j=0; j<costrows; j++)
+ {
+ if (i<scd1.rows && j<scd2.rows)
+ {
+ cv::Mat sig1(scd1.cols,1,CV_32F), sig2(scd2.cols,1,CV_32F);
+ sig1.col(0)=scd1.row(i).t();
+ sig2.col(0)=scd2.row(j).t();
+ costMatrix.at<float>(i,j) = cv::EMDL1(sig1, sig2);
+ }
+ else
+ {
+ costMatrix.at<float>(i,j) = defaultCost;
+ }
+ }
+ }
+}
+
+Ptr <HistogramCostExtractor> createEMDL1HistogramCostExtractor(int nDummies, float defaultCost)
+{
+ return Ptr <HistogramCostExtractor>( new EMDL1HistogramCostExtractorImpl(nDummies, defaultCost) );
+}
+
+} // cv
+
+
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "precomp.hpp"
+
+/* End of file. */
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#ifndef __OPENCV_PRECOMP_H__
+#define __OPENCV_PRECOMP_H__
+
+#include <vector>
+#include <cmath>
+#include <iostream>
+
+#include "opencv2/video/tracking.hpp"
+#include "opencv2/imgproc.hpp"
+#include "opencv2/shape.hpp"
+
+#include "opencv2/core/utility.hpp"
+#include "opencv2/core/private.hpp"
+
+#include "opencv2/opencv_modules.hpp"
+
+#endif
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+/*
+ * Implementation of the paper Shape Matching and Object Recognition Using Shape Contexts
+ * Belongie et al., 2002 by Juan Manuel Perez for GSoC 2013.
+ */
+#include "precomp.hpp"
+//#include "opencv2/highgui.hpp"
+/*
+ * ShapeContextDescriptor class
+ */
+class SCD
+{
+public:
+ //! the full constructor taking all the necessary parameters
+ explicit SCD(int _nAngularBins=12, int _nRadialBins=5,
+ double _innerRadius=0.1, double _outerRadius=1, bool _rotationInvariant=false)
+ {
+ setAngularBins(_nAngularBins);
+ setRadialBins(_nRadialBins);
+ setInnerRadius(_innerRadius);
+ setOuterRadius(_outerRadius);
+ setRotationInvariant(_rotationInvariant);
+ }
+
+ void extractSCD(cv::Mat& contour, cv::Mat& descriptors,
+ const std::vector<int>& queryInliers=std::vector<int>(),
+ const float _meanDistance=-1)
+ {
+ cv::Mat contourMat = contour;
+ cv::Mat disMatrix = cv::Mat::zeros(contourMat.cols, contourMat.cols, CV_32F);
+ cv::Mat angleMatrix = cv::Mat::zeros(contourMat.cols, contourMat.cols, CV_32F);
+
+ std::vector<double> logspaces, angspaces;
+ logarithmicSpaces(logspaces);
+ angularSpaces(angspaces);
+ buildNormalizedDistanceMatrix(contourMat, disMatrix, queryInliers, _meanDistance);
+ buildAngleMatrix(contourMat, angleMatrix);
+
+ // Now, build the descriptor matrix (each row is a point) //
+ descriptors = cv::Mat::zeros(contourMat.cols, descriptorSize(), CV_32F);
+
+ for (int ptidx=0; ptidx<contourMat.cols; ptidx++)
+ {
+ for (int cmp=0; cmp<contourMat.cols; cmp++)
+ {
+ if (ptidx==cmp) continue;
+ if ((int)queryInliers.size()>0)
+ {
+ if (queryInliers[ptidx]==0 || queryInliers[cmp]==0) continue; //avoid outliers
+ }
+
+ int angidx=-1, radidx=-1;
+ for (int i=0; i<nRadialBins; i++)
+ {
+ if (disMatrix.at<float>(ptidx, cmp)<logspaces[i])
+ {
+ radidx=i;
+ break;
+ }
+ }
+ for (int i=0; i<nAngularBins; i++)
+ {
+ if (angleMatrix.at<float>(ptidx, cmp)<angspaces[i])
+ {
+ angidx=i;
+ break;
+ }
+ }
+ if (angidx!=-1 && radidx!=-1)
+ {
+ int idx = angidx+radidx*nAngularBins;
+ descriptors.at<float>(ptidx, idx)++;
+ }
+ }
+ }
+ }
+
+ int descriptorSize() {return nAngularBins*nRadialBins;}
+ void setAngularBins(int angularBins) { nAngularBins=angularBins; }
+ void setRadialBins(int radialBins) { nRadialBins=radialBins; }
+ void setInnerRadius(double _innerRadius) { innerRadius=_innerRadius; }
+ void setOuterRadius(double _outerRadius) { outerRadius=_outerRadius; }
+ void setRotationInvariant(bool _rotationInvariant) { rotationInvariant=_rotationInvariant; }
+ int getAngularBins() const { return nAngularBins; }
+ int getRadialBins() const { return nRadialBins; }
+ double getInnerRadius() const { return innerRadius; }
+ double getOuterRadius() const { return outerRadius; }
+ bool getRotationInvariant() const { return rotationInvariant; }
+ float getMeanDistance() const { return meanDistance; }
+
+private:
+ int nAngularBins;
+ int nRadialBins;
+ double innerRadius;
+ double outerRadius;
+ bool rotationInvariant;
+ float meanDistance;
+
+protected:
+ void logarithmicSpaces(std::vector<double>& vecSpaces) const
+ {
+ double logmin=log10(innerRadius);
+ double logmax=log10(outerRadius);
+ double delta=(logmax-logmin)/(nRadialBins-1);
+ double accdelta=0;
+
+ for (int i=0; i<nRadialBins; i++)
+ {
+ double val = std::pow(10,logmin+accdelta);
+ vecSpaces.push_back(val);
+ accdelta += delta;
+ }
+ }
+
+ void angularSpaces(std::vector<double>& vecSpaces) const
+ {
+ double delta=2*CV_PI/nAngularBins;
+ double val=0;
+
+ for (int i=0; i<nAngularBins; i++)
+ {
+ val += delta;
+ vecSpaces.push_back(val);
+ }
+ }
+
+ void buildNormalizedDistanceMatrix(cv::Mat& contour,
+ cv::Mat& disMatrix, const std::vector<int> &queryInliers,
+ const float _meanDistance=-1)
+ {
+ cv::Mat contourMat = contour;
+ cv::Mat mask(disMatrix.rows, disMatrix.cols, CV_8U);
+
+ for (int i=0; i<contourMat.cols; i++)
+ {
+ for (int j=0; j<contourMat.cols; j++)
+ {
+ disMatrix.at<float>(i,j) = norm( cv::Mat(contourMat.at<cv::Point2f>(0,i)-contourMat.at<cv::Point2f>(0,j)), cv::NORM_L2 );
+ if (_meanDistance<0)
+ {
+ if (queryInliers.size()>0)
+ {
+ mask.at<char>(i,j)=char(queryInliers[j] & queryInliers[i]);
+ }
+ else
+ {
+ mask.at<char>(i,j)=1;
+ }
+ }
+ }
+ }
+
+ if (_meanDistance<0)
+ {
+ meanDistance=mean(disMatrix, mask)[0];
+ }
+ else
+ {
+ meanDistance=_meanDistance;
+ }
+ disMatrix/=meanDistance+FLT_EPSILON;
+ }
+
+ void buildAngleMatrix(cv::Mat& contour,
+ cv::Mat& angleMatrix) const
+ {
+ cv::Mat contourMat = contour;
+
+ // if descriptor is rotationInvariant compute massCenter //
+ cv::Point2f massCenter(0,0);
+ if (rotationInvariant)
+ {
+ for (int i=0; i<contourMat.cols; i++)
+ {
+ massCenter+=contourMat.at<cv::Point2f>(0,i);
+ }
+ massCenter.x=massCenter.x/(float)contourMat.cols;
+ massCenter.y=massCenter.y/(float)contourMat.cols;
+ }
+
+
+ for (int i=0; i<contourMat.cols; i++)
+ {
+ for (int j=0; j<contourMat.cols; j++)
+ {
+ if (i==j)
+ {
+ angleMatrix.at<float>(i,j)=0.0;
+ }
+ else
+ {
+ cv::Point2f dif = contourMat.at<cv::Point2f>(0,i) - contourMat.at<cv::Point2f>(0,j);
+ angleMatrix.at<float>(i,j) = std::atan2(dif.y, dif.x);
+
+ if (rotationInvariant)
+ {
+ cv::Point2f refPt = contourMat.at<cv::Point2f>(0,i) - massCenter;
+ float refAngle = atan2(refPt.y, refPt.x);
+ angleMatrix.at<float>(i,j) -= refAngle;
+ }
+ angleMatrix.at<float>(i,j) = fmod(angleMatrix.at<float>(i,j)+FLT_EPSILON,2*CV_PI)+CV_PI;
+ //angleMatrix.at<float>(i,j) = 1+floor( angleMatrix.at<float>(i,j)*nAngularBins/(2*CV_PI) );
+ }
+ }
+ }
+ }
+};
+
+/*
+ * Matcher
+ */
+class SCDMatcher
+{
+public:
+ // the full constructor
+ SCDMatcher()
+ {
+ }
+
+ // the matcher function using Hungarian method
+ void matchDescriptors(cv::Mat& descriptors1, cv::Mat& descriptors2, std::vector<cv::DMatch>& matches, cv::Ptr<cv::HistogramCostExtractor>& comparer,
+ std::vector<int>& inliers1, std::vector<int> &inliers2)
+ {
+ matches.clear();
+
+ // Build the cost Matrix between descriptors //
+ cv::Mat costMat;
+ buildCostMatrix(descriptors1, descriptors2, costMat, comparer);
+
+ // Solve the matching problem using the hungarian method //
+ hungarian(costMat, matches, inliers1, inliers2, descriptors1.rows, descriptors2.rows);
+ }
+
+ // matching cost
+ float getMatchingCost() const {return minMatchCost;}
+
+private:
+ float minMatchCost;
+ float betaAdditional;
+protected:
+ void buildCostMatrix(const cv::Mat& descriptors1, const cv::Mat& descriptors2,
+ cv::Mat& costMatrix, cv::Ptr<cv::HistogramCostExtractor>& comparer) const
+ {
+ comparer->buildCostMatrix(descriptors1, descriptors2, costMatrix);
+ }
+
+ void hungarian(cv::Mat& costMatrix, std::vector<cv::DMatch>& outMatches, std::vector<int> &inliers1,
+ std::vector<int> &inliers2, int sizeScd1=0, int sizeScd2=0)
+ {
+ std::vector<int> free(costMatrix.rows, 0), collist(costMatrix.rows, 0);
+ std::vector<int> matches(costMatrix.rows, 0), colsol(costMatrix.rows), rowsol(costMatrix.rows);
+ std::vector<float> d(costMatrix.rows), pred(costMatrix.rows), v(costMatrix.rows);
+
+ const float LOWV=1e-10;
+ bool unassignedfound;
+ int i=0, imin=0, numfree=0, prvnumfree=0, f=0, i0=0, k=0, freerow=0;
+ int j=0, j1=0, j2=0, endofpath=0, last=0, low=0, up=0;
+ float min=0, h=0, umin=0, usubmin=0, v2=0;
+
+ // COLUMN REDUCTION //
+ for (j = costMatrix.rows-1; j >= 0; j--)
+ {
+ // find minimum cost over rows.
+ min = costMatrix.at<float>(0,j);
+ imin = 0;
+ for (i = 1; i < costMatrix.rows; i++)
+ if (costMatrix.at<float>(i,j) < min)
+ {
+ min = costMatrix.at<float>(i,j);
+ imin = i;
+ }
+ v[j] = min;
+
+ if (++matches[imin] == 1)
+ {
+ rowsol[imin] = j;
+ colsol[j] = imin;
+ }
+ else
+ {
+ colsol[j]=-1;
+ }
+ }
+
+ // REDUCTION TRANSFER //
+ for (i=0; i<costMatrix.rows; i++)
+ {
+ if (matches[i] == 0)
+ {
+ free[numfree++] = i;
+ }
+ else
+ {
+ if (matches[i] == 1)
+ {
+ j1=rowsol[i];
+ min=std::numeric_limits<float>::max();
+ for (j=0; j<costMatrix.rows; j++)
+ {
+ if (j!=j1)
+ {
+ if (costMatrix.at<float>(i,j)-v[j] < min)
+ {
+ min=costMatrix.at<float>(i,j)-v[j];
+ }
+ }
+ }
+ v[j1] = v[j1]-min;
+ }
+ }
+ }
+ // AUGMENTING ROW REDUCTION //
+ int loopcnt = 0;
+ do
+ {
+ loopcnt++;
+ k=0;
+ prvnumfree=numfree;
+ numfree=0;
+ while (k < prvnumfree)
+ {
+ i=free[k];
+ k++;
+ umin = costMatrix.at<float>(i,0)-v[0];
+ j1=0;
+ usubmin = std::numeric_limits<float>::max();
+ for (j=1; j<costMatrix.rows; j++)
+ {
+ h = costMatrix.at<float>(i,j)-v[j];
+ if (h < usubmin)
+ {
+ if (h >= umin)
+ {
+ usubmin = h;
+ j2 = j;
+ }
+ else
+ {
+ usubmin = umin;
+ umin = h;
+ j2 = j1;
+ j1 = j;
+ }
+ }
+ }
+ i0 = colsol[j1];
+
+ if (fabs(umin-usubmin) > LOWV) //if( umin < usubmin )
+ {
+ v[j1] = v[j1] - (usubmin - umin);
+ }
+ else // minimum and subminimum equal.
+ {
+ if (i0 >= 0) // minimum column j1 is assigned.
+ {
+ j1 = j2;
+ i0 = colsol[j2];
+ }
+ }
+ // (re-)assign i to j1, possibly de-assigning an i0.
+ rowsol[i]=j1;
+ colsol[j1]=i;
+
+ if (i0 >= 0)
+ {
+ //if( umin < usubmin )
+ if (fabs(umin-usubmin) > LOWV)
+ {
+ free[--k] = i0;
+ }
+ else
+ {
+ free[numfree++] = i0;
+ }
+ }
+ }
+ }while (loopcnt<2); // repeat once.
+
+ // AUGMENT SOLUTION for each free row //
+ for (f = 0; f<numfree; f++)
+ {
+ freerow = free[f]; // start row of augmenting path.
+ // Dijkstra shortest path algorithm.
+ // runs until unassigned column added to shortest path tree.
+ for (j = 0; j < costMatrix.rows; j++)
+ {
+ d[j] = costMatrix.at<float>(freerow,j) - v[j];
+ pred[j] = freerow;
+ collist[j] = j; // init column list.
+ }
+
+ low=0; // columns in 0..low-1 are ready, now none.
+ up=0; // columns in low..up-1 are to be scanned for current minimum, now none.
+ unassignedfound = false;
+ do
+ {
+ if (up == low)
+ {
+ last=low-1;
+ min = d[collist[up++]];
+ for (k = up; k < costMatrix.rows; k++)
+ {
+ j = collist[k];
+ h = d[j];
+ if (h <= min)
+ {
+ if (h < min) // new minimum.
+ {
+ up = low; // restart list at index low.
+ min = h;
+ }
+ collist[k] = collist[up];
+ collist[up++] = j;
+ }
+ }
+ for (k=low; k<up; k++)
+ {
+ if (colsol[collist[k]] < 0)
+ {
+ endofpath = collist[k];
+ unassignedfound = true;
+ break;
+ }
+ }
+ }
+
+ if (!unassignedfound)
+ {
+ // update 'distances' between freerow and all unscanned columns, via next scanned column.
+ j1 = collist[low];
+ low++;
+ i = colsol[j1];
+ h = costMatrix.at<float>(i,j1)-v[j1]-min;
+
+ for (k = up; k < costMatrix.rows; k++)
+ {
+ j = collist[k];
+ v2 = costMatrix.at<float>(i,j) - v[j] - h;
+ if (v2 < d[j])
+ {
+ pred[j] = i;
+ if (v2 == min)
+ {
+ if (colsol[j] < 0)
+ {
+ // if unassigned, shortest augmenting path is complete.
+ endofpath = j;
+ unassignedfound = true;
+ break;
+ }
+ else
+ {
+ collist[k] = collist[up];
+ collist[up++] = j;
+ }
+ }
+ d[j] = v2;
+ }
+ }
+ }
+ }while (!unassignedfound);
+
+ // update column prices.
+ for (k = 0; k <= last; k++)
+ {
+ j1 = collist[k];
+ v[j1] = v[j1] + d[j1] - min;
+ }
+
+ // reset row and column assignments along the alternating path.
+ do
+ {
+ i = pred[endofpath];
+ colsol[endofpath] = i;
+ j1 = endofpath;
+ endofpath = rowsol[i];
+ rowsol[i] = j1;
+ }while (i != freerow);
+ }
+
+ // calculate symmetric shape context cost
+ cv::Mat trueCostMatrix(costMatrix, cv::Rect(0,0,sizeScd1, sizeScd2));
+ float leftcost = 0;
+ for (int nrow=0; nrow<trueCostMatrix.rows; nrow++)
+ {
+ double minval;
+ minMaxIdx(trueCostMatrix.row(nrow), &minval);
+ leftcost+=minval;
+ }
+ leftcost /= trueCostMatrix.rows;
+
+ float rightcost = 0;
+ for (int ncol=0; ncol<trueCostMatrix.cols; ncol++)
+ {
+ double minval;
+ minMaxIdx(trueCostMatrix.col(ncol), &minval);
+ rightcost+=minval;
+ }
+ rightcost /= trueCostMatrix.cols;
+
+ minMatchCost = std::max(leftcost,rightcost);
+
+ // Save in a DMatch vector
+ for (i=0;i<costMatrix.cols;i++)
+ {
+ cv::DMatch singleMatch(colsol[i],i,costMatrix.at<float>(colsol[i],i));//queryIdx,trainIdx,distance
+ outMatches.push_back(singleMatch);
+ }
+
+ // Update inliers
+ inliers1.reserve(sizeScd1);
+ for (size_t kc = 0; kc<inliers1.size(); kc++)
+ {
+ if (rowsol[kc]<sizeScd1) // if a real match
+ inliers1[kc]=1;
+ else
+ inliers1[kc]=0;
+ }
+ inliers2.reserve(sizeScd2);
+ for (size_t kc = 0; kc<inliers2.size(); kc++)
+ {
+ if (colsol[kc]<sizeScd2) // if a real match
+ inliers2[kc]=1;
+ else
+ inliers2[kc]=0;
+ }
+ }
+
+};
+
+/*
+ *
+ */
+
+namespace cv
+{
+class ShapeContextDistanceExtractorImpl : public ShapeContextDistanceExtractor
+{
+public:
+ /* Constructors */
+ ShapeContextDistanceExtractorImpl(int _nAngularBins, int _nRadialBins, float _innerRadius, float _outerRadius, int _iterations,
+ const Ptr<HistogramCostExtractor> &_comparer, const Ptr<ShapeTransformer> &_transformer)
+ {
+ nAngularBins=_nAngularBins;
+ nRadialBins=_nRadialBins;
+ innerRadius=_innerRadius;
+ outerRadius=_outerRadius;
+ rotationInvariant=false;
+ comparer=_comparer;
+ iterations=_iterations;
+ transformer=_transformer;
+ bendingEnergyWeight=0.3;
+ imageAppearanceWeight=0.0;
+ shapeContextWeight=1.0;
+ sigma=10;
+ name_ = "ShapeDistanceExtractor.SCD";
+ }
+
+ /* Destructor */
+ ~ShapeContextDistanceExtractorImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operator
+ virtual float computeDistance(InputArray contour1, InputArray contour2);
+
+ //! Setters/Getters
+ virtual void setAngularBins(int _nAngularBins){CV_Assert(_nAngularBins>0); nAngularBins=_nAngularBins;}
+ virtual int getAngularBins() const {return nAngularBins;}
+
+ virtual void setRadialBins(int _nRadialBins){CV_Assert(_nRadialBins>0); nRadialBins=_nRadialBins;}
+ virtual int getRadialBins() const {return nRadialBins;}
+
+ virtual void setInnerRadius(float _innerRadius) {CV_Assert(_innerRadius>0); innerRadius=_innerRadius;}
+ virtual float getInnerRadius() const {return innerRadius;}
+
+ virtual void setOuterRadius(float _outerRadius) {CV_Assert(_outerRadius>0); outerRadius=_outerRadius;}
+ virtual float getOuterRadius() const {return outerRadius;}
+
+ virtual void setRotationInvariant(bool _rotationInvariant) {rotationInvariant=_rotationInvariant;}
+ virtual bool getRotationInvariant() const {return rotationInvariant;}
+
+ virtual void setCostExtractor(Ptr<HistogramCostExtractor> _comparer) { comparer = _comparer; }
+ virtual Ptr<HistogramCostExtractor> getCostExtractor() const { return comparer; }
+
+ virtual void setShapeContextWeight(float _shapeContextWeight) {shapeContextWeight=_shapeContextWeight;}
+ virtual float getShapeContextWeight() const {return shapeContextWeight;}
+
+ virtual void setImageAppearanceWeight(float _imageAppearanceWeight) {imageAppearanceWeight=_imageAppearanceWeight;}
+ virtual float getImageAppearanceWeight() const {return imageAppearanceWeight;}
+
+ virtual void setBendingEnergyWeight(float _bendingEnergyWeight) {bendingEnergyWeight=_bendingEnergyWeight;}
+ virtual float getBendingEnergyWeight() const {return bendingEnergyWeight;}
+
+ virtual void setStdDev(float _sigma) {sigma=_sigma;}
+ virtual float getStdDev() const {return sigma;}
+
+ virtual void setImages(InputArray _image1, InputArray _image2)
+ {
+ Mat image1_=_image1.getMat(), image2_=_image2.getMat();
+ CV_Assert((image1_.depth()==0) & (image2_.depth()==0));
+ image1=image1_;
+ image2=image2_;
+ }
+
+ virtual void getImages(OutputArray _image1, OutputArray _image2) const
+ {
+ CV_Assert((!image1.empty()) & (!image2.empty()));
+ _image1.create(image1.size(), image1.type());
+ _image2.create(image2.size(), image2.type());
+ _image1.getMat()=image1;
+ _image2.getMat()=image2;
+ }
+
+ virtual void setIterations(int _iterations) {CV_Assert(_iterations>0); iterations=_iterations;}
+ virtual int getIterations() const {return iterations;}
+
+ virtual void setTransformAlgorithm(Ptr<ShapeTransformer> _transformer) {transformer=_transformer;}
+ virtual Ptr<ShapeTransformer> getTransformAlgorithm() const {return transformer;}
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "nRads" << nRadialBins
+ << "nAngs" << nAngularBins
+ << "iters" << iterations
+ << "img_1" << image1
+ << "img_2" << image2
+ << "beWei" << bendingEnergyWeight
+ << "scWei" << shapeContextWeight
+ << "iaWei" << imageAppearanceWeight
+ << "costF" << costFlag
+ << "rotIn" << rotationInvariant
+ << "sigma" << sigma;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ nRadialBins = (int)fn["nRads"];
+ nAngularBins = (int)fn["nAngs"];
+ iterations = (int)fn["iters"];
+ bendingEnergyWeight = (float)fn["beWei"];
+ shapeContextWeight = (float)fn["scWei"];
+ imageAppearanceWeight = (float)fn["iaWei"];
+ costFlag = (int)fn["costF"];
+ sigma = (float)fn["sigma"];
+ }
+
+private:
+ int nAngularBins;
+ int nRadialBins;
+ float innerRadius;
+ float outerRadius;
+ bool rotationInvariant;
+ int costFlag;
+ int iterations;
+ Ptr<ShapeTransformer> transformer;
+ Ptr<HistogramCostExtractor> comparer;
+ Mat image1;
+ Mat image2;
+ float bendingEnergyWeight;
+ float imageAppearanceWeight;
+ float shapeContextWeight;
+ float sigma;
+
+protected:
+ String name_;
+};
+
+float ShapeContextDistanceExtractorImpl::computeDistance(InputArray contour1, InputArray contour2)
+{
+ // Checking //
+ Mat sset1=contour1.getMat(), sset2=contour2.getMat(), set1, set2;
+ if (set1.type() != CV_32F)
+ sset1.convertTo(set1, CV_32F);
+ else
+ sset1.copyTo(set1);
+
+ if (set2.type() != CV_32F)
+ sset2.convertTo(set2, CV_32F);
+ else
+ sset1.copyTo(set2);
+
+ CV_Assert((set1.channels()==2) & (set1.cols>0));
+ CV_Assert((set2.channels()==2) & (set2.cols>0));
+ if (imageAppearanceWeight!=0)
+ {
+ CV_Assert((!image1.empty()) & (!image2.empty()));
+ }
+
+ // Initializing Extractor, Descriptor structures and Matcher //
+ SCD set1SCE(nAngularBins, nRadialBins, innerRadius, outerRadius, false);
+ Mat set1SCD;
+ SCD set2SCE(nAngularBins, nRadialBins, innerRadius, outerRadius, false);
+ Mat set2SCD;
+ SCDMatcher matcher;
+ std::vector<DMatch> matches;
+
+ // Distance components (The output is a linear combination of these 3) //
+ float sDistance=0, bEnergy=0, iAppearance=0;
+ float beta;
+
+ // Initializing some variables //
+ std::vector<int> inliers1, inliers2;
+ bool isTPS=false;
+ if ( dynamic_cast<ThinPlateSplineShapeTransformer*>(&*transformer) )
+ isTPS=true;
+ Mat warpedImage;
+ for (int ii=0; ii<iterations; ii++)
+ {
+ // Extract SCD descriptor in the set1 //
+ set1SCE.extractSCD(set1, set1SCD, inliers1);
+
+ // Extract SCD descriptor of the set2 (TARGET) //
+ set2SCE.extractSCD(set2, set2SCD, inliers2, set1SCE.getMeanDistance());
+
+ // regularization parameter with annealing rate annRate //
+ beta=std::pow(set1SCE.getMeanDistance(),2);
+
+ // match //
+ matcher.matchDescriptors(set1SCD, set2SCD, matches, comparer, inliers1, inliers2);
+
+ // apply TPS transform //
+ if ( isTPS )
+ dynamic_cast<ThinPlateSplineShapeTransformer*>(&*transformer)->setRegularizationParameter(beta);
+ transformer->estimateTransformation(set1, set2, matches);
+ bEnergy += transformer->applyTransformation(set1, set1);
+
+ // Image appearance //
+ if (imageAppearanceWeight!=0)
+ {
+ // Have to accumulate the transformation along all the iterations
+ if (ii==0)
+ {
+ if ( isTPS )
+ {
+ image2.copyTo(warpedImage);
+ }
+ else
+ {
+ image1.copyTo(warpedImage);
+ }
+ }
+ transformer->warpImage(warpedImage, warpedImage);
+ }
+ }
+
+ Mat gaussWindow, diffIm;
+ if (imageAppearanceWeight!=0)
+ {
+ // compute appearance cost
+ if ( isTPS )
+ {
+ resize(warpedImage, warpedImage, image1.size());
+ Mat temp=(warpedImage-image1);
+ multiply(temp, temp, diffIm);
+ }
+ else
+ {
+ resize(warpedImage, warpedImage, image2.size());
+ Mat temp=(warpedImage-image2);
+ multiply(temp, temp, diffIm);
+ }
+ gaussWindow = Mat::zeros(warpedImage.rows, warpedImage.cols, CV_32F);
+ for (int pt=0; pt<sset1.cols; pt++)
+ {
+ for (int ii=0; ii<diffIm.rows; ii++)
+ {
+ for (int jj=0; jj<diffIm.cols; jj++)
+ {
+ float xx = sset1.at<Point2f>(0,pt).x;
+ float yy = sset1.at<Point2f>(0,pt).y;
+ float val = std::exp( -( (xx-jj)*(xx-jj) + (yy-ii)*(yy-ii) )/(2*sigma*sigma) ) / (sigma*sigma*2*CV_PI);
+ gaussWindow.at<float>(ii,jj) += val;
+ }
+ }
+ }
+
+ Mat appIm(diffIm.rows, diffIm.cols, CV_32F);
+ for (int ii=0; ii<diffIm.rows; ii++)
+ {
+ for (int jj=0; jj<diffIm.cols; jj++)
+ {
+ float elema=float( diffIm.at<uchar>(ii,jj) )/255;
+ float elemb=gaussWindow.at<float>(ii,jj);
+ appIm.at<float>(ii,jj) = elema*elemb;
+ }
+ }
+ iAppearance = cv::sum(appIm)[0]/sset1.cols;
+ }
+ sDistance = matcher.getMatchingCost();
+
+ return (sDistance*shapeContextWeight+bEnergy*bendingEnergyWeight+iAppearance*imageAppearanceWeight);
+}
+
+Ptr <ShapeContextDistanceExtractor> createShapeContextDistanceExtractor(int nAngularBins, int nRadialBins, float innerRadius, float outerRadius, int iterations,
+ const Ptr<HistogramCostExtractor> &comparer, const Ptr<ShapeTransformer> &transformer)
+{
+ return Ptr <ShapeContextDistanceExtractor> ( new ShapeContextDistanceExtractorImpl(nAngularBins, nRadialBins, innerRadius,
+ outerRadius, iterations, comparer, transformer) );
+}
+
+} // cv
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
+// Copyright (C) 2009, Willow Garage Inc., all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of the copyright holders may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "precomp.hpp"
+
+namespace cv
+{
+
+class ThinPlateSplineShapeTransformerImpl : public ThinPlateSplineShapeTransformer
+{
+public:
+ /* Constructors */
+ ThinPlateSplineShapeTransformerImpl()
+ {
+ regularizationParameter=0;
+ name_ = "ShapeTransformer.TPS";
+ tpsComputed=false;
+ }
+
+ ThinPlateSplineShapeTransformerImpl(double _regularizationParameter)
+ {
+ regularizationParameter=_regularizationParameter;
+ name_ = "ShapeTransformer.TPS";
+ tpsComputed=false;
+ }
+
+ /* Destructor */
+ ~ThinPlateSplineShapeTransformerImpl()
+ {
+ }
+
+ virtual AlgorithmInfo* info() const { return 0; }
+
+ //! the main operators
+ virtual void estimateTransformation(InputArray transformingShape, InputArray targetShape, std::vector<DMatch> &matches);
+ virtual float applyTransformation(InputArray inPts, OutputArray output=noArray());
+ virtual void warpImage(InputArray transformingImage, OutputArray output,
+ int flags, int borderMode, const Scalar& borderValue) const;
+
+ //! Setters/Getters
+ virtual void setRegularizationParameter(double _regularizationParameter) {regularizationParameter=_regularizationParameter;}
+ virtual double getRegularizationParameter() const {return regularizationParameter;}
+
+ //! write/read
+ virtual void write(FileStorage& fs) const
+ {
+ fs << "name" << name_
+ << "regularization" << regularizationParameter;
+ }
+
+ virtual void read(const FileNode& fn)
+ {
+ CV_Assert( (String)fn["name"] == name_ );
+ regularizationParameter = (int)fn["regularization"];
+ }
+
+private:
+ bool tpsComputed;
+ double regularizationParameter;
+ float transformCost;
+ Mat tpsParameters;
+ Mat shapeReference;
+
+protected:
+ String name_;
+};
+
+static double distance(Point2f p, Point2f q)
+{
+ Point2f diff = p - q;
+ float norma = diff.x*diff.x + diff.y*diff.y;// - 2*diff.x*diff.y;
+ if (norma<0) norma=0;
+ //else norma = std::sqrt(norma);
+ norma = norma*std::log(norma+FLT_EPSILON);
+ return norma;
+}
+
+static Point2f _applyTransformation(const Mat &shapeRef, const Point2f point, const Mat &tpsParameters)
+{
+ Point2f out;
+ for (int i=0; i<2; i++)
+ {
+ float a1=tpsParameters.at<float>(tpsParameters.rows-3,i);
+ float ax=tpsParameters.at<float>(tpsParameters.rows-2,i);
+ float ay=tpsParameters.at<float>(tpsParameters.rows-1,i);
+
+ float affine=a1+ax*point.x+ay*point.y;
+ float nonrigid=0;
+ for (int j=0; j<shapeRef.rows; j++)
+ {
+ nonrigid+=tpsParameters.at<float>(j,i)*
+ distance(Point2f(shapeRef.at<float>(j,0),shapeRef.at<float>(j,1)),
+ point);
+ }
+ if (i==0)
+ {
+ out.x=affine+nonrigid;
+ }
+ if (i==1)
+ {
+ out.y=affine+nonrigid;
+ }
+ }
+ return out;
+}
+
+/* public methods */
+void ThinPlateSplineShapeTransformerImpl::warpImage(InputArray transformingImage, OutputArray output,
+ int flags, int borderMode, const Scalar& borderValue) const
+{
+ CV_Assert(tpsComputed==true);
+
+ Mat theinput = transformingImage.getMat();
+ Mat mapX(theinput.rows, theinput.cols, CV_32FC1);
+ Mat mapY(theinput.rows, theinput.cols, CV_32FC1);
+
+ for (int row = 0; row < theinput.rows; row++)
+ {
+ for (int col = 0; col < theinput.cols; col++)
+ {
+ Point2f pt = _applyTransformation(shapeReference, Point2f(float(col), float(row)), tpsParameters);
+ mapX.at<float>(row, col) = pt.x;
+ mapY.at<float>(row, col) = pt.y;
+ }
+ }
+ remap(transformingImage, output, mapX, mapY, flags, borderMode, borderValue);
+}
+
+float ThinPlateSplineShapeTransformerImpl::applyTransformation(InputArray inPts, OutputArray outPts)
+{
+ CV_Assert(tpsComputed);
+ Mat pts1 = inPts.getMat();
+ CV_Assert((pts1.channels()==2) & (pts1.cols>0));
+
+ //Apply transformation in the complete set of points
+ // Ensambling output //
+ if (outPts.needed())
+ {
+ outPts.create(1,pts1.cols, CV_32FC2);
+ Mat outMat = outPts.getMat();
+ for (int i=0; i<pts1.cols; i++)
+ {
+ Point2f pt=pts1.at<Point2f>(0,i);
+ outMat.at<Point2f>(0,i)=_applyTransformation(shapeReference, pt, tpsParameters);
+ }
+ }
+
+ return transformCost;
+}
+
+void ThinPlateSplineShapeTransformerImpl::estimateTransformation(InputArray _pts1, InputArray _pts2,
+ std::vector<DMatch>& _matches )
+{
+ Mat pts1 = _pts1.getMat();
+ Mat pts2 = _pts2.getMat();
+ CV_Assert((pts1.channels()==2) & (pts1.cols>0) & (pts2.channels()==2) & (pts2.cols>0));
+ CV_Assert(_matches.size()>1);
+
+ if (pts1.type() != CV_32F)
+ pts1.convertTo(pts1, CV_32F);
+ if (pts2.type() != CV_32F)
+ pts2.convertTo(pts2, CV_32F);
+
+ // Use only valid matchings //
+ std::vector<DMatch> matches;
+ for (size_t i=0; i<_matches.size(); i++)
+ {
+ if (_matches[i].queryIdx<pts1.cols &&
+ _matches[i].trainIdx<pts2.cols)
+ {
+ matches.push_back(_matches[i]);
+ }
+ }
+
+ // Organizing the correspondent points in matrix style //
+ Mat shape1(matches.size(),2,CV_32F); // transforming shape
+ Mat shape2(matches.size(),2,CV_32F); // target shape
+ for (size_t i=0; i<matches.size(); i++)
+ {
+ Point2f pt1=pts1.at<Point2f>(0,matches[i].queryIdx);
+ shape1.at<float>(i,0) = pt1.x;
+ shape1.at<float>(i,1) = pt1.y;
+
+ Point2f pt2=pts2.at<Point2f>(0,matches[i].trainIdx);
+ shape2.at<float>(i,0) = pt2.x;
+ shape2.at<float>(i,1) = pt2.y;
+ }
+ shape1.copyTo(shapeReference);
+
+ // Building the matrices for solving the L*(w|a)=(v|0) problem with L={[K|P];[P'|0]}
+
+ //Building K and P (Neede to buil L)
+ Mat matK(matches.size(),matches.size(),CV_32F);
+ Mat matP(matches.size(),3,CV_32F);
+ for (size_t i=0; i<matches.size(); i++)
+ {
+ for (size_t j=0; j<matches.size(); j++)
+ {
+ if (i==j)
+ {
+ matK.at<float>(i,j)=regularizationParameter;
+ }
+ else
+ {
+ matK.at<float>(i,j) = distance(Point2f(shape1.at<float>(i,0),shape1.at<float>(i,1)),
+ Point2f(shape1.at<float>(j,0),shape1.at<float>(j,1)));
+ }
+ }
+ matP.at<float>(i,0) = 1;
+ matP.at<float>(i,1) = shape1.at<float>(i,0);
+ matP.at<float>(i,2) = shape1.at<float>(i,1);
+ }
+
+ //Building L
+ Mat matL=Mat::zeros(matches.size()+3,matches.size()+3,CV_32F);
+ Mat matLroi(matL, Rect(0,0,matches.size(),matches.size())); //roi for K
+ matK.copyTo(matLroi);
+ matLroi = Mat(matL,Rect(matches.size(),0,3,matches.size())); //roi for P
+ matP.copyTo(matLroi);
+ Mat matPt;
+ transpose(matP,matPt);
+ matLroi = Mat(matL,Rect(0,matches.size(),matches.size(),3)); //roi for P'
+ matPt.copyTo(matLroi);
+
+ //Building B (v|0)
+ Mat matB = Mat::zeros(matches.size()+3,2,CV_32F);
+ for (size_t i=0; i<matches.size(); i++)
+ {
+ matB.at<float>(i,0) = shape2.at<float>(i,0); //x's
+ matB.at<float>(i,1) = shape2.at<float>(i,1); //y's
+ }
+
+ //Obtaining transformation params (w|a)
+ solve(matL, matB, tpsParameters, DECOMP_LU);
+ //tpsParameters = matL.inv()*matB;
+
+ //Setting transform Cost and Shape reference
+ Mat w(tpsParameters, Rect(0,0,2,tpsParameters.rows-3));
+ Mat Q=w.t()*matK*w;
+ transformCost=fabs(Q.at<float>(0,0)*Q.at<float>(1,1));//fabs(mean(Q.diag(0))[0]);//std::max(Q.at<float>(0,0),Q.at<float>(1,1));
+ tpsComputed=true;
+}
+
+Ptr <ThinPlateSplineShapeTransformer> createThinPlateSplineShapeTransformer(double regularizationParameter)
+{
+ return Ptr<ThinPlateSplineShapeTransformer>( new ThinPlateSplineShapeTransformerImpl(regularizationParameter) );
+}
+
+} // cv
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// Intel License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000, Intel Corporation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of Intel Corporation may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "test_precomp.hpp"
+
+using namespace cv;
+using namespace std;
+
+const int angularBins=12;
+const int radialBins=4;
+const float minRad=0.2;
+const float maxRad=2;
+const int NSN=5;//10;//20; //number of shapes per class
+const int NP=100; //number of points sympliying the contour
+const float outlierWeight=0.1;
+const int numOutliers=20;
+const float CURRENT_MAX_ACCUR=95; //98% and 99% reached in several tests, 95 is fixed as minimum boundary
+
+class CV_ShapeEMDTest : public cvtest::BaseTest
+{
+public:
+ CV_ShapeEMDTest();
+ ~CV_ShapeEMDTest();
+protected:
+ void run(int);
+
+private:
+ void mpegTest();
+ void listShapeNames(vector<string> &listHeaders);
+ vector<Point2f> convertContourType(const Mat &, int n=0 );
+ float computeShapeDistance(vector <Point2f>& queryNormal,
+ vector <Point2f>& queryFlipped1,
+ vector <Point2f>& queryFlipped2,
+ vector<Point2f>& testq);
+ void displayMPEGResults();
+};
+
+CV_ShapeEMDTest::CV_ShapeEMDTest()
+{
+}
+CV_ShapeEMDTest::~CV_ShapeEMDTest()
+{
+}
+
+vector <Point2f> CV_ShapeEMDTest::convertContourType(const Mat& currentQuery, int n)
+{
+ vector<vector<Point> > _contoursQuery;
+ vector <Point2f> contoursQuery;
+ findContours(currentQuery, _contoursQuery, RETR_LIST, CHAIN_APPROX_NONE);
+ for (size_t border=0; border<_contoursQuery.size(); border++)
+ {
+ for (size_t p=0; p<_contoursQuery[border].size(); p++)
+ {
+ contoursQuery.push_back(Point2f((float)_contoursQuery[border][p].x,
+ (float)_contoursQuery[border][p].y));
+ }
+ }
+
+ // In case actual number of points is less than n
+ int dum=0;
+ for (int add=contoursQuery.size()-1; add<n; add++)
+ {
+ contoursQuery.push_back(contoursQuery[dum++]); //adding dummy values
+ }
+
+ // Uniformly sampling
+ random_shuffle(contoursQuery.begin(), contoursQuery.end());
+ int nStart=n;
+ vector<Point2f> cont;
+ for (int i=0; i<nStart; i++)
+ {
+ cont.push_back(contoursQuery[i]);
+ }
+ return cont;
+}
+
+void CV_ShapeEMDTest::listShapeNames( vector<string> &listHeaders)
+{
+ listHeaders.push_back("apple"); //ok
+ listHeaders.push_back("children"); // ok
+ listHeaders.push_back("device7"); // ok
+ listHeaders.push_back("Heart"); // ok
+ listHeaders.push_back("teddy"); // ok
+}
+float CV_ShapeEMDTest::computeShapeDistance(vector <Point2f>& query1, vector <Point2f>& query2,
+ vector <Point2f>& query3, vector <Point2f>& testq)
+{
+ //waitKey(0);
+ Ptr <ShapeContextDistanceExtractor> mysc = createShapeContextDistanceExtractor(angularBins, radialBins, minRad, maxRad);
+ //Ptr <HistogramCostExtractor> cost = createNormHistogramCostExtractor(cv::DIST_L1);
+ //Ptr <HistogramCostExtractor> cost = createChiHistogramCostExtractor(30,0.15);
+ //Ptr <HistogramCostExtractor> cost = createEMDHistogramCostExtractor();
+ // Ptr <HistogramCostExtractor> cost = createEMDL1HistogramCostExtractor();
+ mysc->setIterations(1); //(3)
+ mysc->setCostExtractor( createEMDL1HistogramCostExtractor() );
+ //mysc->setTransformAlgorithm(createAffineTransformer(true));
+ mysc->setTransformAlgorithm( createThinPlateSplineShapeTransformer() );
+ //mysc->setImageAppearanceWeight(1.6);
+ //mysc->setImageAppearanceWeight(0.0);
+ //mysc->setImages(im1,imtest);
+ return ( std::min( mysc->computeDistance(query1, testq),
+ std::min(mysc->computeDistance(query2, testq), mysc->computeDistance(query3, testq) )));
+}
+
+void CV_ShapeEMDTest::mpegTest()
+{
+ string baseTestFolder="shape/mpeg_test/";
+ string path = cvtest::TS::ptr()->get_data_path() + baseTestFolder;
+ vector<string> namesHeaders;
+ listShapeNames(namesHeaders);
+
+ // distance matrix //
+ Mat distanceMat=Mat::zeros(NSN*namesHeaders.size(), NSN*namesHeaders.size(), CV_32F);
+
+ // query contours (normal v flipped, h flipped) and testing contour //
+ vector<Point2f> contoursQuery1, contoursQuery2, contoursQuery3, contoursTesting;
+
+ // reading query and computing its properties //
+ int counter=0;
+ const int loops=NSN*namesHeaders.size()*NSN*namesHeaders.size();
+ for (size_t n=0; n<namesHeaders.size(); n++)
+ {
+ for (int i=1; i<=NSN; i++)
+ {
+ // read current image //
+ stringstream thepathandname;
+ thepathandname<<path+namesHeaders[n]<<"-"<<i<<".png";
+ Mat currentQuery, flippedHQuery, flippedVQuery;
+ currentQuery=imread(thepathandname.str(), IMREAD_GRAYSCALE);
+ Mat currentQueryBuf=currentQuery.clone();
+ flip(currentQuery, flippedHQuery, 0);
+ flip(currentQuery, flippedVQuery, 1);
+ // compute border of the query and its flipped versions //
+ vector<Point2f> origContour;
+ contoursQuery1=convertContourType(currentQuery, NP);
+ origContour=contoursQuery1;
+ contoursQuery2=convertContourType(flippedHQuery, NP);
+ contoursQuery3=convertContourType(flippedVQuery, NP);
+
+ // compare with all the rest of the images: testing //
+ for (size_t nt=0; nt<namesHeaders.size(); nt++)
+ {
+ for (int it=1; it<=NSN; it++)
+ {
+ // skip self-comparisson //
+ counter++;
+ if (nt==n && it==i)
+ {
+ distanceMat.at<float>(NSN*n+i-1,
+ NSN*nt+it-1)=0;
+ continue;
+ }
+ // read testing image //
+ stringstream thetestpathandname;
+ thetestpathandname<<path+namesHeaders[nt]<<"-"<<it<<".png";
+ Mat currentTest;
+ currentTest=imread(thetestpathandname.str().c_str(), 0);
+ // compute border of the testing //
+ contoursTesting=convertContourType(currentTest, NP);
+
+ // compute shape distance //
+ std::cout<<std::endl<<"Progress: "<<counter<<"/"<<loops<<": "<<100*double(counter)/loops<<"% *******"<<std::endl;
+ std::cout<<"Computing shape distance between "<<namesHeaders[n]<<i<<
+ " and "<<namesHeaders[nt]<<it<<": ";
+ distanceMat.at<float>(NSN*n+i-1, NSN*nt+it-1)=
+ computeShapeDistance(contoursQuery1, contoursQuery2, contoursQuery3, contoursTesting);
+ std::cout<<distanceMat.at<float>(NSN*n+i-1, NSN*nt+it-1)<<std::endl;
+ }
+ }
+ }
+ }
+ // save distance matrix //
+ FileStorage fs(cvtest::TS::ptr()->get_data_path() + baseTestFolder + "distanceMatrixMPEGTest.yml", FileStorage::WRITE);
+ fs << "distanceMat" << distanceMat;
+}
+
+const int FIRST_MANY=2*NSN;
+void CV_ShapeEMDTest::displayMPEGResults()
+{
+ string baseTestFolder="shape/mpeg_test/";
+ Mat distanceMat;
+ FileStorage fs(cvtest::TS::ptr()->get_data_path() + baseTestFolder + "distanceMatrixMPEGTest.yml", FileStorage::READ);
+ vector<string> namesHeaders;
+ listShapeNames(namesHeaders);
+
+ // Read generated MAT //
+ fs["distanceMat"]>>distanceMat;
+
+ int corrects=0;
+ int divi=0;
+ for (int row=0; row<distanceMat.rows; row++)
+ {
+ if (row%NSN==0) //another group
+ {
+ divi+=NSN;
+ }
+ for (int col=divi-NSN; col<divi; col++)
+ {
+ int nsmall=0;
+ for (int i=0; i<distanceMat.cols; i++)
+ {
+ if (distanceMat.at<float>(row,col)>distanceMat.at<float>(row,i))
+ {
+ nsmall++;
+ }
+ }
+ if (nsmall<=FIRST_MANY)
+ {
+ corrects++;
+ }
+ }
+ }
+ float porc = 100*float(corrects)/(NSN*distanceMat.rows);
+ std::cout<<"%="<<porc<<std::endl;
+ if (porc >= CURRENT_MAX_ACCUR)
+ ts->set_failed_test_info(cvtest::TS::OK);
+ else
+ ts->set_failed_test_info(cvtest::TS::FAIL_BAD_ACCURACY);
+
+}
+
+void CV_ShapeEMDTest::run( int /*start_from*/ )
+{
+ mpegTest();
+ displayMPEGResults();
+}
+
+TEST(ShapeEMD_SCD, regression) { CV_ShapeEMDTest test; test.safe_run(); }
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// Intel License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000, Intel Corporation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of Intel Corporation may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "test_precomp.hpp"
+#include <stdlib.h>
+
+using namespace cv;
+using namespace std;
+
+const int NSN=5;//10;//20; //number of shapes per class
+const float CURRENT_MAX_ACCUR=85; //90% and 91% reached in several tests, 85 is fixed as minimum boundary
+
+class CV_HaussTest : public cvtest::BaseTest
+{
+public:
+ CV_HaussTest();
+ ~CV_HaussTest();
+protected:
+ void run(int);
+private:
+ float computeShapeDistance(vector<Point> &query1, vector<Point> &query2,
+ vector<Point> &query3, vector<Point> &testq);
+ vector <Point> convertContourType(const Mat& currentQuery, int n=180);
+ vector<Point2f> normalizeContour(const vector <Point>& contour);
+ void listShapeNames( vector<string> &listHeaders);
+ void mpegTest();
+ void displayMPEGResults();
+};
+
+CV_HaussTest::CV_HaussTest()
+{
+}
+CV_HaussTest::~CV_HaussTest()
+{
+}
+
+vector<Point2f> CV_HaussTest::normalizeContour(const vector<Point> &contour)
+{
+ vector<Point2f> output(contour.size());
+ Mat disMat(contour.size(),contour.size(),CV_32F);
+ Point2f meanpt(0,0);
+ float meanVal=1;
+
+ for (size_t ii=0; ii<contour.size(); ii++)
+ {
+ for (size_t jj=0; jj<contour.size(); jj++)
+ {
+ if (ii==jj) disMat.at<float>(ii,jj)=0;
+ else
+ {
+ disMat.at<float>(ii,jj)=
+ fabs(contour[ii].x*contour[jj].x)+fabs(contour[ii].y*contour[jj].y);
+ }
+ }
+ meanpt.x+=contour[ii].x;
+ meanpt.y+=contour[ii].y;
+ }
+ meanpt.x/=contour.size();
+ meanpt.y/=contour.size();
+ meanVal=cv::mean(disMat)[0];
+ for (size_t ii=0; ii<contour.size(); ii++)
+ {
+ output[ii].x = (contour[ii].x-meanpt.x)/meanVal;
+ output[ii].y = (contour[ii].y-meanpt.y)/meanVal;
+ }
+ return output;
+}
+
+void CV_HaussTest::listShapeNames( vector<string> &listHeaders)
+{
+ listHeaders.push_back("apple"); //ok
+ listHeaders.push_back("children"); // ok
+ listHeaders.push_back("device7"); // ok
+ listHeaders.push_back("Heart"); // ok
+ listHeaders.push_back("teddy"); // ok
+}
+
+
+vector <Point> CV_HaussTest::convertContourType(const Mat& currentQuery, int n)
+{
+ vector<vector<Point> > _contoursQuery;
+ vector <Point> contoursQuery;
+ findContours(currentQuery, _contoursQuery, RETR_LIST, CHAIN_APPROX_NONE);
+ for (size_t border=0; border<_contoursQuery.size(); border++)
+ {
+ for (size_t p=0; p<_contoursQuery[border].size(); p++)
+ {
+ contoursQuery.push_back(_contoursQuery[border][p]);
+ }
+ }
+
+ // In case actual number of points is less than n
+ for (int add=contoursQuery.size()-1; add<n; add++)
+ {
+ contoursQuery.push_back(contoursQuery[contoursQuery.size()-add+1]); //adding dummy values
+ }
+
+ // Uniformly sampling
+ random_shuffle(contoursQuery.begin(), contoursQuery.end());
+ int nStart=n;
+ vector<Point> cont;
+ for (int i=0; i<nStart; i++)
+ {
+ cont.push_back(contoursQuery[i]);
+ }
+ return cont;
+}
+
+float CV_HaussTest::computeShapeDistance(vector <Point>& query1, vector <Point>& query2,
+ vector <Point>& query3, vector <Point>& testq)
+{
+ Ptr <HausdorffDistanceExtractor> haus = createHausdorffDistanceExtractor();
+ return std::min(haus->computeDistance(query1,testq), std::min(haus->computeDistance(query2,testq),
+ haus->computeDistance(query3,testq)));
+}
+
+void CV_HaussTest::mpegTest()
+{
+ string baseTestFolder="shape/mpeg_test/";
+ string path = cvtest::TS::ptr()->get_data_path() + baseTestFolder;
+ vector<string> namesHeaders;
+ listShapeNames(namesHeaders);
+
+ // distance matrix //
+ Mat distanceMat=Mat::zeros(NSN*namesHeaders.size(), NSN*namesHeaders.size(), CV_32F);
+
+ // query contours (normal v flipped, h flipped) and testing contour //
+ vector<Point> contoursQuery1, contoursQuery2, contoursQuery3, contoursTesting;
+
+ // reading query and computing its properties //
+ int counter=0;
+ const int loops=NSN*namesHeaders.size()*NSN*namesHeaders.size();
+ for (size_t n=0; n<namesHeaders.size(); n++)
+ {
+ for (int i=1; i<=NSN; i++)
+ {
+ // read current image //
+ stringstream thepathandname;
+ thepathandname<<path+namesHeaders[n]<<"-"<<i<<".png";
+ Mat currentQuery, flippedHQuery, flippedVQuery;
+ currentQuery=imread(thepathandname.str(), IMREAD_GRAYSCALE);
+ flip(currentQuery, flippedHQuery, 0);
+ flip(currentQuery, flippedVQuery, 1);
+ // compute border of the query and its flipped versions //
+ vector<Point> origContour;
+ contoursQuery1=convertContourType(currentQuery);
+ origContour=contoursQuery1;
+ contoursQuery2=convertContourType(flippedHQuery);
+ contoursQuery3=convertContourType(flippedVQuery);
+
+ // compare with all the rest of the images: testing //
+ for (size_t nt=0; nt<namesHeaders.size(); nt++)
+ {
+ for (int it=1; it<=NSN; it++)
+ {
+ /* skip self-comparisson */
+ counter++;
+ if (nt==n && it==i)
+ {
+ distanceMat.at<float>(NSN*n+i-1,
+ NSN*nt+it-1)=0;
+ continue;
+ }
+ // read testing image //
+ stringstream thetestpathandname;
+ thetestpathandname<<path+namesHeaders[nt]<<"-"<<it<<".png";
+ Mat currentTest;
+ currentTest=imread(thetestpathandname.str().c_str(), 0);
+
+ // compute border of the testing //
+ contoursTesting=convertContourType(currentTest);
+
+ // compute shape distance //
+ std::cout<<std::endl<<"Progress: "<<counter<<"/"<<loops<<": "<<100*double(counter)/loops<<"% *******"<<std::endl;
+ std::cout<<"Computing shape distance between "<<namesHeaders[n]<<i<<
+ " and "<<namesHeaders[nt]<<it<<": ";
+ distanceMat.at<float>(NSN*n+i-1, NSN*nt+it-1)=
+ computeShapeDistance(contoursQuery1, contoursQuery2, contoursQuery3, contoursTesting);
+ std::cout<<distanceMat.at<float>(NSN*n+i-1, NSN*nt+it-1)<<std::endl;
+ }
+ }
+ }
+ }
+ // save distance matrix //
+ FileStorage fs(cvtest::TS::ptr()->get_data_path() + baseTestFolder + "distanceMatrixMPEGTest.yml", FileStorage::WRITE);
+ fs << "distanceMat" << distanceMat;
+}
+
+const int FIRST_MANY=2*NSN;
+void CV_HaussTest::displayMPEGResults()
+{
+ string baseTestFolder="shape/mpeg_test/";
+ Mat distanceMat;
+ FileStorage fs(cvtest::TS::ptr()->get_data_path() + baseTestFolder + "distanceMatrixMPEGTest.yml", FileStorage::READ);
+ vector<string> namesHeaders;
+ listShapeNames(namesHeaders);
+
+ // Read generated MAT //
+ fs["distanceMat"]>>distanceMat;
+
+ int corrects=0;
+ int divi=0;
+ for (int row=0; row<distanceMat.rows; row++)
+ {
+ if (row%NSN==0) //another group
+ {
+ divi+=NSN;
+ }
+ for (int col=divi-NSN; col<divi; col++)
+ {
+ int nsmall=0;
+ for (int i=0; i<distanceMat.cols; i++)
+ {
+ if (distanceMat.at<float>(row,col)>distanceMat.at<float>(row,i))
+ {
+ nsmall++;
+ }
+ }
+ if (nsmall<=FIRST_MANY)
+ {
+ corrects++;
+ }
+ }
+ }
+ float porc = 100*float(corrects)/(NSN*distanceMat.rows);
+ std::cout<<"%="<<porc<<std::endl;
+ if (porc >= CURRENT_MAX_ACCUR)
+ ts->set_failed_test_info(cvtest::TS::OK);
+ else
+ ts->set_failed_test_info(cvtest::TS::FAIL_BAD_ACCURACY);
+
+}
+
+
+void CV_HaussTest::run(int /* */)
+{
+ mpegTest();
+ displayMPEGResults();
+ ts->set_failed_test_info(cvtest::TS::OK);
+}
+
+TEST(Hauss, regression) { CV_HaussTest test; test.safe_run(); }
--- /dev/null
+#include "test_precomp.hpp"
+
+CV_TEST_MAIN("cv")
--- /dev/null
+#include "test_precomp.hpp"
--- /dev/null
+#ifdef __GNUC__
+# pragma GCC diagnostic ignored "-Wmissing-declarations"
+# if defined __clang__ || defined __APPLE__
+# pragma GCC diagnostic ignored "-Wmissing-prototypes"
+# pragma GCC diagnostic ignored "-Wextra"
+# endif
+#endif
+
+#ifndef __OPENCV_TEST_PRECOMP_HPP__
+#define __OPENCV_TEST_PRECOMP_HPP__
+
+#include <iostream>
+#include "opencv2/ts.hpp"
+#include "opencv2/core.hpp"
+#include "opencv2/imgproc.hpp"
+#include "opencv2/highgui.hpp"
+#include "opencv2/shape.hpp"
+
+#include "opencv2/opencv_modules.hpp"
+
+#endif
--- /dev/null
+/*M///////////////////////////////////////////////////////////////////////////////////////
+//
+// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
+//
+// By downloading, copying, installing or using the software you agree to this license.
+// If you do not agree to this license, do not download, install,
+// copy or use the software.
+//
+//
+// Intel License Agreement
+// For Open Source Computer Vision Library
+//
+// Copyright (C) 2000, Intel Corporation, all rights reserved.
+// Third party copyrights are property of their respective owners.
+//
+// Redistribution and use in source and binary forms, with or without modification,
+// are permitted provided that the following conditions are met:
+//
+// * Redistribution's of source code must retain the above copyright notice,
+// this list of conditions and the following disclaimer.
+//
+// * Redistribution's in binary form must reproduce the above copyright notice,
+// this list of conditions and the following disclaimer in the documentation
+// and/or other materials provided with the distribution.
+//
+// * The name of Intel Corporation may not be used to endorse or promote products
+// derived from this software without specific prior written permission.
+//
+// This software is provided by the copyright holders and contributors "as is" and
+// any express or implied warranties, including, but not limited to, the implied
+// warranties of merchantability and fitness for a particular purpose are disclaimed.
+// In no event shall the Intel Corporation or contributors be liable for any direct,
+// indirect, incidental, special, exemplary, or consequential damages
+// (including, but not limited to, procurement of substitute goods or services;
+// loss of use, data, or profits; or business interruption) however caused
+// and on any theory of liability, whether in contract, strict liability,
+// or tort (including negligence or otherwise) arising in any way out of
+// the use of this software, even if advised of the possibility of such damage.
+//
+//M*/
+
+#include "test_precomp.hpp"
+
+using namespace cv;
+using namespace std;
+
+const int angularBins=12;
+const int radialBins=4;
+const float minRad=0.2;
+const float maxRad=2;
+const int NSN=5;//10;//20; //number of shapes per class
+const int NP=120; //number of points sympliying the contour
+const float outlierWeight=0.1;
+const int numOutliers=20;
+const float CURRENT_MAX_ACCUR=95.0; //99% and 100% reached in several tests, 95 is fixed as minimum boundary
+
+class CV_ShapeTest : public cvtest::BaseTest
+{
+public:
+ CV_ShapeTest();
+ ~CV_ShapeTest();
+protected:
+ void run(int);
+
+private:
+ void mpegTest();
+ void listShapeNames(vector<string> &listHeaders);
+ vector<Point2f> convertContourType(const Mat &, int n=0 );
+ float computeShapeDistance(vector <Point2f>& queryNormal,
+ vector <Point2f>& queryFlipped1,
+ vector <Point2f>& queryFlipped2,
+ vector<Point2f>& testq);
+ void displayMPEGResults();
+};
+
+CV_ShapeTest::CV_ShapeTest()
+{
+}
+CV_ShapeTest::~CV_ShapeTest()
+{
+}
+
+vector <Point2f> CV_ShapeTest::convertContourType(const Mat& currentQuery, int n)
+{
+ vector<vector<Point> > _contoursQuery;
+ vector <Point2f> contoursQuery;
+ findContours(currentQuery, _contoursQuery, RETR_LIST, CHAIN_APPROX_NONE);
+ for (size_t border=0; border<_contoursQuery.size(); border++)
+ {
+ for (size_t p=0; p<_contoursQuery[border].size(); p++)
+ {
+ contoursQuery.push_back(Point2f((float)_contoursQuery[border][p].x,
+ (float)_contoursQuery[border][p].y));
+ }
+ }
+
+ // In case actual number of points is less than n
+ for (int add=contoursQuery.size()-1; add<n; add++)
+ {
+ contoursQuery.push_back(contoursQuery[contoursQuery.size()-add+1]); //adding dummy values
+ }
+
+ // Uniformly sampling
+ random_shuffle(contoursQuery.begin(), contoursQuery.end());
+ int nStart=n;
+ vector<Point2f> cont;
+ for (int i=0; i<nStart; i++)
+ {
+ cont.push_back(contoursQuery[i]);
+ }
+ return cont;
+}
+
+void CV_ShapeTest::listShapeNames( vector<string> &listHeaders)
+{
+ listHeaders.push_back("apple"); //ok
+ listHeaders.push_back("children"); // ok
+ listHeaders.push_back("device7"); // ok
+ listHeaders.push_back("Heart"); // ok
+ listHeaders.push_back("teddy"); // ok
+}
+
+float CV_ShapeTest::computeShapeDistance(vector <Point2f>& query1, vector <Point2f>& query2,
+ vector <Point2f>& query3, vector <Point2f>& testq)
+{
+ //waitKey(0);
+ Ptr <ShapeContextDistanceExtractor> mysc = createShapeContextDistanceExtractor(angularBins, radialBins, minRad, maxRad);
+ //Ptr <HistogramCostExtractor> cost = createNormHistogramCostExtractor(cv::DIST_L1);
+ Ptr <HistogramCostExtractor> cost = createChiHistogramCostExtractor(30,0.15);
+ //Ptr <HistogramCostExtractor> cost = createEMDHistogramCostExtractor();
+ //Ptr <HistogramCostExtractor> cost = createEMDL1HistogramCostExtractor();
+ mysc->setIterations(1);
+ mysc->setCostExtractor( cost );
+ //mysc->setTransformAlgorithm(createAffineTransformer(true));
+ mysc->setTransformAlgorithm( createThinPlateSplineShapeTransformer() );
+ //mysc->setImageAppearanceWeight(1.6);
+ //mysc->setImageAppearanceWeight(0.0);
+ //mysc->setImages(im1,imtest);
+ return ( std::min( mysc->computeDistance(query1, testq),
+ std::min(mysc->computeDistance(query2, testq), mysc->computeDistance(query3, testq) )));
+}
+
+void CV_ShapeTest::mpegTest()
+{
+ string baseTestFolder="shape/mpeg_test/";
+ string path = cvtest::TS::ptr()->get_data_path() + baseTestFolder;
+ vector<string> namesHeaders;
+ listShapeNames(namesHeaders);
+
+ // distance matrix //
+ Mat distanceMat=Mat::zeros(NSN*namesHeaders.size(), NSN*namesHeaders.size(), CV_32F);
+
+ // query contours (normal v flipped, h flipped) and testing contour //
+ vector<Point2f> contoursQuery1, contoursQuery2, contoursQuery3, contoursTesting;
+
+ // reading query and computing its properties //
+ int counter=0;
+ const int loops=NSN*namesHeaders.size()*NSN*namesHeaders.size();
+ for (size_t n=0; n<namesHeaders.size(); n++)
+ {
+ for (int i=1; i<=NSN; i++)
+ {
+ // read current image //
+ stringstream thepathandname;
+ thepathandname<<path+namesHeaders[n]<<"-"<<i<<".png";
+ Mat currentQuery, flippedHQuery, flippedVQuery;
+ currentQuery=imread(thepathandname.str(), IMREAD_GRAYSCALE);
+ Mat currentQueryBuf=currentQuery.clone();
+ flip(currentQuery, flippedHQuery, 0);
+ flip(currentQuery, flippedVQuery, 1);
+ // compute border of the query and its flipped versions //
+ vector<Point2f> origContour;
+ contoursQuery1=convertContourType(currentQuery, NP);
+ origContour=contoursQuery1;
+ contoursQuery2=convertContourType(flippedHQuery, NP);
+ contoursQuery3=convertContourType(flippedVQuery, NP);
+
+ // compare with all the rest of the images: testing //
+ for (size_t nt=0; nt<namesHeaders.size(); nt++)
+ {
+ for (int it=1; it<=NSN; it++)
+ {
+ // skip self-comparisson //
+ counter++;
+ if (nt==n && it==i)
+ {
+ distanceMat.at<float>(NSN*n+i-1,
+ NSN*nt+it-1)=0;
+ continue;
+ }
+ // read testing image //
+ stringstream thetestpathandname;
+ thetestpathandname<<path+namesHeaders[nt]<<"-"<<it<<".png";
+ Mat currentTest;
+ currentTest=imread(thetestpathandname.str().c_str(), 0);
+ // compute border of the testing //
+ contoursTesting=convertContourType(currentTest, NP);
+
+ // compute shape distance //
+ std::cout<<std::endl<<"Progress: "<<counter<<"/"<<loops<<": "<<100*double(counter)/loops<<"% *******"<<std::endl;
+ std::cout<<"Computing shape distance between "<<namesHeaders[n]<<i<<
+ " and "<<namesHeaders[nt]<<it<<": ";
+ distanceMat.at<float>(NSN*n+i-1, NSN*nt+it-1)=
+ computeShapeDistance(contoursQuery1, contoursQuery2, contoursQuery3, contoursTesting);
+ std::cout<<distanceMat.at<float>(NSN*n+i-1, NSN*nt+it-1)<<std::endl;
+ }
+ }
+ }
+ }
+ // save distance matrix //
+ FileStorage fs(cvtest::TS::ptr()->get_data_path() + baseTestFolder + "distanceMatrixMPEGTest.yml", FileStorage::WRITE);
+ fs << "distanceMat" << distanceMat;
+}
+
+const int FIRST_MANY=2*NSN;
+void CV_ShapeTest::displayMPEGResults()
+{
+ string baseTestFolder="shape/mpeg_test/";
+ Mat distanceMat;
+ FileStorage fs(cvtest::TS::ptr()->get_data_path() + baseTestFolder + "distanceMatrixMPEGTest.yml", FileStorage::READ);
+ vector<string> namesHeaders;
+ listShapeNames(namesHeaders);
+
+ // Read generated MAT //
+ fs["distanceMat"]>>distanceMat;
+
+ int corrects=0;
+ int divi=0;
+ for (int row=0; row<distanceMat.rows; row++)
+ {
+ if (row%NSN==0) //another group
+ {
+ divi+=NSN;
+ }
+ for (int col=divi-NSN; col<divi; col++)
+ {
+ int nsmall=0;
+ for (int i=0; i<distanceMat.cols; i++)
+ {
+ if (distanceMat.at<float>(row,col)>distanceMat.at<float>(row,i))
+ {
+ nsmall++;
+ }
+ }
+ if (nsmall<=FIRST_MANY)
+ {
+ corrects++;
+ }
+ }
+ }
+ float porc = 100*float(corrects)/(NSN*distanceMat.rows);
+ std::cout<<"%="<<porc<<std::endl;
+ if (porc >= CURRENT_MAX_ACCUR)
+ ts->set_failed_test_info(cvtest::TS::OK);
+ else
+ ts->set_failed_test_info(cvtest::TS::FAIL_BAD_ACCURACY);
+ //done
+}
+
+void CV_ShapeTest::run( int /*start_from*/ )
+{
+ mpegTest();
+ displayMPEGResults();
+ ts->set_failed_test_info(cvtest::TS::OK);
+}
+
+TEST(Shape_SCD, regression) { CV_ShapeTest test; test.safe_run(); }
SET(OPENCV_CPP_SAMPLES_REQUIRED_DEPS opencv_core opencv_flann opencv_imgproc
opencv_highgui opencv_ml opencv_video opencv_objdetect opencv_photo opencv_nonfree opencv_softcascade
- opencv_features2d opencv_calib3d opencv_legacy opencv_contrib opencv_stitching opencv_videostab opencv_bioinspired)
+ opencv_features2d opencv_calib3d opencv_legacy opencv_contrib opencv_stitching opencv_videostab opencv_bioinspired opencv_shape)
ocv_check_dependencies(${OPENCV_CPP_SAMPLES_REQUIRED_DEPS})
--- /dev/null
+/*
+ * shape_context.cpp -- Shape context demo for shape matching
+ */
+
+#include "opencv2/shape.hpp"
+#include "opencv2/highgui.hpp"
+#include "opencv2/imgproc.hpp"
+#include <opencv2/core/utility.hpp>
+#include <iostream>
+#include <string>
+
+using namespace std;
+using namespace cv;
+
+static void help()
+{
+ printf("\n"
+ "This program demonstrates a method for shape comparisson based on Shape Context\n"
+ "You should run the program providing a number between 1 and 20 for selecting an image in the folder shape_sample.\n"
+ "Call\n"
+ "./shape_example [number between 1 and 20]\n\n");
+}
+
+static vector<Point> simpleContour( const Mat& currentQuery, int n=300 )
+{
+ vector<vector<Point> > _contoursQuery;
+ vector <Point> contoursQuery;
+ findContours(currentQuery, _contoursQuery, RETR_LIST, CHAIN_APPROX_NONE);
+ for (size_t border=0; border<_contoursQuery.size(); border++)
+ {
+ for (size_t p=0; p<_contoursQuery[border].size(); p++)
+ {
+ contoursQuery.push_back( _contoursQuery[border][p] );
+ }
+ }
+
+ // In case actual number of points is less than n
+ int dummy=0;
+ for (int add=contoursQuery.size()-1; add<n; add++)
+ {
+ contoursQuery.push_back(contoursQuery[dummy++]); //adding dummy values
+ }
+
+ // Uniformly sampling
+ random_shuffle(contoursQuery.begin(), contoursQuery.end());
+ vector<Point> cont;
+ for (int i=0; i<n; i++)
+ {
+ cont.push_back(contoursQuery[i]);
+ }
+ return cont;
+}
+
+int main(int argc, char** argv)
+{
+ help();
+ string path = "./shape_sample/";
+ int indexQuery = 1;
+ if( argc < 2 )
+ {
+ std::cout<<"Using first image as query."<<std::endl;
+ }
+ else
+ {
+ sscanf( argv[1], "%i", &indexQuery );
+ }
+ cv::Ptr <cv::ShapeContextDistanceExtractor> mysc = cv::createShapeContextDistanceExtractor();
+
+ Size sz2Sh(300,300);
+ stringstream queryName;
+ queryName<<path<<indexQuery<<".png";
+ Mat query=imread(queryName.str(), IMREAD_GRAYSCALE);
+ Mat queryToShow;
+ resize(query, queryToShow, sz2Sh);
+ imshow("QUERY", queryToShow);
+ moveWindow("TEST", 0,0);
+ vector<Point> contQuery = simpleContour(query);
+ int bestMatch;
+ float bestDis=FLT_MAX;
+ for ( int ii=1; ii<=20; ii++ )
+ {
+ if (ii==indexQuery) continue;
+ waitKey(30);
+ stringstream iiname;
+ iiname<<path<<ii<<".png";
+ cout<<"name: "<<iiname.str()<<endl;
+ Mat iiIm=imread(iiname.str(), 0);
+ Mat iiToShow;
+ resize(iiIm, iiToShow, sz2Sh);
+ imshow("TEST", iiToShow);
+ moveWindow("TEST", sz2Sh.width+50,0);
+ vector<Point> contii = simpleContour(iiIm);
+ float dis = mysc->computeDistance( contQuery, contii );
+ if ( dis<bestDis )
+ {
+ bestMatch = ii;
+ bestDis = dis;
+ }
+ std::cout<<" distance between "<<queryName.str()<<" and "<<iiname.str()<<" is: "<<dis<<std::endl;
+ }
+ destroyWindow("TEST");
+ stringstream bestname;
+ bestname<<path<<bestMatch<<".png";
+ Mat iiIm=imread(bestname.str(), 0);
+ Mat bestToShow;
+ resize(iiIm, bestToShow, sz2Sh);
+ imshow("BEST MATCH", bestToShow);
+ moveWindow("BEST MATCH", sz2Sh.width+50,0);
+
+ return 0;
+}
--- /dev/null
+/*
+ * shape_context.cpp -- Shape context demo for shape matching
+ */
+
+#include "opencv2/shape.hpp"
+#include "opencv2/highgui.hpp"
+#include "opencv2/imgproc.hpp"
+#include "opencv2/features2d/features2d.hpp"
+#include "opencv2/nonfree/nonfree.hpp"
+#include <opencv2/core/utility.hpp>
+#include <iostream>
+#include <string>
+
+using namespace std;
+using namespace cv;
+
+static void help()
+{
+ printf("\nThis program demonstrates how to use common interface for shape transformers\n"
+ "Call\n"
+ "shape_transformation [image1] [image2]\n");
+}
+
+int main(int argc, char** argv)
+{
+ help();
+ Mat img1 = imread(argv[1], IMREAD_GRAYSCALE);
+ Mat img2 = imread(argv[2], IMREAD_GRAYSCALE);
+ if(img1.empty() || img2.empty() || argc<2)
+ {
+ printf("Can't read one of the images\n");
+ return -1;
+ }
+
+ // detecting keypoints
+ SurfFeatureDetector detector(5000);
+ vector<KeyPoint> keypoints1, keypoints2;
+ detector.detect(img1, keypoints1);
+ detector.detect(img2, keypoints2);
+
+ // computing descriptors
+ SurfDescriptorExtractor extractor;
+ Mat descriptors1, descriptors2;
+ extractor.compute(img1, keypoints1, descriptors1);
+ extractor.compute(img2, keypoints2, descriptors2);
+
+ // matching descriptors
+ BFMatcher matcher(NORM_L2);
+ vector<DMatch> matches;
+ matcher.match(descriptors1, descriptors2, matches);
+
+ // drawing the results
+ namedWindow("matches", 1);
+ Mat img_matches;
+ drawMatches(img1, keypoints1, img2, keypoints2, matches, img_matches);
+ imshow("matches", img_matches);
+
+ // extract points
+ vector<Point2f> pts1, pts2;
+ for (size_t ii=0; ii<keypoints1.size(); ii++)
+ pts1.push_back( keypoints1[ii].pt );
+ for (size_t ii=0; ii<keypoints2.size(); ii++)
+ pts2.push_back( keypoints2[ii].pt );
+
+ // Apply TPS
+ Ptr<ThinPlateSplineShapeTransformer> mytps = createThinPlateSplineShapeTransformer(25000); //TPS with a relaxed constraint
+ mytps->estimateTransformation(pts1, pts2, matches);
+ mytps->warpImage(img2, img2);
+
+ imshow("Tranformed", img2);
+ waitKey(0);
+
+ return 0;
+}
projects / profile / ivi / opencv.git / commitdiff
