在MATLAB自适应椭圆形的结构化元素(adaptive elliptical structurin

我试图创建一个图像扩张或侵蚀它的自适应椭圆形结构元素。我编写此代码但不幸的是所有的结构元素是ones(2*M+1)

I = input('Enter the input image: ');
M = input('Enter the maximum allowed semi-major axes length: ');

% determining ellipse parameteres from eigen value decomposition of LST

row = size(I,1);
col = size(I,2);
SE = cell(row,col);
padI = padarray(I,[M M],'replicate','both');
padrow = size(padI,1);
padcol = size(padI,2);

for m = M+1:padrow-M
   for n = M+1:padcol-M

      a = (l2(m-M,n-M)+eps/l1(m-M,n-M)+l2(m-M,n-M)+2*eps)*M;
      b = (l1(m-M,n-M)+eps/l1(m-M,n-M)+l2(m-M,n-M)+2*eps)*M;

      if e1(m-M,n-M,1)==0
         phi = pi/2;
      else
         phi = atan(e1(m-M,n-M,2)/e1(m-M,n-M,1));
      end

      % defining structuring element for each pixel of image

      x0 = m;  
      y0 = n;
      se = zeros(2*M+1);
      row_se = 0;
      for i = x0-M:x0+M
         row_se = row_se+1;
         col_se = 0;
         for j = y0-M:y0+M
            col_se = col_se+1;
            x = j-y0;
            y = x0-i;
            if ((x*cos(phi)+y*sin(phi))^2)/a^2+((x*sin(phi)-y*cos(phi))^2)/b^2 <= 1
               se(row_se,col_se) = 1;
            end
         end
      end

      SE{m-M,n-M} = se;
   end
end

a ， b和phi是半长轴和半短轴长度和phi之间角度a和x轴。

我用2个MATLAB函数来计算所述图像的局部结构张量，然后将其特征向量对于每个像素。这些都是矩阵l1 ， l2 ， e1和e2 。

这是你的代码我不明白的一点：

 a = (l2(mM,nM)+eps/l1(mM,nM)+l2(mM,nM)+2*eps)*M; b = (l1(mM,nM)+eps/l1(mM,nM)+l2(mM,nM)+2*eps)*M;

我简化的表达式b至（只删除索引）：

b = (l1+eps/l1+l2+2*eps)*M;

对于l1和l2在正常范围内，我们得到：

b =(approx)= (l1+0/l1+l2+2*0)*M = (l1+l2)*M;

因此， b可以很容易地大于M ，我不认为这是你的意图。的eps在这种情况下也不能防止被零除，其通常添加的目的eps ：如果l1是零， eps/l1是Inf 。

看着这个表情，似乎对我说，你打算这个代替：

b = (l1+eps)/(l1+l2+2*eps)*M;

在这里，您要添加eps到每个特征值，使他们保证非零（结构张量是对称的，半正定）。然后，你将l1的特征值的总和，乘以M ，导致之间的值0和M对于每个轴的。

所以，这似乎是放错地方的括号的情况。

只是为了记录在案，这是你在你的代码需要：

a = (l2(m-M,n-M)+eps ) / ( l1(m-M,n-M)+l2(m-M,n-M)+2*eps)*M;
b = (l1(m-M,n-M)+eps ) / ( l1(m-M,n-M)+l2(m-M,n-M)+2*eps)*M;
                     ^   ^
                added parentheses

请注意，您可以通过定义，环外的简化代码：

[se_x,se_y] = meshgrid(-M:M,-M:M);

内两个环路，经i和j ，构造se然后可简单地写为：

se = ((se_x.*cos(phi)+se_y.*sin(phi)).^2)./a.^2 + ...
     ((se_x.*sin(phi)-se_y.*cos(phi)).^2)./b.^2 <= 1;

（请注意.*及.^运营商，这些做逐元素乘法和力量。）

进一步的轻微改进来自意识到phi首先从计算e1(m,n,1)和e1(m,n,2)然后在呼叫中使用cos和sin 。如果我们假设特征向量是正确的标准化，然后

cos(phi) == e1(m,n,1)
sin(phi) == e1(m,n,2)

但你总是可以确保它们是标准化的：

cos_phi = e1(m-M,n-M,1);
sin_phi = e1(m-M,n-M,2);
len = hypot(cos_phi,sin_phi);
cos_phi = cos_phi / len;
sin_phi = sin_phi / len;
se = ((se_x.*cos_phi+se_y.*sin_phi).^2)./a.^2 + ...
     ((se_x.*sin_phi-se_y.*cos_phi).^2)./b.^2 <= 1;

考虑三角运算是相当昂贵的，这应该加快你的代码位。