5

我正在努力理解在此处找到的 LBP 算法的 Matlab 实现。我试图找出它如何计算每个像素的二进制文件?它只是计算相邻像素大于实际中心像素大小的位置。我想计算每个像素的二进制文件,以便使用局部直方图来计算图像的特征。

[ysize, xsize] = size(image);

miny=min(spoints(:,1));
maxy=max(spoints(:,1));
minx=min(spoints(:,2));
maxx=max(spoints(:,2));

% Block size, each LBP code is computed within a block of size bsizey*bsizex
bsizey=ceil(max(maxy,0))-floor(min(miny,0))+1;
bsizex=ceil(max(maxx,0))-floor(min(minx,0))+1;

% Coordinates of origin (0,0) in the block
origy=1-floor(min(miny,0));
origx=1-floor(min(minx,0));

% Minimum allowed size for the input image depends
% on the radius of the used LBP operator.
if(xsize < bsizex || ysize < bsizey)
   error('Too small input image. Should be at least (2*radius+1) x (2*radius+1)');
end

% Calculate dx and dy;
dx = xsize - bsizex;
dy = ysize - bsizey;

% Fill the center pixel matrix C.
C = image(origy:origy+dy,origx:origx+dx);
d_C = double(C);

bins = 2^neighbors;

% Initialize the result matrix with zeros.
result=zeros(dy+1,dx+1);

%Compute the LBP code image
% the whole process here
for i = 1:neighbors
  y = spoints(i,1)+origy;
  x = spoints(i,2)+origx;
  % Calculate floors, ceils and rounds for the x and y.
  fy = floor(y); cy = ceil(y); ry = round(y);
  fx = floor(x); cx = ceil(x); rx = round(x);
  % Check if interpolation is needed.
  if (abs(x - rx) < 1e-6) && (abs(y - ry) < 1e-6)
    % Interpolation is not needed, use original datatypes
    N = image(ry:ry+dy,rx:rx+dx);
    D = N >= C; 
  else
  % Interpolation needed, use double type images 
  ty = y - fy;
  tx = x - fx;

  % Calculate the interpolation weights.
  w1 = roundn((1 - tx) * (1 - ty),-6);
  w2 = roundn(tx * (1 - ty),-6);
  w3 = roundn((1 - tx) * ty,-6) ;
  % w4 = roundn(tx * ty,-6) ;
  w4 = roundn(1 - w1 - w2 - w3, -6);

  % Compute interpolated pixel values
  N = w1*d_image(fy:fy+dy,fx:fx+dx) + w2*d_image(fy:fy+dy,cx:cx+dx) + ...
 w3*d_image(cy:cy+dy,fx:fx+dx) + w4*d_image(cy:cy+dy,cx:cx+dx);
  N = roundn(N,-4);
  D = N >= d_C; 
 end  
   % Update the result matrix.
  v = 2^(i-1);
  result = result + v*D;
end

 %Apply mapping if it is defined
 if isstruct(mapping)
 bins = mapping.num;
 for i = 1:size(result,1)
    for j = 1:size(result,2)
        result(i,j) = mapping.table(result(i,j)+1);
    end
  end
 end

 if (strcmp(mode,'h') || strcmp(mode,'hist') || strcmp(mode,'nh'))
  % Return with LBP histogram if mode equals 'hist'.
  result=hist(result(:),0:(bins-1));
  if (strcmp(mode,'nh'))
    result=result/sum(result);
  end
 else
 %Otherwise return a matrix of unsigned integers
 if ((bins-1)<=intmax('uint8'))
     result=uint8(result);
 elseif ((bins-1)<=intmax('uint16'))
     result=uint16(result);
 else
     result=uint32(result);
 end
end
 size(result)
end

迭代地,它为每个像素的所有 8 个邻居在结果中添加一些值。但它如何与 LBP 二进制文件相关联?它与以下 c++ LBP 方法的以下代码有何关联:

 uchar lbp(const Mat_<uchar> & img, int x, int y)
 {
  // this is pretty much the same what you already got..
  uchar v = 0;
  uchar c = img(y,x);
  v += (img(y-1,x  ) > c) << 0;
  v += (img(y-1,x+1) > c) << 1;
  v += (img(y  ,x+1) > c) << 2;
  v += (img(y+1,x+1) > c) << 3;
  v += (img(y+1,x  ) > c) << 4;
  v += (img(y+1,x-1) > c) << 5;
  v += (img(y  ,x-1) > c) << 6;
  v += (img(y-1,x-1) > c) << 7;
  return v;

}

4

2 回答 2

8

它是 LBP 的矢量化实现,非常适合 Matlab。

在初始化说明之后,让我们看看主循环,从“ for i = 1:neighbors”行开始。循环非常清晰:它计算一个邻居与中心像素的比较,循环遍历所有邻居。你已经明白了这一点,所以现在深入循环以了解它是如何累积所有结果的。

循环的核心实际上过于复杂,因为它考虑了实圆而不是近似整数圆。所以指令的主要部分的目的是计算相邻像素的插值强度。在这里,它与您作为参考的 C++ 代码不同,它仅采用整数、1 像素宽的半径圆。请记住,使用 lbp.m 代码,您可以(理论上,稍后我将讨论)沿半径为 R 的圆计算 LBP,该圆具有 N 个采样点,因此 C++ 将对应于半径为 1 且具有 8 个采样的圆点,只要没有插值。但是当邻居不适合图像的像素网格时有一个插值,当(abs(x - rx) < 1e-6) && (abs(y - ry) < 1e-6)为假时)。

如果(abs(x - rx) < 1e-6) && (abs(y - ry) < 1e-6)为真,则没有插值,因此中心像素和当前邻居之间所有比较的计算直接存储到D. 否则,它会在整个图像上计算采样相邻点的强度的双线性插值N = w1*d_image(fy:fy+dy,fx:fx+dx) + w2*d_image(fy:fy+dy,cx:cx+dx) + w3*d_image(cy:cy+dy,fx:fx+dx) + w4*d_image(cy:cy+dy,cx:cx+dx);

最后,转到更新部分:v = 2^(i-1); result = result + v*D;. v等效于移位:对于第 i 个邻居,您将比较值i-1向左移动,或等效地乘以 be 2^(i-1)。然后与 相加result。因此,在循环结束时,计算实际上等同于您的 C++ 代码,只是它是在整个图像上完成的,而不是在单个像素上完成的。C++ 代码可以看作是 matlab 循环的展开版本,具有半径为 1 的相邻圆和 8 个采样点。至此,计算LBP图,接下来的块是对LBP图的附加处理(通过映射表重新映射,并且可选地计算LBP图像的直方图而不是LBP图像本身)。

现在,对整个脚本进行一点讨论。这里有一个隐藏在脚本末尾的缺陷。实际上,通过代码,您被限制为 32 个邻居,仅此而已,因为最后 LBP 图像被转换为int32​​ . 缺陷是变量result被分配为双矩阵而不是整数矩阵,所以我真的希望在更新result时没有近似问题,然后在转换为整数时,导致 LBP 中的位发生变化。通常不应该有,因为至少有 52 个精度位(根据IEEE 754 规范的维基百科)。我认为这里有风险......相反,我不知道用于长固定大小、有效位向量的 matlab 类型。我会用int64而不是int32,但限制将是 64 个采样邻居。

编辑

现在,如果您希望计算一些限制在 3*3 邻域内的局部二进制模式,那么这个 Matlab 函数对您来说太通用了,最好的办法是展开这个邻域的循环,从而真正接近C++ 代码。这是一段代码(我使用按位或代替加法,但它是等价的):

result = uint8(ysize, xsize);
result = (image(1:end-2,2:end-1) > image(2:end-1,2:end-1));                                 % <=> v += (img(y-1,x  ) > c) << 0;
result = result|bitshift((image(1:end-2,3:end) > image(2:end-1,2:end-1)), 1, 'uint8');      % <=> v += (img(y-1,x+1) > c) << 1;
result = result|bitshift((image(2:end-1,3:end) > image(2:end-1,2:end-1)), 2, 'uint8');      % <=> v += (img(y  ,x+1) > c) << 2;
result = result|bitshift((image(3:end,3:end) > image(2:end-1,2:end-1)), 3, 'uint8');        % <=> v += (img(y+1,x+1) > c) << 3;
result = result|bitshift((image(3:end,2:end-1) > image(2:end-1,2:end-1)), 4, 'uint8');      % <=> v += (img(y+1,x  ) > c) << 4;
result = result|bitshift((image(3:end,1:end-2) > image(2:end-1,2:end-1)), 5, 'uint8');      % <=> v += (img(y+1,x-1) > c) << 5;
result = result|bitshift((image(2:end-1,3:end) > image(2:end-1,2:end-1)), 6, 'uint8');      % <=> v += (img(y  ,x-1) > c) << 6;
result = result|bitshift((image(1:end-2,1:end-2) > image(2:end-1,2:end-1)), 7, 'uint8');    % <=> v += (img(y-1,x-1) > c) << 7;

它是使用强大的矢量化将 C 代码精确转换为 Matlab 脚本。有了这个,就可以很容易地改变这个社区的另一个订单或不同的测试。我也提到了这一点,因为在这种情况下,Matlab 脚本中有一个错误,第 53 行有一个错误的标志:neighobrhood 更好,spoints=[-1 -1; -1 0; -1 1; 0 -1; 0 -1; 1 -1; 1 0; 1 1];而不是spoints=[-1 -1; -1 0; -1 1; 0 -1; -0 1; 1 -1; 1 0; 1 1];.

于 2014-12-01T10:08:04.413 回答
8

我完成了关于本地二进制模式的最后一年项目。我看到了您指向的代码,但我决定编写自己的代码。

这是我的代码:

function [ LBP ] = LBP( I2)
m=size(I2,1);
n=size(I2,2);
for i=2:m-1
    for j=2:n-1
        J0=I2(i,j);
        I3(i-1,j-1)=I2(i-1,j-1)>J0;
        I3(i-1,j)=I2(i-1,j)>J0;
        I3(i-1,j+1)=I2(i-1,j+1)>J0; 
        I3(i,j+1)=I2(i,j+1)>J0;
        I3(i+1,j+1)=I2(i+1,j+1)>J0; 
        I3(i+1,j)=I2(i+1,j)>J0; 
        I3(i+1,j-1)=I2(i+1,j-1)>J0; 
        I3(i,j-1)=I2(i,j-1)>J0;
        LBP(i,j)=I3(i-1,j-1)*2^7+I3(i-1,j)*2^6+I3(i-1,j+1)*2^5+I3(i,j+1)*2^4+I3(i+1,j+1)*2^3+I3(i+1,j)*2^2+I3(i+1,j-1)*2^1+I3(i,j-1)*2^0;
    end
end
end

I2 是您传递的图像,LBP 是输出。为本地二进制模式写了这个:http: //quantgreeks.com/local-binary-pattern-in-matlab/。我知道我可以以更有效的形式编写代码。但是这样写,可以清楚地说明本地二进制模式是如何工作的。

于 2014-12-01T14:14:03.337 回答