y(n)=y(n)+W(m,n)*X(n-m+1);％滤波器输出 end end end
E=D-y;％计算误差 a=a+1;
lmsMSE(a)=mean(E.^2);％根据目标函数计算MSE for n=1:L-1
W(i,n+1)=W(i,n)+u*E(n)*X(n);％更新权值 end end
a=linspace(1,L,L);%设置画MSE变化曲线的横坐标 f=reshape(E,h,k)+H;％重构图像矩阵 figure,subplot(2,2,1),imshow(F); title('原图像') D=reshape(D,h,k); subplot(2,2,2),imshow(D); title('退化图像') myMSE=mean2((F-f).^2) subplot(2,2,3),imshow(f); title('复原图像') subplot(2,2,4) plot(a,lmsMSE,'r-');
title('当u＝0.000001时的 MSE 变化曲线') xlabel('迭代次数'),ylabel('MSE'); legend('lmsMSE(a)',0);
附录B
function [f,noise] = wiener2(varargin)
%WIENER2 Perform 2-D adaptive noise-removal filtering. % WIENER2 lowpass filters an intensity image that has been % degraded by constant power additive noise. WIENER2 uses a % pixel-wise adaptive Wiener method based on statistics % estimated from a local neighborhood of each pixel. %
% J = WIENER2(I,[M N],NOISE) filters the image I using % pixel-wise adaptive Wiener filtering, using neighborhoods of % size M-by-N to estimate the local image mean and standard % deviation. If you omit the [M N] argument, M and N default to % 3. The additive noise (Gaussian white noise) power is assumed % to be NOISE. %
% [J,NOISE] = WIENER2(I,[M N]) also estimates the additive % noise power before doing the filtering. WIENER2 returns this % estimate as NOISE. %
% Class Support % -------------
% The input image I can be of class uint8, uint16, or double. % The output image J is of the same class as I. %
% Example % -------
% I = imread('saturn.tif');
% J = imnoise(I,'gaussian',0,0.005); % K = wiener2(J,[5 5]); % imshow(J), figure, imshow(K) %
% See also FILTER2, MEDFILT2.
% The following syntax is grandfathered: %
% J = WIENER2(I,[M N],[MBLOCK NBLOCK],NOISE) or [J,NOISE] = % WIENER2(I,[M N],[MBLOCK NBLOCK]) processes the intensity % image I as above but in blocks of size MBLOCK-by-NBLOCK. Use % J = WIENER2(I,[M N],SIZE(I),NOISE) to process the matrix all % at once.
% Copyright 1993-2002 The MathWorks, Inc. % $Revision: 5.18 $ $Date: 2002/03/15 15:29:28 $
% Uses algorithm developed by Lee (1980).
% Reference: \ % Jae S. Lim, pp.536-540.
[g, nhood, block, noise, msg] = ParseInputs(varargin{:});%调用ParseInput()处理不同的输入参数的情况,且将输入参数以数组的形式表示 if (~isempty(msg)) error(msg); end
classin = class(g); classChanged = 0; if ~isa(g, 'double') classChanged = 1; g = im2double(g); end
% Estimate the local mean of f.
localMean = filter2(ones(nhood), g) / prod(nhood);
%用一个所有的值全相等的小矩阵中的二维FIR滤波器对二维数据g进行滤波
%其实质是利用小矩阵对g进行二维卷积运算，最终结果的维数大小
%等于g的维数，取的是卷积结果的中心部分。类似于求均值
% Estimate of the local variance of f.
localVar = filter2(ones(nhood), g.^2) / prod(nhood) - localMean.^2;
%求取类似于方差的localVar。这样的算法正好体现了图像矩阵中 相关性
% Estimate the noise power if necessary. if (isempty(noise)) noise = mean2(localVar);
%若没有输入noise，则取localVar的均值作为noise end
% Compute result
% f = localMean + (max(0, localVar - noise) ./ ... % max(localVar, noise)) .* (g - localMean); %
% Computation is split up to minimize use of memory % for temp arrays. f = g - localMean;
%去除图像信号中的直流分量，使得图像信号成为零均值的信号 g = localVar - noise; g = max(g, 0);
localVar = max(localVar, noise); f = f ./ localVar; f = f .* g; f = f + localMean;
%对处理后的信号再重新加上直流分量
%元素间的