0

我在互联网上找到了这个 Matlab 代码。它实际上是整个代码的一部分，可以在这里找到。有人可以逐行解释发生了什么。我真的很绝望...

% ------------- % This is code to make the edge detecting filter % ----%
function filter=gaussfilt(N)

% calculate alpha so the filter fills all N points
alpha=N;
first=-(1-N/2)*exp(-(1-N/2)^2/alpha);
count=0;
while first<.1*(-(1530/4000*N-N/2)*exp(-(1530/4000*N-N/2)^2/alpha))
    count=count+1;
    alpha=N*500*count;
    first=-(1-N/2)*exp(-(1-N/2)^2/alpha);
end

for n=1:N
     filter(n)=-(n-N/2)*exp(-(n-N/2)^2/alpha);   % d/dt of a gaussian
end
filter=filter/sum(abs(filter));     % normalization

return
4

2 回答 2

7

此函数尝试返回高斯梯度 - 当您将其与一组 o 数据点进行卷积时，它将检测数据中的边缘，同时平滑平均水平没有突然变化的区域中的点。这是 N = 15 时的输出：

在此处输入图像描述

然而，如果你将 N 增加到 20，就会发生一些疯狂的事情，因为代码坦率地说是错误的（也很丑陋）。你最终得到的不是一条漂亮的曲线，而是一条直线——一个糟糕的过滤器。这是因为计算新值的尝试alpha是可怕的。

Robert P 已经提供了代码正在执行的逐步描述。让我向您展示编写此函数的“正确方法”（您将看到我使用了 Robert 提到的一些技术，以及其他一些我会解释）...

function myFilter = gaussfilt2(N, alpha)
% myFilter = gaussfilt(N, alpha)
% returns an N point normalized array of filter coefficients
% corresponding to the gradient of a Gaussian over the interval [-1 1]
% with a standard deviation of 1/alpha 
% in other words, the higher alpha, the sharper the filter
% default value for alpha is 3

if ~exist('alpha', 'var')
  alpha = 3; 
end 

x = linspace( -1, 1, N); % create a vector of N values between -1 and 1 inclusive
sigma = 1.0 / alpha; % convert from alpha to sigma as used in Gaussian formula

% compute first derivative, but leave constants out
% we will normalize later by summing over the coefficients
myFilter = -x .* exp( -(x.^2)/(2*sigma.^2)); % using .* for element-by-element operation

% normalize:
myFilter = myFilter / sum( abs( myFilter ) ); % absolute sum of coefficients is now one

当您使用 and 运行此函数时N=150，alpha = 3曲线如下所示：

在此处输入图像描述

一些需要解释的技术：

当%您在函数声明正下方开始注释块（用F1)。总是个好主意

矢量化Matlab 在做“显式循环”方面很糟糕，而在“隐式循环”方面非常擅长。只要有可能，您就想一次对一大堆值进行“相同的计算”。在这种情况下，我计算一个向量 x using thelinspace( x1, x2, n ) function. This returns an equally spaced array ofn values fromx1 tox2` 包括在内。现在我可以用一条语句计算整个函数

逐元素乘法Matlab 真正用于矩阵操作（即“Matlab”中的“Mat”——它代表矩阵，而不是数学）。如果您有两个向量a和b，并且您想将它们逐个元素相乘（因此结果为[a(1)8b(1) a(2)*b(2) ... a(n)*b(n)]），则使用.*运算符。

灵活性不是“硬连线” alpha 的值，而是将其作为第二个参数允许您重复使用相同的功能并更改滤镜的清晰度。允许用户省略变量并提供默认值意味着当函数只有一个参数时编写的程序仍然可以使用它。当具有名称的变量 ( )存在时，exist('alpha', 'var')调用返回。在函数前面添加否定结果 - 就像在其他一些语言中一样。true'var''alpha'~If Not

于 2013-06-16T17:05:38.353 回答
5

我完全同意弗洛里斯的观点，但我会给你一些关于正在发生的事情的线索。根据您的问题，我假设您以前从未使用过 MATLAB，因此我强烈建议您尝试在线提供的许多教程中的一些。

function filter=gaussfilt(N)

这定义了函数的名称gaussfilt、输入变量N和输出变量filter。您的 m 文件必须另存为gaussfilt.m.

filter是 MATLAB 中的许多内置函数之一，因此不是一个好的变量名。我建议你使用一个名称，例如gauss_filter，或基本上其他任何东西。在您的情况下，这可能并不重要，因为您没有使用该filter函数，但是，使用这样的名称是一个坏习惯。这也适用于size,length等名称max。

alpha=N;
first=-(1-N/2)*exp(-(1-N/2)^2/alpha);
count=0;

这些行只不过是为变量名赋值。exp()是自然指数。count用于跟踪while循环循环的次数。

while first<.1*(-(1530/4000*N-N/2)*exp(-(1530/4000*N-N/2)^2/alpha))
    count=count+1;
    alpha=N*500*count;
    first=-(1-N/2)*exp(-(1-N/2)^2/alpha);
end

只要满足以下条件while，就会计算介于两者之间的所有内容。end条件应该是不言自明的。

first < .1*(-(1530/4000*N-N/2)*exp(-(1530/4000*N-N/2)^2/alpha))

count递增，以便您可以跟踪循环执行了多少次。但是，它从未使用过，因此在这种情况下是不必要的。但是，我建议您保留计数器，并在while循环中包含另一个条件，即它应该停止运行 ifcount > 1e6或其他一些较大的数字。这样，您将避免循环在不满足其他条件的情况下永远运行。

for n=1:N
     filter(n)=-(n-N/2)*exp(-(n-N/2)^2/alpha);   % d/dt of a gaussian
end

for n = 1:N是一个将运行N多次的循环，其中n将是1第一次，2第二次等等。filter(n) = ..将值分配给变量中的n第 th 个位置filter，从而创建一个长度为 N 的向量。这是创建向量的一种不好的坏方法。您应该始终为向量分配内存，以避免向量在循环内增长。“增长”的向量非常非常慢。因此，在开始循环之前，您应该执行以下操作：

filter = zeros(1,N);

这会创建一个零向量。这可能有点太多了，但分配值的更好方法filter是使用arrayfun：

filter = arrayfun(@(n) (-(n-N/2) * exp(-(n-N/2)^2 / alpha)), 1:N);

在这个答案中查看第 6 点以了解原因。

最后一行：

filter=filter/sum(abs(filter));

好吧，您将filter值除以值的总和，从而创建filter一个总和等于的新值1。这是您将从函数中获得的输出。

使用该函数时，必须编写如下内容：

filter_vector = gaussfilt(N) % where N is an integer

最后一点，使用空格！它更容易阅读！

再次，我推荐一些 MATLAB 教程......

于 2013-06-16T15:23:06.950 回答