双边滤波原理理解
2017-08-26 19:04
302 查看
bilateral filter
双边滤波是一种保边缘滤波,输出像素值由邻域像素的加权组合得到,公式如下:当前像素点p的值,由邻域N(p)中的所有像素点q∈N(p)加权求和得到,权重系数为Wpq,
Wpq=Gσs(||p−q||)Gσr(||Ip−Iq||
Gσ为高斯函数,均值为0,标准差为σ.Gσs,Gσr分别为定义域核,值域核,:
Gσs(||p−q||)=exp(−(px−qx)2+(py−qy)22σ2s)
px,py,qx,qy表示像素点p,q在图像中的位置坐标点,Gσs表示当前的滤波输出值与像素点的位置距离有关,离p点越近的点,对p的滤波输出影响越大.
Gσr(||p−q||)=exp(−(Ip−Iq)22σ2r)
Ip,Iq表示像素点p,q的亮度值,表示亮度值越大的像素点对p的滤波输出影响越大.
之所有能够起到保边缘的效果,是因为在权重中引入了值域核Gσr.在边缘处,相似点p,q的亮度值相差较大,因此此时的Gσr在边缘处的跳变也很大,所以输出值在边缘处也有明显的跳变,为了更好说明,引用http://blog.csdn.net/jfuck/article/details/8932978的示意图如下:
由图像可知,原始图像在边缘处,亮度值存在明显的跳变,从而使得值域核Gσr也在边缘处存在明显跳变,从而使得滤波输出能够保留边缘信息.
Bilateral Grid
bilateral grid算法matlab代码如下:function output = bilateralFilter( data, edge, sigmaSpatial, sigmaRange, ... samplingSpatial, samplingRange ) if ~exist( 'edge', 'var' ), edge = data; end inputHeight = size( data, 1 ); inputWidth = size( data, 2 ); if ~exist( 'sigmaSpatial', 'var' ), sigmaSpatial = min( inputWidth, inputHeight ) / 16; end edgeMin = min( edge( : ) ); edgeMax = max( edge( : ) ); edgeDelta = edgeMax - edgeMin; if ~exist( 'sigmaRange', 'var' ), sigmaRange = 0.1 * edgeDelta; end if ~exist( 'samplingSpatial', 'var' ), samplingSpatial = sigmaSpatial; end if ~exist( 'samplingRange', 'var' ), samplingRange = sigmaRange; end if size( data ) ~= size( edge ), error( 'data and edge must be of the same size' ); end % parameters derivedSigmaSpatial = sigmaSpatial / samplingSpatial; derivedSigmaRange = sigmaRange / samplingRange; paddingXY = floor( 2 * derivedSigmaSpatial ) + 1; paddingZ = floor( 2 * derivedSigmaRange ) + 1; % allocate 3D grid downsampledWidth = floor( ( inputWidth - 1 ) / samplingSpatial ) + 1 + 2 * paddingXY; downsampledHeight = floor( ( inputHeight - 1 ) / samplingSpatial ) + 1 + 2 * paddingXY; downsampledDepth = floor( edgeDelta / samplingRange ) + 1 + 2 * paddingZ; gridData = zeros( downsampledHeight, downsampledWidth, downsampledDepth ); gridWeights = zeros( downsampledHeight, downsampledWidth, downsampledDepth ); % compute downsampled indices [ jj, ii ] = meshgrid( 0 : inputWidth - 1, 0 : inputHeight - 1 ); % ii = % 0 0 0 0 0 % 1 1 1 1 1 % 2 2 2 2 2 % jj = % 0 1 2 3 4 % 0 1 2 3 4 % 0 1 2 3 4 % so when iterating over ii( k ), jj( k ) % get: ( 0, 0 ), ( 1, 0 ), ( 2, 0 ), ... (down columns first) di = round( ii / samplingSpatial ) + paddingXY + 1; dj = round( jj / samplingSpatial ) + paddingXY + 1; dz = round( ( edge - edgeMin ) / samplingRange ) + paddingZ + 1; % perform scatter (there's probably a faster way than this) % normally would do downsampledWeights( di, dj, dk ) = 1, but we have to % perform a summation to do box downsampling for k = 1 : numel( dz ), dataZ = data( k ); % traverses the image column wise, same as di( k ) if ~isnan( dataZ ), dik = di( k ); djk = dj( k ); dzk = dz( k ); gridData( dik, djk, dzk ) = gridData( dik, djk, dzk ) + dataZ; gridWeights( dik, djk, dzk ) = gridWeights( dik, djk, dzk ) + 1; end end % make gaussian kernel kernelWidth = 2 * derivedSigmaSpatial + 1; kernelHeight = kernelWidth; kernelDepth = 2 * derivedSigmaRange + 1; halfKernelWidth = floor( kernelWidth / 2 ); halfKernelHeight = floor( kernelHeight / 2 ); halfKernelDepth = floor( kernelDepth / 2 ); [gridX, gridY, gridZ] = meshgrid( 0 : kernelWidth - 1, 0 : kernelHeight - 1, 0 : kernelDepth - 1 ); gridX = gridX - halfKernelWidth; gridY = gridY - halfKernelHeight; gridZ = gridZ - halfKernelDepth; gridRSquared = ( gridX .* gridX + gridY .* gridY ) / ( derivedSigmaSpatial * derivedSigmaSpatial ) + ( gridZ .* gridZ ) / ( derivedSigmaRange * derivedSigmaRange ); kernel = exp( -0.5 * gridRSquared ); % convolve blurredGridData = convn( gridData, kernel, 'same' ); blurredGridWeights = convn( gridWeights, kernel, 'same' ); % divide blurredGridWeights( blurredGridWeights == 0 ) = -2; % avoid divide by 0, won't read there anyway normalizedBlurredGrid = blurredGridData ./ blurredGridWeights; normalizedBlurredGrid( blurredGridWeights < -1 ) = 0; % put 0s where it's undefined blurredGridWeights( blurredGridWeights < -1 ) = 0; % put zeros back % upsample [ jj, ii ] = meshgrid( 0 : inputWidth - 1, 0 : inputHeight - 1 ); % meshgrid does x, then y, so output arguments need to be reversed % no rounding di = ( ii / samplingSpatial ) + paddingXY + 1; dj = ( jj / samplingSpatial ) + paddingXY + 1; dz = ( edge - edgeMin ) / samplingRange + paddingZ + 1; % interpn takes rows, then cols, etc % i.e. size(v,1), then size(v,2), ... output = interpn( normalizedBlurredGrid, di, dj, dz );
相关文章推荐
- 双边滤波和引导滤波的原理
- 图像平滑技术之核算子、均值滤波、中值滤波、高斯滤波、双边滤滤、导向滤波的原理概要及OpenCV代码实现
- 双边滤波算法原理及实现
- 双边滤波原理(Bilateral Filtering)
- 双边滤波的白话理解
- 双边滤波(BilateralFilter)原理
- 2017年10月13日笔记双边滤波的理解及结合另一篇博客的实现
- 双边滤波算法原理
- 双边滤波原理与C++实现
- 双边滤波原理
- 高斯滤波与双边滤波原理
- 图像处理-双边滤波原理
- 双边滤波算法原理
- 深入理解Spark 2.1 Core (六):Standalone模式运行的原理与源码分析
- 理解 OpenStack Swift (2):架构、原理及功能 [Architecture, Implementation and Features]
- 深入理解PHP原理之Session Gc的一个小概率Notice【转自(风雪之隅)】
- 通过Xen理解Oracle VM原理
- 深入理解Spark 2.1 Core (九):迭代计算和Shuffle的原理与源码分析
- 【OpenCV】邻域滤波:方框、高斯、中值、双边滤波
- Unix/Linux fork() 函数一次调用2次返回实现原理(个人理解仅供参考)