canny算子源码（matlab）

function CVHomework1
        I = imread（‘001.bmp‘）;
        I=rgb2gray（I）

        BW = canny_edge（I2.5）;
        figure（1）
        imshow（I）
        figure（2）
        imshow（BW）


% Canny边缘检测的函数
% Input:
%   a: input image
%   sigma: Gaussian的均方差
function e=canny_edge（asigma）

a = double（a）;     % 将图像像素数据转换为double类型
[mn] = size（a）;
e = repmat（logical（uint8（0））mn）;  % 生成初始矩阵

OffGate = 0.0001;  
Perc = 0.7;
Th = 0.4;

pw = 1:30;             % hardcoding. But it‘s large enough if sigma isn‘t too large
sig2 = sigma * sigma;  % 方差

width = max（find（exp（-（pw.*pw）/（2*sig2）） > OffGate））;  % 寻找截断点

t = （-width:width）;
len = 2*width+1;
t3=[t-0.5;t;t+0.5];

dgau = （-t.*exp（-（t.*t）/（2*sig2））/sig2）.‘;              % 一阶高斯函数的导数

ra = size（a1）;  % 图像行数
ca = size（a2）;  % 图像列数
ay = 255*a;
ax = 255*a‘;

ax = conv2（axdgau‘same‘）.‘;   % 高斯平滑滤波后的图像的x方向梯度
ay = conv2（aydgau‘same‘）;     % 高斯平滑滤波后的图像的y方向梯度
mag = sqrt（（ax.*ax）+（ay.*ay））;  % 每个像素点的梯度强度

magmax = max（mag（:））;
if magmax>0
    mag = mag/magmax;  % 归一化
end

[countsx] = imhist（mag64）;    % 直方图统计
high = min（find（cumsum（counts）>Perc*m*n））/64;
low = Th*high;

thresh = [low high];   % 根据直方图统计确定上下限

% Nonmax-suppression
idxStrong = [];
for dir = 1:4
    Localmax = Findlocalmax（diraxaymag）;
    idxWeak = Localmax（mag（Localmax）>low）;
    e（idxWeak） = 1;
    idxStrong = [idxStrong; idxWeak（mag（idxWeak）>high）];
end

rstrong = rem（idxStrong-1m）+1;