Tekkotsu Homepage
Demos
Overview
Downloads
Dev. Resources
Reference
Credits

sift.ipp

Go 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