From ac177b849c2860d7469fcd0dda6edf6a18c997b6 Mon Sep 17 00:00:00 2001 From: danielenricocahall Date: Tue, 4 Aug 2020 22:18:32 -0400 Subject: [PATCH] Improve initialization performance of Brisk reformatting Improve initialization performance of Brisk fix formatting Improve initialization performance of Brisk formatting Improve initialization performance of Brisk make a lookup table for ring use cosine/sine lookup table for theta in brisk and utilize trig identity fix ring lookup table use cosine/sine lookup table for theta in brisk and utilize trig identity formatting use cosine/sine lookup table for theta in brisk and utilize trig identity move scale radius product to ring loop to ensure it's not recomputed for each rot revert change move scale radius product to ring loop to ensure it's not recomputed for each rot remove rings lookup table move scale radius product to ring loop to ensure it's not recomputed for each rot fix formatting of for loop move scale radius product to ring loop to ensure it's not recomputed for each rot use sine/cosine approximations for brisk lookup table. add documentation for sine/cosine lookup tables Improve initialization performance of BRISK --- modules/features2d/src/brisk.cpp | 96 ++++++++++++++++++++------------ 1 file changed, 59 insertions(+), 37 deletions(-) diff --git a/modules/features2d/src/brisk.cpp b/modules/features2d/src/brisk.cpp index 4038279d75..98ae1c4f4e 100644 --- a/modules/features2d/src/brisk.cpp +++ b/modules/features2d/src/brisk.cpp @@ -353,13 +353,30 @@ BRISK_Impl::generateKernel(const std::vector &radiusList, const int rings = (int)radiusList.size(); CV_Assert(radiusList.size() != 0 && radiusList.size() == numberList.size()); points_ = 0; // remember the total number of points + double sineThetaLookupTable[n_rot_]; + double cosThetaLookupTable[n_rot_]; for (int ring = 0; ring < rings; ring++) { points_ += numberList[ring]; } + + // using a sine/cosine approximation for the lookup table + // utilizes the trig identities: + // sin(a + b) = sin(a)cos(b) + cos(a)sin(b) + // cos(a + b) = cos(a)cos(b) - sin(a)sin(b) + // and the fact that sin(0) = 0, cos(0) = 1 + double cosval = 1., sinval = 0.; + double dcos = cos(2*CV_PI/double(n_rot_)), dsin = sin(2*CV_PI/double(n_rot_)); + for( size_t rot = 0; rot < n_rot_; ++rot) + { + sineThetaLookupTable[rot] = sinval; + cosThetaLookupTable[rot] = cosval; + double t = sinval*dcos + cosval*dsin; + cosval = cosval*dcos - sinval*dsin; + sinval = t; + } // set up the patterns patternPoints_ = new BriskPatternPoint[points_ * scales_ * n_rot_]; - BriskPatternPoint* patternIterator = patternPoints_; // define the scale discretization: static const float lb_scale = (float)(std::log(scalerange_) / std::log(2.0)); @@ -370,46 +387,51 @@ BRISK_Impl::generateKernel(const std::vector &radiusList, const float sigma_scale = 1.3f; - for (unsigned int scale = 0; scale < scales_; ++scale) - { - scaleList_[scale] = (float)std::pow((double) 2.0, (double) (scale * lb_scale_step)); - sizeList_[scale] = 0; - - // generate the pattern points look-up - double alpha, theta; - for (size_t rot = 0; rot < n_rot_; ++rot) - { - theta = double(rot) * 2 * CV_PI / double(n_rot_); // this is the rotation of the feature - for (int ring = 0; ring < rings; ++ring) - { - for (int num = 0; num < numberList[ring]; ++num) - { - // the actual coordinates on the circle - alpha = (double(num)) * 2 * CV_PI / double(numberList[ring]); - patternIterator->x = (float)(scaleList_[scale] * radiusList[ring] * cos(alpha + theta)); // feature rotation plus angle of the point - patternIterator->y = (float)(scaleList_[scale] * radiusList[ring] * sin(alpha + theta)); - // and the gaussian kernel sigma - if (ring == 0) - { - patternIterator->sigma = sigma_scale * scaleList_[scale] * 0.5f; - } - else - { - patternIterator->sigma = (float)(sigma_scale * scaleList_[scale] * (double(radiusList[ring])) - * sin(CV_PI / numberList[ring])); + for (unsigned int scale = 0; scale < scales_; ++scale) { + scaleList_[scale] = (float) std::pow((double) 2.0, (double) (scale * lb_scale_step)); + sizeList_[scale] = 0; + BriskPatternPoint *patternIteratorOuter = patternPoints_ + (scale * n_rot_ * points_); + // generate the pattern points look-up + for (int ring = 0; ring < rings; ++ring) { + double scaleRadiusProduct = scaleList_[scale] * radiusList[ring]; + float patternSigma = 0.0f; + if (ring == 0) { + patternSigma = sigma_scale * scaleList_[scale] * 0.5f; + } else { + patternSigma = (float) (sigma_scale * scaleList_[scale] * (double(radiusList[ring])) + * sin(CV_PI / numberList[ring])); } // adapt the sizeList if necessary - const unsigned int size = cvCeil(((scaleList_[scale] * radiusList[ring]) + patternIterator->sigma)) + 1; - if (sizeList_[scale] < size) - { - sizeList_[scale] = size; + const unsigned int size = cvCeil(((scaleList_[scale] * radiusList[ring]) + patternSigma)) + 1; + if (sizeList_[scale] < size) { + sizeList_[scale] = size; + } + for (int num = 0; num < numberList[ring]; ++num) { + BriskPatternPoint *patternIterator = patternIteratorOuter; + double alpha = (double(num)) * 2 * CV_PI / double(numberList[ring]); + double sine_alpha = sin(alpha); + double cosine_alpha = cos(alpha); + + for (size_t rot = 0; rot < n_rot_; ++rot) { + double cosine_theta = cosThetaLookupTable[rot]; + double sine_theta = sineThetaLookupTable[rot]; + + // the actual coordinates on the circle + // sin(a + b) = sin(a) cos(b) + cos(a) sin(b) + // cos(a + b) = cos(a) cos(b) - sin(a) sin(b) + patternIterator->x = (float) (scaleRadiusProduct * + (cosine_theta * cosine_alpha - + sine_theta * sine_alpha)); // feature rotation plus angle of the point + patternIterator->y = (float) (scaleRadiusProduct * + (sine_theta * cosine_alpha + cosine_theta * sine_alpha)); + patternIterator->sigma = patternSigma; + // and the gaussian kernel sigma + // increment the iterator + patternIterator += points_; + } + ++patternIteratorOuter; } - - // increment the iterator - ++patternIterator; - } } - } } // now also generate pairings