((FREAK*)this)->buildPattern();
+ // Convert to gray if not already
+ Mat grayImage = image;
+ if( image.channels() > 1 )
+ cvtColor( image, grayImage, COLOR_BGR2GRAY );
+
+ // Use 32-bit integers if we won't overflow in the integral image
+ if ((image.depth() == CV_8U || image.depth() == CV_8S) &&
+ (image.rows * image.cols) < 8388608 ) // 8388608 = 2 ^ (32 - 8(bit depth) - 1(sign bit))
+ {
+ // Create the integral image appropriate for our type & usage
+ if (image.depth() == CV_8U)
+ computeDescriptors<uchar, int>(grayImage, keypoints, descriptors);
+ else if (image.depth() == CV_8S)
+ computeDescriptors<char, int>(grayImage, keypoints, descriptors);
+ else
+ CV_Error( Error::StsUnsupportedFormat, "" );
+ } else {
+ // Create the integral image appropriate for our type & usage
+ if ( image.depth() == CV_8U )
+ computeDescriptors<uchar, double>(grayImage, keypoints, descriptors);
+ else if ( image.depth() == CV_8S )
+ computeDescriptors<char, double>(grayImage, keypoints, descriptors);
+ else if ( image.depth() == CV_16U )
+ computeDescriptors<ushort, double>(grayImage, keypoints, descriptors);
+ else if ( image.depth() == CV_16S )
+ computeDescriptors<short, double>(grayImage, keypoints, descriptors);
+ else
+ CV_Error( Error::StsUnsupportedFormat, "" );
+ }
+}
+
+template <typename srcMatType>
+void FREAK::extractDescriptor(srcMatType *pointsValue, void ** ptr) const
+{
+ std::bitset<FREAK_NB_PAIRS>** ptrScalar = (std::bitset<FREAK_NB_PAIRS>**) ptr;
+
+ // extracting descriptor preserving the order of SSE version
+ int cnt = 0;
+ for( int n = 7; n < FREAK_NB_PAIRS; n += 128)
+ {
+ for( int m = 8; m--; )
+ {
+ int nm = n-m;
+ for(int kk = nm+15*8; kk >= nm; kk-=8, ++cnt)
+ {
+ (*ptrScalar)->set(kk, pointsValue[descriptionPairs[cnt].i] >= pointsValue[descriptionPairs[cnt].j]);
+ }
+ }
+ }
+ --(*ptrScalar);
+}
+
+#if CV_SSE2
+template <>
+void FREAK::extractDescriptor(uchar *pointsValue, void ** ptr) const
+{
+ __m128i** ptrSSE = (__m128i**) ptr;
+
+ // note that comparisons order is modified in each block (but first 128 comparisons remain globally the same-->does not affect the 128,384 bits segmanted matching strategy)
+ int cnt = 0;
+ for( int n = FREAK_NB_PAIRS/128; n-- ; )
+ {
+ __m128i result128 = _mm_setzero_si128();
+ for( int m = 128/16; m--; cnt += 16 )
+ {
+ __m128i operand1 = _mm_set_epi8(pointsValue[descriptionPairs[cnt+0].i],
+ pointsValue[descriptionPairs[cnt+1].i],
+ pointsValue[descriptionPairs[cnt+2].i],
+ pointsValue[descriptionPairs[cnt+3].i],
+ pointsValue[descriptionPairs[cnt+4].i],
+ pointsValue[descriptionPairs[cnt+5].i],
+ pointsValue[descriptionPairs[cnt+6].i],
+ pointsValue[descriptionPairs[cnt+7].i],
+ pointsValue[descriptionPairs[cnt+8].i],
+ pointsValue[descriptionPairs[cnt+9].i],
+ pointsValue[descriptionPairs[cnt+10].i],
+ pointsValue[descriptionPairs[cnt+11].i],
+ pointsValue[descriptionPairs[cnt+12].i],
+ pointsValue[descriptionPairs[cnt+13].i],
+ pointsValue[descriptionPairs[cnt+14].i],
+ pointsValue[descriptionPairs[cnt+15].i]);
+
+ __m128i operand2 = _mm_set_epi8(pointsValue[descriptionPairs[cnt+0].j],
+ pointsValue[descriptionPairs[cnt+1].j],
+ pointsValue[descriptionPairs[cnt+2].j],
+ pointsValue[descriptionPairs[cnt+3].j],
+ pointsValue[descriptionPairs[cnt+4].j],
+ pointsValue[descriptionPairs[cnt+5].j],
+ pointsValue[descriptionPairs[cnt+6].j],
+ pointsValue[descriptionPairs[cnt+7].j],
+ pointsValue[descriptionPairs[cnt+8].j],
+ pointsValue[descriptionPairs[cnt+9].j],
+ pointsValue[descriptionPairs[cnt+10].j],
+ pointsValue[descriptionPairs[cnt+11].j],
+ pointsValue[descriptionPairs[cnt+12].j],
+ pointsValue[descriptionPairs[cnt+13].j],
+ pointsValue[descriptionPairs[cnt+14].j],
+ pointsValue[descriptionPairs[cnt+15].j]);
+
+ __m128i workReg = _mm_min_epu8(operand1, operand2); // emulated "not less than" for 8-bit UNSIGNED integers
+ workReg = _mm_cmpeq_epi8(workReg, operand2); // emulated "not less than" for 8-bit UNSIGNED integers
+
+ workReg = _mm_and_si128(_mm_set1_epi16(short(0x8080 >> m)), workReg); // merge the last 16 bits with the 128bits std::vector until full
+ result128 = _mm_or_si128(result128, workReg);
+ }
+ (**ptrSSE) = result128;
+ ++(*ptrSSE);
+ }
+ (*ptrSSE) -= 8;
+}
+#endif
+
+template <typename srcMatType, typename iiMatType>
+void FREAK::computeDescriptors( const Mat& image, std::vector<KeyPoint>& keypoints, Mat& descriptors ) const {
+
Mat imgIntegral;
- integral(image, imgIntegral);
+ integral(image, imgIntegral, DataType<iiMatType>::type);
std::vector<int> kpScaleIdx(keypoints.size()); // used to save pattern scale index corresponding to each keypoints
const std::vector<int>::iterator ScaleIdxBegin = kpScaleIdx.begin(); // used in std::vector erase function
const std::vector<cv::KeyPoint>::iterator kpBegin = keypoints.begin(); // used in std::vector erase function
const float sizeCst = static_cast<float>(FREAK_NB_SCALES/(FREAK_LOG2* nOctaves));
- uchar pointsValue[FREAK_NB_POINTS];
+ srcMatType pointsValue[FREAK_NB_POINTS];
int thetaIdx = 0;
int direction0;
int direction1;
if( !extAll ) {
// extract the best comparisons only
descriptors = cv::Mat::zeros((int)keypoints.size(), FREAK_NB_PAIRS/8, CV_8U);
-#if CV_SSE2
- __m128i* ptr= (__m128i*) (descriptors.data+(keypoints.size()-1)*descriptors.step[0]);
-#else
- std::bitset<FREAK_NB_PAIRS>* ptr = (std::bitset<FREAK_NB_PAIRS>*) (descriptors.data+(keypoints.size()-1)*descriptors.step[0]);
-#endif
+ void *ptr = descriptors.data+(keypoints.size()-1)*descriptors.step[0];
+
for( size_t k = keypoints.size(); k--; ) {
// estimate orientation (gradient)
if( !orientationNormalized ) {
else {
// get the points intensity value in the un-rotated pattern
for( int i = FREAK_NB_POINTS; i--; ) {
- pointsValue[i] = meanIntensity(image, imgIntegral, keypoints[k].pt.x,keypoints[k].pt.y, kpScaleIdx[k], 0, i);
+ pointsValue[i] = meanIntensity<srcMatType, iiMatType>(image, imgIntegral,
+ keypoints[k].pt.x, keypoints[k].pt.y,
+ kpScaleIdx[k], 0, i);
}
direction0 = 0;
direction1 = 0;
}
// extract descriptor at the computed orientation
for( int i = FREAK_NB_POINTS; i--; ) {
- pointsValue[i] = meanIntensity(image, imgIntegral, keypoints[k].pt.x,keypoints[k].pt.y, kpScaleIdx[k], thetaIdx, i);
+ pointsValue[i] = meanIntensity<srcMatType, iiMatType>(image, imgIntegral,
+ keypoints[k].pt.x, keypoints[k].pt.y,
+ kpScaleIdx[k], thetaIdx, i);
}
-#if CV_SSE2
- // note that comparisons order is modified in each block (but first 128 comparisons remain globally the same-->does not affect the 128,384 bits segmanted matching strategy)
- int cnt = 0;
- for( int n = FREAK_NB_PAIRS/128; n-- ; )
- {
- __m128i result128 = _mm_setzero_si128();
- for( int m = 128/16; m--; cnt += 16 )
- {
- __m128i operand1 = _mm_set_epi8(
- pointsValue[descriptionPairs[cnt+0].i],
- pointsValue[descriptionPairs[cnt+1].i],
- pointsValue[descriptionPairs[cnt+2].i],
- pointsValue[descriptionPairs[cnt+3].i],
- pointsValue[descriptionPairs[cnt+4].i],
- pointsValue[descriptionPairs[cnt+5].i],
- pointsValue[descriptionPairs[cnt+6].i],
- pointsValue[descriptionPairs[cnt+7].i],
- pointsValue[descriptionPairs[cnt+8].i],
- pointsValue[descriptionPairs[cnt+9].i],
- pointsValue[descriptionPairs[cnt+10].i],
- pointsValue[descriptionPairs[cnt+11].i],
- pointsValue[descriptionPairs[cnt+12].i],
- pointsValue[descriptionPairs[cnt+13].i],
- pointsValue[descriptionPairs[cnt+14].i],
- pointsValue[descriptionPairs[cnt+15].i]);
-
- __m128i operand2 = _mm_set_epi8(
- pointsValue[descriptionPairs[cnt+0].j],
- pointsValue[descriptionPairs[cnt+1].j],
- pointsValue[descriptionPairs[cnt+2].j],
- pointsValue[descriptionPairs[cnt+3].j],
- pointsValue[descriptionPairs[cnt+4].j],
- pointsValue[descriptionPairs[cnt+5].j],
- pointsValue[descriptionPairs[cnt+6].j],
- pointsValue[descriptionPairs[cnt+7].j],
- pointsValue[descriptionPairs[cnt+8].j],
- pointsValue[descriptionPairs[cnt+9].j],
- pointsValue[descriptionPairs[cnt+10].j],
- pointsValue[descriptionPairs[cnt+11].j],
- pointsValue[descriptionPairs[cnt+12].j],
- pointsValue[descriptionPairs[cnt+13].j],
- pointsValue[descriptionPairs[cnt+14].j],
- pointsValue[descriptionPairs[cnt+15].j]);
-
- __m128i workReg = _mm_min_epu8(operand1, operand2); // emulated "not less than" for 8-bit UNSIGNED integers
- workReg = _mm_cmpeq_epi8(workReg, operand2); // emulated "not less than" for 8-bit UNSIGNED integers
-
- workReg = _mm_and_si128(_mm_set1_epi16(short(0x8080 >> m)), workReg); // merge the last 16 bits with the 128bits std::vector until full
- result128 = _mm_or_si128(result128, workReg);
- }
- (*ptr) = result128;
- ++ptr;
- }
- ptr -= 8;
-#else
- // extracting descriptor preserving the order of SSE version
- int cnt = 0;
- for( int n = 7; n < FREAK_NB_PAIRS; n += 128)
- {
- for( int m = 8; m--; )
- {
- int nm = n-m;
- for(int kk = nm+15*8; kk >= nm; kk-=8, ++cnt)
- {
- ptr->set(kk, pointsValue[descriptionPairs[cnt].i] >= pointsValue[descriptionPairs[cnt].j]);
- }
- }
- }
- --ptr;
-#endif
+
+ // Extract descriptor
+ extractDescriptor<srcMatType>(pointsValue, &ptr);
}
}
else { // extract all possible comparisons for selection
else {
//get the points intensity value in the un-rotated pattern
for( int i = FREAK_NB_POINTS;i--; )
- pointsValue[i] = meanIntensity(image, imgIntegral, keypoints[k].pt.x,keypoints[k].pt.y, kpScaleIdx[k], 0, i);
+ pointsValue[i] = meanIntensity<srcMatType, iiMatType>(image, imgIntegral,
+ keypoints[k].pt.x,keypoints[k].pt.y,
+ kpScaleIdx[k], 0, i);
direction0 = 0;
direction1 = 0;
}
// get the points intensity value in the rotated pattern
for( int i = FREAK_NB_POINTS; i--; ) {
- pointsValue[i] = meanIntensity(image, imgIntegral, keypoints[k].pt.x,
- keypoints[k].pt.y, kpScaleIdx[k], thetaIdx, i);
+ pointsValue[i] = meanIntensity<srcMatType, iiMatType>(image, imgIntegral,
+ keypoints[k].pt.x, keypoints[k].pt.y,
+ kpScaleIdx[k], thetaIdx, i);
}
int cnt(0);
}
// simply take average on a square patch, not even gaussian approx
-uchar FREAK::meanIntensity( const cv::Mat& image, const cv::Mat& integral,
- const float kp_x,
- const float kp_y,
- const unsigned int scale,
- const unsigned int rot,
- const unsigned int point) const {
+template <typename imgType, typename iiType>
+imgType FREAK::meanIntensity( const cv::Mat& image, const cv::Mat& integral,
+ const float kp_x,
+ const float kp_y,
+ const unsigned int scale,
+ const unsigned int rot,
+ const unsigned int point) const {
// get point position in image
const PatternPoint& FreakPoint = patternLookup[scale*FREAK_NB_ORIENTATION*FREAK_NB_POINTS + rot*FREAK_NB_POINTS + point];
const float xf = FreakPoint.x+kp_x;
const float yf = FreakPoint.y+kp_y;
const int x = int(xf);
const int y = int(yf);
- const int& imagecols = image.cols;
// get the sigma:
const float radius = FreakPoint.sigma;
const int r_y = static_cast<int>((yf-y)*1024);
const int r_x_1 = (1024-r_x);
const int r_y_1 = (1024-r_y);
- uchar* ptr = image.data+x+y*imagecols;
unsigned int ret_val;
// linear interpolation:
- ret_val = (r_x_1*r_y_1*int(*ptr));
- ptr++;
- ret_val += (r_x*r_y_1*int(*ptr));
- ptr += imagecols;
- ret_val += (r_x*r_y*int(*ptr));
- ptr--;
- ret_val += (r_x_1*r_y*int(*ptr));
+ ret_val = r_x_1*r_y_1*int(image.at<imgType>(y , x ))
+ + r_x *r_y_1*int(image.at<imgType>(y , x+1))
+ + r_x_1*r_y *int(image.at<imgType>(y+1, x ))
+ + r_x *r_y *int(image.at<imgType>(y+1, x+1));
//return the rounded mean
ret_val += 2 * 1024 * 1024;
- return static_cast<uchar>(ret_val / (4 * 1024 * 1024));
+ return static_cast<imgType>(ret_val / (4 * 1024 * 1024));
}
// expected case:
const int y_top = int(yf-radius+0.5);
const int x_right = int(xf+radius+1.5);//integral image is 1px wider
const int y_bottom = int(yf+radius+1.5);//integral image is 1px higher
- int ret_val;
+ iiType ret_val;
- ret_val = integral.at<int>(y_bottom,x_right);//bottom right corner
- ret_val -= integral.at<int>(y_bottom,x_left);
- ret_val += integral.at<int>(y_top,x_left);
- ret_val -= integral.at<int>(y_top,x_right);
+ ret_val = integral.at<iiType>(y_bottom,x_right);//bottom right corner
+ ret_val -= integral.at<iiType>(y_bottom,x_left);
+ ret_val += integral.at<iiType>(y_top,x_left);
+ ret_val -= integral.at<iiType>(y_top,x_right);
ret_val = ret_val/( (x_right-x_left)* (y_bottom-y_top) );
//~ std::cout<<integral.step[1]<<std::endl;
- return static_cast<uchar>(ret_val);
+ return static_cast<imgType>(ret_val);
}
// pair selection algorithm from a set of training images and corresponding keypoints
namespace cv
{
-static void
-computeIntegralImages( const Mat& matI, Mat& matS, Mat& matT, Mat& _FT )
+template <typename inMatType, typename outMatType> static void
+computeIntegralImages( const Mat& matI, Mat& matS, Mat& matT, Mat& _FT,
+ int iiType )
{
- CV_Assert( matI.type() == CV_8U );
-
int x, y, rows = matI.rows, cols = matI.cols;
- matS.create(rows + 1, cols + 1, CV_32S);
- matT.create(rows + 1, cols + 1, CV_32S);
- _FT.create(rows + 1, cols + 1, CV_32S);
+ matS.create(rows + 1, cols + 1, iiType );
+ matT.create(rows + 1, cols + 1, iiType );
+ _FT.create(rows + 1, cols + 1, iiType );
+
+ const inMatType* I = matI.ptr<inMatType>();
+
+ outMatType *S = matS.ptr<outMatType>();
+ outMatType *T = matT.ptr<outMatType>();
+ outMatType *FT = _FT.ptr<outMatType>();
- const uchar* I = matI.ptr<uchar>();
- int *S = matS.ptr<int>(), *T = matT.ptr<int>(), *FT = _FT.ptr<int>();
- int istep = (int)matI.step, step = (int)(matS.step/sizeof(S[0]));
+ int istep = (int)(matI.step/matI.elemSize());
+ int step = (int)(matS.step/matS.elemSize());
for( x = 0; x <= cols; x++ )
S[x] = T[x] = FT[x] = 0;
}
}
-struct StarFeature
-{
- int area;
- int* p[8];
-};
-
-static int
-StarDetectorComputeResponses( const Mat& img, Mat& responses, Mat& sizes, int maxSize )
+template <typename iiMatType> static int
+StarDetectorComputeResponses( const Mat& img, Mat& responses, Mat& sizes,
+ int maxSize, int iiType )
{
const int MAX_PATTERN = 17;
static const int sizes0[] = {1, 2, 3, 4, 6, 8, 11, 12, 16, 22, 23, 32, 45, 46, 64, 90, 128, -1};
__m128 sizes1_4[MAX_PATTERN];
union { int i; float f; } absmask;
absmask.i = 0x7fffffff;
- volatile bool useSIMD = cv::checkHardwareSupport(CV_CPU_SSE2);
+ volatile bool useSIMD = cv::checkHardwareSupport(CV_CPU_SSE2) && iiType == CV_32S;
#endif
+
+ struct StarFeature
+ {
+ int area;
+ iiMatType* p[8];
+ };
+
StarFeature f[MAX_PATTERN];
Mat sum, tilted, flatTilted;
int y, rows = img.rows, cols = img.cols;
int border, npatterns=0, maxIdx=0;
- CV_Assert( img.type() == CV_8UC1 );
-
responses.create( img.size(), CV_32F );
sizes.create( img.size(), CV_16S );
npatterns += (pairs[npatterns-1][0] >= 0);
maxIdx = pairs[npatterns-1][0];
- computeIntegralImages( img, sum, tilted, flatTilted );
+ // Create the integral image appropriate for our type & usage
+ if ( img.type() == CV_8U )
+ computeIntegralImages<uchar, iiMatType>( img, sum, tilted, flatTilted, iiType );
+ else if ( img.type() == CV_8S )
+ computeIntegralImages<char, iiMatType>( img, sum, tilted, flatTilted, iiType );
+ else if ( img.type() == CV_16U )
+ computeIntegralImages<ushort, iiMatType>( img, sum, tilted, flatTilted, iiType );
+ else if ( img.type() == CV_16S )
+ computeIntegralImages<short, iiMatType>( img, sum, tilted, flatTilted, iiType );
+ else
+ CV_Error( Error::StsUnsupportedFormat, "" );
+
int step = (int)(sum.step/sum.elemSize());
for(int i = 0; i <= maxIdx; i++ )
int ur_area = (2*ur_size + 1)*(2*ur_size + 1);
int t_area = t_size*t_size + (t_size + 1)*(t_size + 1);
- f[i].p[0] = sum.ptr<int>() + (ur_size + 1)*step + ur_size + 1;
- f[i].p[1] = sum.ptr<int>() - ur_size*step + ur_size + 1;
- f[i].p[2] = sum.ptr<int>() + (ur_size + 1)*step - ur_size;
- f[i].p[3] = sum.ptr<int>() - ur_size*step - ur_size;
+ f[i].p[0] = sum.ptr<iiMatType>() + (ur_size + 1)*step + ur_size + 1;
+ f[i].p[1] = sum.ptr<iiMatType>() - ur_size*step + ur_size + 1;
+ f[i].p[2] = sum.ptr<iiMatType>() + (ur_size + 1)*step - ur_size;
+ f[i].p[3] = sum.ptr<iiMatType>() - ur_size*step - ur_size;
- f[i].p[4] = tilted.ptr<int>() + (t_size + 1)*step + 1;
- f[i].p[5] = flatTilted.ptr<int>() - t_size;
- f[i].p[6] = flatTilted.ptr<int>() + t_size + 1;
- f[i].p[7] = tilted.ptr<int>() - t_size*step + 1;
+ f[i].p[4] = tilted.ptr<iiMatType>() + (t_size + 1)*step + 1;
+ f[i].p[5] = flatTilted.ptr<iiMatType>() - t_size;
+ f[i].p[6] = flatTilted.ptr<iiMatType>() + t_size + 1;
+ f[i].p[7] = tilted.ptr<iiMatType>() - t_size*step + 1;
f[i].area = ur_area + t_area;
sizes1[i] = sizes0[i];
for(int i = 0; i <= maxIdx; i++ )
{
- const int** p = (const int**)&f[i].p[0];
+ const iiMatType** p = (const iiMatType**)&f[i].p[0];
__m128i r0 = _mm_sub_epi32(_mm_loadu_si128((const __m128i*)(p[0]+ofs)),
_mm_loadu_si128((const __m128i*)(p[1]+ofs)));
__m128i r1 = _mm_sub_epi32(_mm_loadu_si128((const __m128i*)(p[3]+ofs)),
for(int i = 0; i <= maxIdx; i++ )
{
- const int** p = (const int**)&f[i].p[0];
- vals[i] = p[0][ofs] - p[1][ofs] - p[2][ofs] + p[3][ofs] +
- p[4][ofs] - p[5][ofs] - p[6][ofs] + p[7][ofs];
+ const iiMatType** p = (const iiMatType**)&f[i].p[0];
+ vals[i] = (int)(p[0][ofs] - p[1][ofs] - p[2][ofs] + p[3][ofs] +
+ p[4][ofs] - p[5][ofs] - p[6][ofs] + p[7][ofs]);
}
for(int i = 0; i < npatterns; i++ )
{
void StarDetector::detectImpl( const Mat& image, std::vector<KeyPoint>& keypoints, const Mat& mask ) const
{
Mat grayImage = image;
- if( image.type() != CV_8U ) cvtColor( image, grayImage, COLOR_BGR2GRAY );
+ if( image.channels() > 1 ) cvtColor( image, grayImage, COLOR_BGR2GRAY );
(*this)(grayImage, keypoints);
KeyPointsFilter::runByPixelsMask( keypoints, mask );
void StarDetector::operator()(const Mat& img, std::vector<KeyPoint>& keypoints) const
{
Mat responses, sizes;
- int border = StarDetectorComputeResponses( img, responses, sizes, maxSize );
+ int border;
+
+ // Use 32-bit integers if we won't overflow in the integral image
+ if ((img.depth() == CV_8U || img.depth() == CV_8S) &&
+ (img.rows * img.cols) < 8388608 ) // 8388608 = 2 ^ (32 - 8(bit depth) - 1(sign bit))
+ border = StarDetectorComputeResponses<int>( img, responses, sizes, maxSize, CV_32S );
+ else
+ border = StarDetectorComputeResponses<double>( img, responses, sizes, maxSize, CV_64F );
+
keypoints.clear();
if( border >= 0 )
StarDetectorSuppressNonmax( responses, sizes, keypoints, border,
projects / profile / ivi / opencv.git / commitdiff