一、簡介

基於SVD（奇異值分解）的去噪聲技術屬於子空間演算法的一種。簡單的來說我們希望將帶噪訊號向量空間分解為分別由純淨訊號主導和噪聲訊號主導的兩個子空間，然後通過簡單地去除落在“噪聲空間”中的帶噪訊號向量分量來估計純淨訊號。要將帶噪訊號向量空間分解為“訊號子空間”和“噪聲子空間”，可以採用線性代數中的正交矩陣分解技術，特別是奇異值分解（SVD）和特徵值分解（EVD）。

二、原始碼

%------------------------------中值濾波------------------------------
% 首先通過中值濾波對原資料進行處理
clear all;close all;clc;
 load data.mat
figure(1);
orig=Data(5000:9000)';
L=length(orig);
plot(1:L,orig);             %原始資料
N=2:12;
Y1=zeros(length(N),L);
Var=zeros(length(N),L);
mse=zeros(1,length(N));
for i=1:length(N)    %通過計算均方誤差判斷中值濾波的階數最優值
    Y1(i,:)=medfilt1(orig,N(i));
    Var(i,:)=(Y1(i,:)-orig).^2;
    mse(i)=sum(Var(i,:));
end
clear i;
mse_best=min(mse);
n_best=N(find(mse==mse_best));
x1_best=Y1(find(N==n_best),:);

%----------------------奇異值分解的相關操作-----------------------------
%首先是對映到相空間
% 相空間矩陣行數
%m=floor(L/2);   % 最佳維數基本在P=N／2處的一個鄰域內產生，並且在此鄰域所取的維數值的降噪效果比較理想，而且能滿足工程要求。
                % 因此，重構矩陣的結構可以根據振動訊號長度N來確定，工程應用中可以取P=N/2。
% 相空間間隔： 通過對時間序列計算自相關函式，取其小於某值域最下值（一般取0.1）
% r_k=autocorr(orig,length(orig)-1);
% m_k=find(r_k<0.1);
tao=2;
m=floor((L+tao)/(tao+1));
n=L+tao-m*tao;
%n=L+1-m;
y=zeros(1,n);
X=zeros(m,n);
% for i=1:m
%    y=x1_best(i:n+i-1); 
%    X(i,:)=y;
% end
for i=1:m
   y=x1_best(1+(i-1)*tao:n+(i-1)*tao); 
   X(i,:)=y;
end
clear i;
%矩陣的奇異值分解
[U,S,V] = svd(X);
sv = diag(S);
newSv=sv;
[a,b] = size(S);
%startPosToZero=18;    % 重構階次（其確定較難）
%------------重構階次的確定方法-奇異值能力差分譜法---------------------
% for i=1:length(sv)-1
%     pho(i)=(sv(i)^2-sv(i+1)^2)/(max(sv)^2-min(sv)^2);
% end
% figure(2)
% plot(pho);
startPosToZero=11;
newSv(startPosToZero:end) = 0;
newS = [diag(newSv),zeros(a,b-a)];
newX = U*newS*V';  % 重構相空間矩陣

%降噪後的恢復
 

三、執行結果

四、備註

版本：2014a
完整程式碼或代寫加1564658423