Tekkotsu Homepage | Demos | Overview | Downloads | Dev. Resources | Reference | Credits |
sift.ippGo to the documentation of this file.00001 // file: sift.ipp 00002 // author: Andrea Vedaldi 00003 // description: Sift inline members definition 00004 00005 // AUTORIGHTS 00006 // Copyright (c) 2006 The Regents of the University of California 00007 // All Rights Reserved. 00008 // 00009 // Created by Andrea Vedaldi (UCLA VisionLab) 00010 // 00011 // Permission to use, copy, modify, and distribute this software and its 00012 // documentation for educational, research and non-profit purposes, 00013 // without fee, and without a written agreement is hereby granted, 00014 // provided that the above copyright notice, this paragraph and the 00015 // following three paragraphs appear in all copies. 00016 // 00017 // This software program and documentation are copyrighted by The Regents 00018 // of the University of California. The software program and 00019 // documentation are supplied "as is", without any accompanying services 00020 // from The Regents. The Regents does not warrant that the operation of 00021 // the program will be uninterrupted or error-free. The end-user 00022 // understands that the program was developed for research purposes and 00023 // is advised not to rely exclusively on the program for any reason. 00024 // 00025 // This software embodies a method for which the following patent has 00026 // been issued: "Method and apparatus for identifying scale invariant 00027 // features in an image and use of same for locating an object in an 00028 // image," David G. Lowe, US Patent 6,711,293 (March 23, 00029 // 2004). Provisional application filed March 8, 1999. Asignee: The 00030 // University of British Columbia. 00031 // 00032 // IN NO EVENT SHALL THE UNIVERSITY OF CALIFORNIA BE LIABLE TO ANY PARTY 00033 // FOR DIRECT, INDIRECT, SPECIAL, INCIDENTAL, OR CONSEQUENTIAL DAMAGES, 00034 // INCLUDING LOST PROFITS, ARISING OUT OF THE USE OF THIS SOFTWARE AND 00035 // ITS DOCUMENTATION, EVEN IF THE UNIVERSITY OF CALIFORNIA HAS BEEN 00036 // ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. THE UNIVERSITY OF 00037 // CALIFORNIA SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT 00038 // LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 00039 // A PARTICULAR PURPOSE. THE SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS" 00040 // BASIS, AND THE UNIVERSITY OF CALIFORNIA HAS NO OBLIGATIONS TO PROVIDE 00041 // MAINTENANCE, SUPPORT, UPDATES, ENHANCEMENTS, OR MODIFICATIONS. 00042 00043 00044 /** 00045 ** @file 00046 ** @brief SIFT class - inline functions and members 00047 **/ 00048 #include<iostream> 00049 #include<cassert> 00050 00051 namespace VL 00052 { 00053 00054 namespace Detail 00055 { 00056 extern int const expnTableSize ; 00057 extern VL::float_t const expnTableMax ; 00058 extern VL::float_t expnTable [] ; 00059 } 00060 00061 /** @brief Get width of source image 00062 ** @result width. 00063 **/ 00064 inline 00065 int 00066 Sift::getWidth() const 00067 { 00068 return width ; 00069 } 00070 00071 /** @brief Get height of source image 00072 ** @result height. 00073 **/ 00074 inline 00075 int 00076 Sift::getHeight() const 00077 { 00078 return height ; 00079 } 00080 00081 /** @brief Get width of an octave 00082 ** @param o octave index. 00083 ** @result width of octave @a o. 00084 **/ 00085 inline 00086 int 00087 Sift::getOctaveWidth(int o) const 00088 { 00089 assert( omin <= o && o < omin + O ) ; 00090 return (o >= 0) ? (width >> o) : (width << -o) ; 00091 } 00092 00093 /** @brief Get height of an octave 00094 ** @param o octave index. 00095 ** @result height of octave @a o. 00096 **/ 00097 inline 00098 int 00099 Sift::getOctaveHeight(int o) const 00100 { 00101 assert( omin <= o && o < omin + O ) ; 00102 return (o >= 0) ? (height >> o) : (height << -o) ; 00103 } 00104 00105 /** @brief Get octave 00106 ** @param o octave index. 00107 ** @return pointer to octave @a o. 00108 **/ 00109 inline 00110 VL::pixel_t * 00111 Sift::getOctave(int o) 00112 { 00113 assert( omin <= o && o < omin + O ) ; 00114 return octaves[o-omin] ; 00115 } 00116 00117 /** @brief Get level 00118 ** @param o octave index. 00119 ** @param s level index. 00120 ** @result pointer to level @c (o,s). 00121 **/ 00122 inline 00123 VL::pixel_t * 00124 Sift::getLevel(int o, int s) 00125 { 00126 assert( omin <= o && o < omin + O ) ; 00127 assert( smin <= s && s <= smax ) ; 00128 return octaves[o - omin] + 00129 getOctaveWidth(o)*getOctaveHeight(o) * (s-smin) ; 00130 } 00131 00132 /** @brief Get octave sampling period 00133 ** @param o octave index. 00134 ** @result Octave sampling period (in pixels). 00135 **/ 00136 inline 00137 VL::float_t 00138 Sift::getOctaveSamplingPeriod(int o) const 00139 { 00140 return (o >= 0) ? (1 << o) : 1.0f / (1 << -o) ; 00141 } 00142 00143 /** @brief Convert index into scale 00144 ** @param o octave index. 00145 ** @param s scale index. 00146 ** @return scale. 00147 **/ 00148 inline 00149 VL::float_t 00150 Sift::getScaleFromIndex(VL::float_t o, VL::float_t s) const 00151 { 00152 return sigma0 * powf( 2.0f, o + s / S ) ; 00153 } 00154 00155 /** @brief Get keypoint list begin 00156 ** @return iterator to the beginning. 00157 **/ 00158 inline 00159 Sift::KeypointsIter 00160 Sift::keypointsBegin() 00161 { 00162 return keypoints.begin() ; 00163 } 00164 00165 /** @brief Get keypoint list end 00166 ** @return iterator to the end. 00167 **/ 00168 inline 00169 Sift::KeypointsIter 00170 Sift::keypointsEnd() 00171 { 00172 return keypoints.end() ; 00173 } 00174 00175 /** @brief Set normalize descriptor flag */ 00176 inline 00177 void 00178 Sift::setNormalizeDescriptor(bool flag) 00179 { 00180 normalizeDescriptor = flag ; 00181 } 00182 00183 /** @brief Get normalize descriptor flag */ 00184 inline 00185 bool 00186 Sift::getNormalizeDescriptor() const 00187 { 00188 return normalizeDescriptor ; 00189 } 00190 00191 /** @brief Set descriptor magnification */ 00192 inline 00193 void 00194 Sift::setMagnification(VL::float_t _magnif) 00195 { 00196 magnif = _magnif ; 00197 } 00198 00199 /** @brief Get descriptor magnification */ 00200 inline 00201 VL::float_t 00202 Sift::getMagnification() const 00203 { 00204 return magnif ; 00205 } 00206 00207 /** @brief Fast @ exp(-x) 00208 ** 00209 ** The argument must be in the range 0-25.0 (bigger arguments may be 00210 ** truncated to zero). 00211 ** 00212 ** @param x argument. 00213 ** @return @c exp(-x) 00214 **/ 00215 inline 00216 VL::float_t 00217 fast_expn(VL::float_t x) 00218 { 00219 assert(VL::float_t(0) <= x && x <= Detail::expnTableMax) ; 00220 #ifdef VL_USEFASTMATH 00221 x *= Detail::expnTableSize / Detail::expnTableMax ; 00222 VL::int32_t i = fast_floor(x) ; 00223 VL::float_t r = x - i ; 00224 VL::float_t a = VL::Detail::expnTable[i] ; 00225 VL::float_t b = VL::Detail::expnTable[i+1] ; 00226 return a + r * (b - a) ; 00227 #else 00228 return std::exp(-x) ; 00229 #endif 00230 } 00231 00232 /** @brief Fast @c mod(x,2pi) 00233 ** 00234 ** The function quickly computes the value @c mod(x,2pi). 00235 ** 00236 ** @remark The computation is fast only for arguments @a x which are 00237 ** small in modulus. 00238 ** 00239 ** @remark For negative arguments, the semantic of the function is 00240 ** not equivalent to the standard library @c fmod function. 00241 ** 00242 ** @param x function argument. 00243 ** @return @c mod(x,2pi) 00244 **/ 00245 inline 00246 VL::float_t 00247 fast_mod_2pi(VL::float_t x) 00248 { 00249 const VL::float_t TWO_PI = static_cast<VL::float_t>(2 * M_PI); 00250 #ifdef VL_USEFASTMATH 00251 while(x < VL::float_t(0) ) x += VL::float_t(TWO_PI) ; 00252 while(x > VL::float_t(2*M_PI) ) x -= VL::float_t(TWO_PI) ; 00253 return x ; 00254 #else 00255 return (x>=0) ? std::fmod(x, TWO_PI) 00256 : TWO_PI + std::fmod(x, TWO_PI) ; 00257 #endif 00258 } 00259 00260 /** @brief Fast @c (int) floor(x) 00261 ** @param x argument. 00262 ** @return @c float(x) 00263 **/ 00264 inline 00265 int32_t 00266 fast_floor(VL::float_t x) 00267 { 00268 #ifdef VL_USEFASTMATH 00269 return (x>=0)? int32_t(x) : std::floor(x) ; 00270 // return int32_t( x - ((x>=0)?0:1) ) ; 00271 #else 00272 return int32_t( std::floor(x) ) ; 00273 #endif 00274 } 00275 00276 /** @brief Fast @c abs(x) 00277 ** @param x argument. 00278 ** @return @c abs(x) 00279 **/ 00280 inline 00281 VL::float_t 00282 fast_abs(VL::float_t x) 00283 { 00284 #ifdef VL_USEFASTMATH 00285 return (x >= 0) ? x : -x ; 00286 #else 00287 return std::fabs(x) ; 00288 #endif 00289 } 00290 00291 /** @brief Fast @c atan2 00292 ** @param x argument. 00293 ** @param y argument. 00294 ** @return Approximation of @c atan2(x). 00295 **/ 00296 inline 00297 VL::float_t 00298 fast_atan2(VL::float_t y, VL::float_t x) 00299 { 00300 #ifdef VL_USEFASTMATH 00301 00302 /* 00303 The function f(r)=atan((1-r)/(1+r)) for r in [-1,1] is easier to 00304 approximate than atan(z) for z in [0,inf]. To approximate f(r) to 00305 the third degree we may solve the system 00306 00307 f(+1) = c0 + c1 + c2 + c3 = atan(0) = 0 00308 f(-1) = c0 - c1 + c2 - c3 = atan(inf) = pi/2 00309 f(0) = c0 = atan(1) = pi/4 00310 00311 which constrains the polynomial to go through the end points and 00312 the middle point. 00313 00314 We still miss a constrain, which might be simply a constarint on 00315 the derivative in 0. Instead we minimize the Linf error in the 00316 range [0,1] by searching for an optimal value of the free 00317 parameter. This turns out to correspond to the solution 00318 00319 c0=pi/4, c1=-0.9675, c2=0, c3=0.1821 00320 00321 which has maxerr = 0.0061 rad = 0.35 grad. 00322 */ 00323 00324 VL::float_t angle, r ; 00325 VL::float_t const c3 = 0.1821 ; 00326 VL::float_t const c1 = 0.9675 ; 00327 VL::float_t abs_y = fast_abs(y) + VL::float_t(1e-10) ; 00328 00329 if (x >= 0) { 00330 r = (x - abs_y) / (x + abs_y) ; 00331 angle = VL::float_t(M_PI/4.0) ; 00332 } else { 00333 r = (x + abs_y) / (abs_y - x) ; 00334 angle = VL::float_t(3*M_PI/4.0) ; 00335 } 00336 angle += (c3*r*r - c1) * r ; 00337 return (y < 0) ? -angle : angle ; 00338 #else 00339 return std::atan2(y,x) ; 00340 #endif 00341 } 00342 00343 /** @brief Fast @c resqrt 00344 ** @param x argument. 00345 ** @return Approximation to @c resqrt(x). 00346 **/ 00347 inline 00348 float 00349 fast_resqrt(float x) 00350 { 00351 #ifdef VL_USEFASTMATH 00352 // Works if VL::float_t is 32 bit ... 00353 union { 00354 float x ; 00355 VL::int32_t i ; 00356 } u ; 00357 float xhalf = float(0.5) * x ; 00358 u.x = x ; // get bits for floating value 00359 u.i = 0x5f3759df - (u.i>>1); // gives initial guess y0 00360 //u.i = 0xdf59375f - (u.i>>1); // gives initial guess y0 00361 u.x = u.x*(float(1.5) - xhalf*u.x*u.x); // Newton step (may repeat) 00362 u.x = u.x*(float(1.5) - xhalf*u.x*u.x); // Newton step (may repeat) 00363 return u.x ; 00364 #else 00365 return float(1.0) / std::sqrt(x) ; 00366 #endif 00367 } 00368 00369 /** @brief Fast @c resqrt 00370 ** @param x argument. 00371 ** @return Approximation to @c resqrt(x). 00372 **/ 00373 inline 00374 double 00375 fast_resqrt(double x) 00376 { 00377 #ifdef VL_USEFASTMATH 00378 // Works if double is 64 bit ... 00379 union { 00380 double x ; 00381 VL::int64_t i ; 00382 } u ; 00383 double xhalf = double(0.5) * x ; 00384 u.x = x ; // get bits for floating value 00385 u.i = 0x5fe6ec85e7de30daLL - (u.i>>1); // gives initial guess y0 00386 u.x = u.x*(double(1.5) - xhalf*u.x*u.x); // Newton step (may repeat) 00387 u.x = u.x*(double(1.5) - xhalf*u.x*u.x); // Newton step (may repeat) 00388 return u.x ; 00389 #else 00390 return double(1.0) / std::sqrt(x) ; 00391 #endif 00392 } 00393 00394 /** @brief Fast @c sqrt 00395 ** @param x argument. 00396 ** @return Approximation to @c sqrt(x). 00397 **/ 00398 inline 00399 VL::float_t 00400 fast_sqrt(VL::float_t x) 00401 { 00402 #ifdef VL_USEFASTMATH 00403 return (x < 1e-8) ? 0 : x * fast_resqrt(x) ; 00404 #else 00405 return std::sqrt(x) ; 00406 #endif 00407 } 00408 00409 } 00410 00411 // Emacs: 00412 // Local Variables: 00413 // mode: C++ 00414 // End: |
Tekkotsu v5.1CVS |
Generated Mon May 9 04:58:51 2016 by Doxygen 1.6.3 |