帮我写一个图像去噪的matlab代码，图像中主要为大量高斯噪声，评价指标需要包括CNR，SNR，ENL以及边缘保护指标

在MATLAB中，我们可以使用中值滤波器或者基于小波变换的方法来进行图像去噪，特别是针对高斯噪声。下面是一个简单的例子，使用中值滤波器去除图像噪声，并计算一些常见的评价指标：

% 加载或创建测试图像
img = imread('your_image.jpg'); % 替换为实际图像文件名

% 对图像进行中值滤波去噪
filtered_img = medfilt2(img, [3 3]); % 使用3x3窗口大小的中值滤波器

% 计算原始图像和滤波后的图像的均值和标准差
mean_img = mean(img(:));
std_dev_img = std(img(:));
mean_filtered_img = mean(filtered_img(:));
std_dev_filtered_img = std(filtered_img(:));

% SNR (Signal-to-Noise Ratio) 计算
SNR = 10 * log10(mean_img^2 / std_dev_img^2);

% CNR (Contrast to Noise Ratio) 计算
if exist('edge_img', 'var')
    % 如果已存在边缘图像（例如通过Canny边检测），则使用边缘信息
    CNR = 10 * log10(mean(edge_img.^2)) / std_dev_filtered_img;
else
    % 否则，假设图像本身作为对比度参考
    CNR = 10 * log10(mean_img^2) / std_dev_filtered_img;
end

% ENL (Edge Number to Leakage) 是一种专门评估边缘保留的量，这里简化的处理方式是比较过滤前后的边缘数量差异
edge_detection = edge(img, 'canny'); % Canny 边缘检测
ENL = sum(edge_detection(:)) - sum(filtered_img == edge_detection(:));

% 边缘保护指标可以简单地通过对比边缘检测前后像素的变化率来估算
edge_change_rate = sum(abs(edge_detection(:) - filtered_img(:))) / numel(filtered_img);
edge_protection_index = 1 - edge_change_rate; % 边缘越少变化，保护越好

% 打印结果
fprintf('SNR: %.2f dB\n', SNR);
fprintf('CNR: %.2f dB\n', CNR);
fprintf('ENL: %d\n', ENL);
fprintf('Edge Protection Index: %.4f\n', edge_protection_index);