|
TuttleOFX
1
|
00001 /** 00002 * @brief This file provides a set of reimplementations of CImg functions by means of the generic image library (gil). 00003 * Functions for image blurring 00004 */ 00005 #ifndef _BLUR_FILTERS_HPP_ 00006 #define _BLUR_FILTERS_HPP_ 00007 00008 #include <tuttle/plugin/exceptions.hpp> 00009 00010 #include <cmath> 00011 #include "boost/scoped_array.hpp" 00012 #include "boost/gil/gil_all.hpp" 00013 00014 00015 namespace tuttle { 00016 namespace imageUtils { 00017 using namespace tuttle::plugin; 00018 namespace bgil = boost::gil; 00019 00020 /** 00021 * @brief The Canny-Deriche filter is a recursive algorithm allowing 00022 * to compute blurred derivatives of order 0,1 or 2 of an image. 00023 * 00024 * @param[in, out] src Source view 00025 * @param[out] dst Destination view 00026 * @param[in] sigma Standard variation of the gaussian distribution 00027 * @param[in] order Derivative order (0 = kind of blur, 00028 * 1 = kind of blurred edge detection, 00029 * 2 = kind of very blurred edge detection) 00030 * @param[in] axe Direction of the blur 00031 * @param[in] cond Miscallenous 00032 * @return Return the result of the Deriche filter 00033 * @see blur 00034 */ 00035 template <typename S_VIEW, typename D_VIEW> 00036 void 00037 deriche( const S_VIEW& src, D_VIEW& dst, 00038 const float sigma, const int order = 0, 00039 const char axe = 'x', const bool cond = true ) 00040 { 00041 int w = src.width( ); 00042 int h = src.height( ); 00043 00044 // Concept checking 00045 bgil::gil_function_requires < bgil::ColorSpacesCompatibleConcept< 00046 typename bgil::color_space_type<S_VIEW>::type, 00047 typename bgil::color_space_type<D_VIEW>::type > > ( ); 00048 00049 if( sigma >= 0.1f || order ) 00050 { 00051 const float 00052 nsigma = sigma < 0.1f ? 0.1f : sigma, 00053 alpha = 1.695f / nsigma, 00054 ema = expf( -alpha ), 00055 ema2 = expf( -2.0f * alpha ), 00056 b1 = -2.0f * ema, 00057 b2 = ema2; 00058 float a0 = 0.0f, a1 = 0.0f, a2 = 0.0f, a3 = 0.0f, coefp = 0.0f, coefn = 0.0f; 00059 switch( order ) 00060 { 00061 case 0: 00062 { 00063 const float k = ( 1.0f - ema ) * ( 1.0f - ema ) / ( 1.0f + 2.0f * alpha * ema - ema2 ); 00064 a0 = k; 00065 a1 = k * ( alpha - 1.0f ) * ema; 00066 a2 = k * ( alpha + 1.0f ) * ema; 00067 a3 = -k * ema2; 00068 break; 00069 } 00070 case 1: 00071 { 00072 const float k = ( 1.0f - ema ) * ( 1.0f - ema ) / ema; 00073 a0 = k * ema; 00074 a1 = a3 = 0; 00075 a2 = -a0; 00076 break; 00077 } 00078 case 2: 00079 { 00080 const float 00081 ea = expf( -alpha ), 00082 k = -( ema2 - 1 ) / ( 2 * alpha * ema ), 00083 kn = ( -2.0f * ( -1.0f + 3.0f * ea - 3.0f * ea * ea + ea * ea * ea ) / 00084 ( 3.0f * ea + 1.0f + 3.0f * ea * ea + ea * ea * ea ) ); 00085 a0 = kn; 00086 a1 = -kn * ( 1.0f + k * alpha ) * ema; 00087 a2 = kn * ( 1.0f - k * alpha ) * ema; 00088 a3 = -kn*ema2; 00089 break; 00090 } 00091 default: { 00092 std::cerr << "Given filter order (order = " << order << ") must be 0, 1 or 2." << std::endl; 00093 return; 00094 } 00095 } 00096 typedef typename bgil::channel_type<S_VIEW>::type sPix_t; 00097 typedef typename bgil::channel_type<D_VIEW>::type dPix_t; 00098 const int nc = src.num_channels( ) > dst.num_channels( ) ? 00099 dst.num_channels( ) : src.num_channels( ); 00100 coefp = ( a0 + a1 ) / ( 1.0f + b1 + b2 ); 00101 coefn = ( a2 + a3 ) / ( 1.0f + b1 + b2 ); 00102 switch( axe ) 00103 { 00104 case 'x': 00105 { 00106 typedef typename S_VIEW::x_iterator sIterator; 00107 typedef typename D_VIEW::x_iterator dIterator; 00108 const int N = w; 00109 boost::scoped_array<sPix_t> Y_scptr( new sPix_t[N] ); 00110 sPix_t *Y = Y_scptr.get( ); 00111 // Recursive based algorithm 00112 for( int v = 0; v < nc; ++v ) 00113 { 00114 for( int y = 0; y < h; ++y ) 00115 { 00116 sIterator iterSrcX = src.row_begin( y ); 00117 dIterator iterDstX = dst.row_end( y ); 00118 sPix_t *ptrY = Y, yb = 0, yp = 0; 00119 sPix_t xp = ( sPix_t ) 0; 00120 if( cond ) 00121 { 00122 xp = ( *iterSrcX )[v]; 00123 yb = yp = ( sPix_t ) ( coefp * xp ); 00124 } 00125 for( int m = 0; m < N; ++m ) 00126 { 00127 const sPix_t xc = ( *iterSrcX )[v]; 00128 iterSrcX++; 00129 const sPix_t yc = *( ptrY++ ) = ( sPix_t ) ( a0 * xc + a1 * xp - b1 * yp - b2 * yb ); 00130 xp = xc; 00131 yb = yp; 00132 yp = yc; 00133 } 00134 sPix_t xn = ( sPix_t ) 0, xa = ( sPix_t ) 0; 00135 sPix_t yn = 0, ya = 0; 00136 if( cond ) 00137 { 00138 xn = xa = ( *( iterSrcX - 1 ) )[v]; 00139 yn = ya = ( sPix_t ) coefn * xn; 00140 } 00141 for( int n = N - 1; n >= 0; --n ) 00142 { 00143 const sPix_t xc = ( *( --iterSrcX ) )[v]; 00144 const sPix_t yc = ( sPix_t ) ( a2 * xn + a3 * xa - b1 * yn - b2 * ya ); 00145 xa = xn; 00146 xn = xc; 00147 ya = yn; 00148 yn = yc; 00149 ( *--iterDstX )[v] = ( dPix_t ) ( *( --ptrY ) + yc ); 00150 } 00151 } 00152 } 00153 break; 00154 } 00155 case 'y': 00156 { 00157 typedef typename S_VIEW::y_iterator sIterator; 00158 typedef typename D_VIEW::y_iterator dIterator; 00159 const int N = h; 00160 boost::scoped_array<sPix_t> Y_scptr( new sPix_t[N] ); 00161 sPix_t *Y = Y_scptr.get( ); 00162 // Recursive based algorithm 00163 for( int v = 0; v < nc; ++v ) 00164 { 00165 for( int x = 0; x < w; ++x ) 00166 { 00167 sIterator ptrX = src.col_begin( x ); 00168 dIterator dPtrX = dst.col_end( x ); 00169 sPix_t *ptrY = Y, yb = 0, yp = 0; 00170 sPix_t xp = ( sPix_t ) 0; 00171 if( cond ) 00172 { 00173 xp = ( *ptrX )[v]; 00174 yb = yp = ( sPix_t ) ( coefp * xp ); 00175 } 00176 for( int m = 0; m < N; ++m ) 00177 { 00178 const sPix_t xc = ( *ptrX )[v]; 00179 ptrX++; 00180 const sPix_t yc = *( ptrY++ ) = ( sPix_t ) ( a0 * xc + a1 * xp - b1 * yp - b2 * yb ); 00181 xp = xc; 00182 yb = yp; 00183 yp = yc; 00184 } 00185 sPix_t xn = ( sPix_t ) 0, xa = ( sPix_t ) 0; 00186 sPix_t yn = 0, ya = 0; 00187 if( cond ) 00188 { 00189 xn = xa = ( *( ptrX - 1 ) )[v]; 00190 yn = ya = ( sPix_t ) coefn*xn; 00191 } 00192 for( int n = N - 1; n >= 0; --n ) 00193 { 00194 const sPix_t xc = ( *( --ptrX ) )[v]; 00195 const sPix_t yc = ( sPix_t ) ( a2 * xn + a3 * xa - b1 * yn - b2 * ya ); 00196 xa = xn; 00197 xn = xc; 00198 ya = yn; 00199 yn = yc; 00200 ( *--dPtrX )[v] = ( dPix_t ) ( *( --ptrY ) + yc ); 00201 } 00202 } 00203 } 00204 break; 00205 } 00206 default: 00207 { 00208 std::cerr << "Invalid axe value: '" << axe << "'" << std::endl; 00209 return; 00210 } 00211 } 00212 } 00213 else 00214 bgil::copy_pixels( src, dst ); 00215 } 00216 00217 /** 00218 * @brief Use deriche filter amon x and y direction to apply blur. 00219 * 00220 * @param[in, out] src Source view 00221 * @param[out] dst Destination view 00222 * @param[in] sigma Standard variation of the gaussian distribution 00223 * 00224 */ 00225 template <typename S_VIEW, typename D_VIEW> 00226 void dericheFilter( const S_VIEW& src, D_VIEW& dst, 00227 const float sigma, int order = 0, 00228 float threshold = 0.0f ) 00229 { 00230 if( order > 0 ) 00231 { 00232 if( src.num_channels( ) < 3 ) { 00233 std::cerr << "Error: destination must have at least 3 channels" << std::endl; 00234 return; 00235 } 00236 int w = src.width( ); 00237 int h = src.height( ); 00238 int xf, c; 00239 const int mc = 5; 00240 const int mc2 = mc / 2; 00241 float Ix, Iy, Ixm1, Iym1, chaos, mean, q; 00242 // Concept checking 00243 bgil::gil_function_requires < bgil::ColorSpacesCompatibleConcept< 00244 typename bgil::color_space_type<S_VIEW>::type, 00245 typename bgil::color_space_type<D_VIEW>::type > > ( ); 00246 typedef typename S_VIEW::locator Loc; 00247 typedef typename S_VIEW::x_iterator sxIterator; 00248 typedef typename D_VIEW::x_iterator dxIterator; 00249 typedef typename bgil::channel_type<D_VIEW>::type pix_t; 00250 00251 typename terry::image_from_view<S_VIEW>::type tmpX( src.dimensions( ) ); 00252 typename terry::image_from_view<S_VIEW>::type tmpY( src.dimensions( ) ); 00253 00254 S_VIEW vTmpX( bgil::view( tmpX ) ); 00255 S_VIEW vTmpY( bgil::view( tmpY ) ); 00256 00257 deriche( src, vTmpX, sigma, order, 'x' ); 00258 deriche( src, vTmpY, sigma, order, 'y' ); 00259 if( threshold == 0.0f ) 00260 { 00261 for( int y = 0; y < h; ++y ) 00262 { 00263 sxIterator svx_it = vTmpX.row_begin( y ); 00264 sxIterator svy_it = vTmpY.row_begin( y ); 00265 dxIterator dst_it = dst.row_begin( y ); 00266 for( int x = 0; x < w; ++x ) 00267 { 00268 Ix = ( ( *svx_it )[0] + ( *svx_it )[1] + ( *svx_it )[2] ) / 3.0f; 00269 Iy = ( ( *svy_it )[0] + ( *svy_it )[1] + ( *svy_it )[2] ) / 3.0f; 00270 ( *dst_it )[0] = ( pix_t ) ( Ix * Ix ); 00271 ( *dst_it )[1] = ( pix_t ) ( Ix * Iy ); 00272 ( *dst_it )[2] = ( pix_t ) ( Iy * Iy ); 00273 ++svx_it; 00274 ++svy_it; 00275 ++dst_it; 00276 } 00277 } 00278 } 00279 else 00280 { 00281 for( int y = 0; y < h; ++y ) 00282 { 00283 sxIterator src_it = src.row_begin( y ); 00284 sxIterator svx_it = vTmpX.row_begin( y ); 00285 sxIterator svy_it = vTmpY.row_begin( y ); 00286 dxIterator dst_it = dst.row_begin( y ); 00287 chaos = Iym1 = Ixm1 = 0.0f; 00288 for( int x = 0; x < w; ++x ) 00289 { 00290 chaos = mean = 0.0f; 00291 Loc loc = src.xy_at( x, y ); 00292 for( c = 0; c < 3; ++c ) 00293 { 00294 if( x > mc2 && x < w - mc2 ) 00295 for( xf = 0; xf < mc; ++xf ) 00296 { 00297 mean += loc( xf - mc2, 0 )[c]; 00298 chaos += powf( loc( xf - mc2, 0 )[c], 2.0f ); 00299 } 00300 if( y > mc2 && y < h - mc2 ) 00301 for( xf = 0; xf < mc; ++xf ) 00302 { 00303 mean += loc( 0, xf - mc2 )[c]; 00304 chaos += powf( loc( 0, xf - mc2 )[c], 2.0f ); 00305 } 00306 } 00307 mean = powf( mean / ( mc * 2.0f * 3 ), 2.0f ); 00308 chaos = fabs( ( chaos / ( mc * 2.0f * 3 ) ) - mean ); 00309 00310 q = threshold * chaos; 00311 00312 Ix = ( ( *svx_it )[0] + ( *svx_it )[1] + ( *svx_it )[2] ) / 3.0f; 00313 Iy = ( ( *svy_it )[0] + ( *svy_it )[1] + ( *svy_it )[2] ) / 3.0f; 00314 if( q > 0.0f ) 00315 { 00316 ( *dst_it )[0] = ( pix_t ) ( floorf( ( Ix * Ix ) / q ) * q ); 00317 ( *dst_it )[1] = ( pix_t ) ( Ix * Iy ); 00318 ( *dst_it )[2] = ( pix_t ) ( floorf( ( Iy * Iy ) / q ) * q ); 00319 } 00320 else 00321 { 00322 ( *dst_it )[0] = ( pix_t ) ( Ix * Ix ); 00323 ( *dst_it )[1] = ( pix_t ) ( Ix * Iy ); 00324 ( *dst_it )[2] = ( pix_t ) ( Iy * Iy ); 00325 } 00326 Iym1 = Iy; 00327 Ixm1 = Ix; 00328 ++svx_it; 00329 ++svy_it; 00330 ++dst_it; 00331 ++src_it; 00332 } 00333 } 00334 } 00335 } 00336 else if( sigma > 0.2f ) 00337 { 00338 typename terry::image_from_view<S_VIEW>::type itmp( src.dimensions( ) ); 00339 00340 S_VIEW vTmp( view( itmp ) ); 00341 00342 deriche( src, vTmp, sigma, order, 'x' ); 00343 deriche( vTmp, dst, sigma, order, 'y' ); 00344 } 00345 else 00346 { 00347 bgil::copy_and_convert_pixels( src, dst ); 00348 } 00349 } 00350 00351 00352 } 00353 } 00354 00355 #endif 00356
1.7.5.1 
