--- /dev/null
+// This is based on Rainer Lienhart contribution. Below is the original copyright:
+//
+/*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.
+//
+//
+// University of Augsburg License Agreement
+// For Open Source MultiMedia Computing (MMC) Library
+//
+// Copyright (C) 2007, University of Augsburg, Germany, 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 University of Augsburg, Germany 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 University of Augsburg, Germany 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*/
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
+// Author: Rainer Lienhart
+// email: Rainer.Lienhart@informatik.uni-augsburg.de
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
+
+// Please cite the following two papers:
+// 1. Shechtman, E., Irani, M.:
+// Matching local self-similarities across images and videos.
+// CVPR, (2007)
+// 2. Eva Horster, Thomas Greif, Rainer Lienhart, Malcolm Slaney.
+// Comparing Local Feature Descriptors in pLSA-Based Image Models.
+// 30th Annual Symposium of the German Association for
+// Pattern Recognition (DAGM) 2008, Munich, Germany, June 2008.
+
+#include "_cvaux.h"
+
+namespace cv
+{
+
+SelfSimDescriptor::SelfSimDescriptor()
+{
+ smallSize = DEFAULT_SMALL_SIZE;
+ largeSize = DEFAULT_LARGE_SIZE;
+ numberOfAngles = DEFAULT_NUM_ANGLES;
+ startDistanceBucket = DEFAULT_START_DISTANCE_BUCKET;
+ numberOfDistanceBuckets = DEFAULT_NUM_DISTANCE_BUCKETS;
+}
+
+SelfSimDescriptor::SelfSimDescriptor(int _ssize, int _lsize,
+ int _startDistanceBucket,
+ int _numberOfDistanceBuckets, int _numberOfAngles)
+{
+ smallSize = _ssize;
+ largeSize = _lsize;
+ startDistanceBucket = _startDistanceBucket;
+ numberOfDistanceBuckets = _numberOfDistanceBuckets;
+ numberOfAngles = _numberOfAngles;
+}
+
+SelfSimDescriptor::SelfSimDescriptor(const SelfSimDescriptor& ss)
+{
+ smallSize = ss.smallSize;
+ largeSize = ss.largeSize;
+ startDistanceBucket = ss.startDistanceBucket;
+ numberOfDistanceBuckets = ss.numberOfDistanceBuckets;
+ numberOfAngles = ss.numberOfAngles;
+}
+
+SelfSimDescriptor::~SelfSimDescriptor()
+{
+}
+
+SelfSimDescriptor& SelfSimDescriptor::operator = (const SelfSimDescriptor& ss)
+{
+ if( this != &ss )
+ {
+ smallSize = ss.smallSize;
+ largeSize = ss.largeSize;
+ startDistanceBucket = ss.startDistanceBucket;
+ numberOfDistanceBuckets = ss.numberOfDistanceBuckets;
+ numberOfAngles = ss.numberOfAngles;
+ }
+ return *this;
+}
+
+size_t SelfSimDescriptor::getDescriptorSize() const
+{
+ return numberOfAngles*(numberOfDistanceBuckets - startDistanceBucket);
+}
+
+Size SelfSimDescriptor::getGridSize( Size imgSize, Size winStride ) const
+{
+ winStride.width = std::max(winStride.width, 1);
+ winStride.height = std::max(winStride.height, 1);
+ int border = largeSize/2 + smallSize/2;
+ return Size(std::max(imgSize.width - border*2 + winStride.width - 1, 0)/winStride.width,
+ std::max(imgSize.height - border*2 + winStride.height - 1, 0)/winStride.height);
+}
+
+// TODO: optimized with SSE2
+void SelfSimDescriptor::SSD(const Mat& img, Point pt, Mat& ssd) const
+{
+ int x, y, dx, dy, r0 = largeSize/2, r1 = smallSize/2;
+ int step = img.step;
+ for( y = -r0; y <= r0; y++ )
+ {
+ float* sptr = ssd.ptr<float>(y+r0) + r0;
+ for( x = -r0; x <= r0; x++ )
+ {
+ int sum = 0;
+ const uchar* src0 = img.ptr<uchar>(y + pt.y - r1) + x + pt.x;
+ const uchar* src1 = img.ptr<uchar>(pt.y - r1) + pt.x;
+ for( dy = -r1; dy <= r1; dy++, src0 += step, src1 += step )
+ for( dx = -r1; dx <= r1; dx++ )
+ {
+ int t = src0[dx] - src1[dx];
+ sum += t*t;
+ }
+ sptr[x] = (float)sum;
+ }
+ }
+}
+
+
+void SelfSimDescriptor::compute(const Mat& img, vector<float>& descriptors, Size winStride,
+ const vector<Point>& locations) const
+{
+ CV_Assert( img.depth() == CV_8U );
+
+ winStride.width = std::max(winStride.width, 1);
+ winStride.height = std::max(winStride.height, 1);
+ Size gridSize = getGridSize(img.size(), winStride);
+ int i, nwindows = locations.empty() ? gridSize.width*gridSize.height : (int)locations.size();
+ int border = largeSize/2 + smallSize/2;
+ int fsize = (int)getDescriptorSize();
+ vector<float> tempFeature(fsize+1);
+ descriptors.resize(fsize*nwindows + 1);
+ Mat ssd(largeSize, largeSize, CV_32F), mappingMask;
+ computeLogPolarMapping(mappingMask);
+
+#if 0 //def _OPENMP
+ int nthreads = cvGetNumThreads();
+ #pragma omp parallel for num_threads(nthreads)
+#endif
+ for( i = 0; i < nwindows; i++ )
+ {
+ Point pt;
+ float* feature0 = &descriptors[fsize*i];
+ float* feature = &tempFeature[0];
+ int x, y, j;
+
+ if( !locations.empty() )
+ {
+ pt = locations[i];
+ if( pt.x < border || pt.x >= img.cols - border ||
+ pt.y < border || pt.y >= img.rows - border )
+ {
+ for( j = 0; j < fsize; j++ )
+ feature0[j] = 0.f;
+ continue;
+ }
+ }
+ else
+ pt = Point((i % gridSize.width)*winStride.width + border,
+ (i / gridSize.width)*winStride.height + border);
+
+ SSD(img, pt, ssd);
+
+ // Determine in the local neighborhood the largest difference and use for normalization
+ float var_noise = 1000.f;
+ for( y = -1; y <= 1 ; y++ )
+ for( x = -1 ; x <= 1 ; x++ )
+ var_noise = std::max(var_noise, ssd.at<float>(largeSize/2+y, largeSize/2+x));
+
+ for( j = 0; j <= fsize; j++ )
+ feature[j] = FLT_MAX;
+
+ // Derive feature vector before exp(-x) computation
+ // Idea: for all x,a >= 0, a=const. we have:
+ // max [ exp( -x / a) ] = exp ( -min(x) / a )
+ // Thus, determine min(ssd) and store in feature[...]
+ for( y = 0; y < ssd.rows; y++ )
+ {
+ const schar *mappingMaskPtr = mappingMask.ptr<schar>(y);
+ const float *ssdPtr = ssd.ptr<float>(y);
+ for( x = 0 ; x < ssd.cols; x++ )
+ {
+ int index = mappingMaskPtr[x];
+ feature[index] = std::max(feature[index], ssdPtr[x]);
+ }
+ }
+
+ var_noise = -1.f/var_noise;
+ for( j = 0; j < fsize; j++ )
+ feature0[j] = feature[j]*var_noise;
+ Mat _f(1, fsize, CV_32F, feature0);
+ cv::exp(_f, _f);
+ }
+}
+
+void SelfSimDescriptor::computeLogPolarMapping(Mat& mappingMask) const
+{
+ mappingMask.create(largeSize, largeSize, CV_8S);
+
+ // What we want is
+ // log_m (radius) = numberOfDistanceBuckets
+ // <==> log_10 (radius) / log_10 (m) = numberOfDistanceBuckets
+ // <==> log_10 (radius) / numberOfDistanceBuckets = log_10 (m)
+ // <==> m = 10 ^ log_10(m) = 10 ^ [log_10 (radius) / numberOfDistanceBuckets]
+ //
+ int radius = largeSize/2, angleBucketSize = 360 / numberOfAngles;
+ int fsize = (int)getDescriptorSize();
+ double inv_log10m = (double)numberOfDistanceBuckets/log10((double)radius);
+
+ for (int y=-radius ; y<=radius ; y++)
+ {
+ schar* mrow = mappingMask.ptr<schar>(y+radius);
+ for (int x=-radius ; x<=radius ; x++)
+ {
+ int index = fsize;
+ float dist = (float)std::sqrt((float)x*x + (float)y*y);
+ int distNo = dist > 0 ? cvRound(log10(dist)*inv_log10m) : 0;
+ if( startDistanceBucket <= distNo && distNo < numberOfDistanceBuckets )
+ {
+ float angle = std::atan2( (float)y, (float)x ) / (float)CV_PI * 180.0f;
+ if (angle < 0) angle += 360.0f;
+ int angleInt = (cvRound(angle) + angleBucketSize/2) % 360;
+ int angleIndex = angleInt / angleBucketSize;
+ index = (distNo-startDistanceBucket)*numberOfAngles + angleIndex;
+ }
+ mrow[x + radius] = saturate_cast<schar>(index);
+ }
+ }
+}
+
+}
projects / opencv / blobdiff