MATLAB 实用案例三：图像边缘检测、数据拟合与可视化、信号处理

MATLAB 实用案例：图像边缘检测

使用 MATLAB 实现经典的 Canny 边缘检测算法，完整代码如下：

% 读取图像并转换为灰度图
img = imread('lena.png');
if size(img, 3) == 3
    gray_img = rgb2gray(img);
else
    gray_img = img;
end

% 高斯滤波降噪
sigma = 1.5;
filter_size = 2*ceil(3*sigma)+1;
gauss_filter = fspecial('gaussian', filter_size, sigma);
smoothed_img = imfilter(gray_img, gauss_filter, 'replicate');

% Sobel算子计算梯度
[Gx, Gy] = imgradientxy(smoothed_img, 'sobel');
[Gmag, Gdir] = imgradient(Gx, Gy);
Gdir = atan2d(Gy, Gx);

% 非极大值抑制
[nms_img] = non_max_suppression(Gmag, Gdir);

% 双阈值检测和边缘连接
low_threshold = 0.05 * max(nms_img(:));
high_threshold = 0.2 * max(nms_img(:));
edge_img = hysteresis_threshold(nms_img, low_threshold, high_threshold);

% 显示结果
figure;
subplot(1,2,1); imshow(gray_img); title('原始图像');
subplot(1,2,2); imshow(edge_img); title('边缘检测结果');

% 非极大值抑制函数
function [nms_img] = non_max_suppression(mag, dir)
    [rows, cols] = size(mag);
    nms_img = zeros(rows, cols);
    dir = mod(dir + 180, 180);
    
    for i = 2:rows-1
        for j = 2:cols-1
            angle = dir(i,j);
            if (angle >= 0 && angle < 22.5) || (angle >= 157.5 && angle <= 180)
                neighbor1 = mag(i, j+1);
                neighbor2 = mag(i, j-1);
            elseif angle >= 22.5 && angle < 67.5
                neighbor1 = mag(i+1, j-1);
                neighbor2 = mag(i-1, j+1);
            elseif angle >= 67.5 && angle < 112.5
                neighbor1 = mag(i+1, j);
                neighbor2 = mag(i-1, j);
            else
                neighbor1 = mag(i-1, j-1);
                neighbor2 = mag(i+1, j+1);
            end
            
            if mag(i,j) >= neighbor1 && mag(i,j) >= neighbor2
                nms_img(i,j) = mag(i,j);
            end
        end
    end
end

% 滞后阈值函数
function [edge_img] = hysteresis_threshold(img, low, high)
    [rows, cols] = size(img);
    edge_img = zeros(rows, cols);
    
    strong_edges = (img >= high);
    weak_edges = (img >= low) & (img < high);
    
    % 8-连通区域分析
    [L, num] = bwlabel(strong_edges, 8);
    stats = regionprops(L, 'PixelIdxList');
    
    for k = 1:num
        weak_pixels = weak_edges(stats(k).PixelIdxList);
        if any(weak_pixels)
            edge_img(stats(k).PixelIdxList) = 1;
        end
    end
end

MATLAB 实用案例：数据拟合与可视化

使用 MATLAB 进行多项式拟合和三维曲面绘制：

% 生成示例数据
x = linspace(0, 10, 100);
y = 2*x.^3 - 3*x.^2 + 5*x - 1 + 10*randn(size(x));

% 多项式拟合
p = polyfit(x, y, 3);
y_fit = polyval(p, x);

% 计算拟合指标
R2 = 1 - sum((y - y_fit).^2)/sum((y - mean(y)).^2);
RMSE = sqrt(mean((y - y_fit).^2));

% 绘制结果
figure;
plot(x, y, 'o', 'MarkerSize', 6, 'MarkerEdgeColor', 'b');
hold on;
plot(x, y_fit, 'r-', 'LineWidth', 2);
xlabel('X'); ylabel('Y');
title(sprintf('多项式拟合: R²=%.4f, RMSE=%.4f', R2, RMSE));
legend('原始数据', '拟合曲线');

% 三维曲面绘制
[X,Y] = meshgrid(-10:0.5:10);
Z = sin(sqrt(X.^2 + Y.^2))./sqrt(X.^2 + Y.^2);

figure;
surf(X, Y, Z);
shading interp;
colormap(jet);
xlabel('X'); ylabel('Y'); zlabel('Z');
title('三维曲面图');
colorbar;

MATLAB 实用案例：信号处理

实现傅里叶变换和滤波器设计：

% 生成含噪声信号
fs = 1000; % 采样频率
t = 0:1/fs:1-1/fs;
f1 = 50;  % 信号频率1
f2 = 120; % 信号频率2
x = 0.7*sin(2*pi*f1*t) + sin(2*pi*f2*t);
noise = 2*randn(size(t));
x_noisy = x + noise;

% 傅里叶变换
N = length(x_noisy);
f = (0:N-1)*(fs/N);
X = fft(x_noisy);

% 设计带阻滤波器
order = 6;
cutoff = [40 60; 110 130];
[b, a] = butter(order, cutoff/(fs/2), 'stop');
y_filtered = filtfilt(b, a, x_noisy);

% 绘制结果
figure;
subplot(3,1,1);
plot(t, x_noisy);
title('含噪声信号');
xlabel('时间 (s)'); ylabel('幅值');

subplot(3,1,2);
plot(f(1:N/2), abs(X(1:N/2)));
title('频域分析');
xlabel('频率 (Hz)'); ylabel('幅值');

subplot(3,1,3);
plot(t, y_filtered);
title('滤波后信号');
xlabel('时间 (s)'); ylabel('幅值');