移动窗口最小二乘多项式平滑（sg smooth）

移动窗口最小二乘多项式平滑，即S-G数据平滑源代码。可滤掉光谱数据中的高频噪声，较好地保留原信号的真实性。

function y=sgsmooth（xwidthorder）

%   A smooth method first mentioned by Savitzky and Golay in the paper
%   Abraham Savitzky & Marcel J. E. Golay. Smoothing and Differentiation of Data by 
%   Simplified Least Squares Procedures. Analytical Chemistry.196436（8）:1627-1639.
%   
%   Then corrected by Jeans Teinier Yves Termonia & Jules Deltour in 1972
%   published by Analytical Chemistry as well. And Peter A. Gorry in his
%   article: General Least-Squares Smoothing and Differentiation by the Convolution 
%            （Savitzky-Golay） Method （Anal.Chem.199062（6）: 570-573） showed us the
%   generalized implementation of SG least square smoothing method.
%
%   The algorithm showed in this program mainly from Jeans Ternier etc.for calculating 
%   the weighting parameters also from Peter A. Gorry for calculating the leading and 
%   trailing value.Here we set the derivative order s=0. 
%   
%   y=SavGo（xwidthorder）:
%   x ---------- Data you must input and it must be a vector.
%
%   width ------ The window width for smoothing.
%
%   order ------ The order of polynomials. It must be smaller than the window width you 
%                choose for smooothing.
%
%   DNP. 2007.12.16

y=[];
[mn]=size（x）;

if ~any（[mn]==1）
    error（‘The data you want to be smoothed should be a vector!‘）;
    return
elseif m==1
    x=x‘;              % Make x be a column vector.
end

if width<=order
    error（‘The order of polynomials must be smaller then the window width you choose!‘）;
    return
elseif all（[width order]<=0）
    error（‘The window width and order must be positive!‘）;
    return
end

if nargin < 3
    error（‘You must input data for smoothing and choose the window width and polynomial order.‘）;
    return
elseif nargin > 3
    error（‘Too many inputs!‘）;
    return
end

if mod（width2）==0
    error（‘The window width should be odd.‘）;
    return
end

m=（width-1）/2;
B=[];
h=zeros（width）;

% Coefficient matrix ................
for i=-m:m
%     for t=-m:m
        for o=0:order
%             coef=（（2*o+1）*factorial（2*m）/factorial（2*m-o））/（factorial（2*m+o+1）/factorial（2*m））;
%             p_i=grampoly（iom）;
%             p_t=grampoly（tom）;
%             h（i+m+1t+m+1）=h（i+m+1t+m+1）+coef*p_i*p_t;
            B（i+m+1o+1）=i^o;
%         end
    end
end

% ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
% ......... Smoothing data ..............

% y=zeros（length（x）1）;
for i=1:length（x）
    if i-1-m < 0
        y（i）=h（:i）‘*x（1:2*m+1）;
    elseif i > length（x）-m
        y（i）=B（2*m+1-（length（x）-i）:）*（br\bq‘）*x（length（x）-2*m:end）;
    else
        [bq br]=qr（B0）;         % QR decomposition of B for Least Square Fitting.
        y（i）=B（m+1:）*（br\bq‘）*x（i-m:i+m）;
    end
end
        
        
                                                                                                                                                
  

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----
     文件        3591  2007-12-26 21:11  sgsmooth.m