sgolay (Signal Processing Toolbox)

Savitzky-Golay filter design

Syntax

b = sgolay(k,f)
b = sgolay(k,f,w)
[b,g] = sgolay(...)

Description

b = sgolay(k,f) designs a Savitzky-Golay FIR smoothing filter b. The polynomial order k must be less than the frame size, f, which must be odd. If k = f-1, the designed filter produces no smoothing. The output, b, is an f-by-f matrix whose rows represent the time-varying FIR filter coefficients. In a smoothing filter implementation (for example, sgolayfilt), the last (f-1)/2 rows (each an FIR filter) are applied to the signal during the startup transient, and the first (f-1)/2 rows are applied to the signal during the terminal transient. The center row is applied to the signal in the steady state.

b = sgolay(k,f,w) specifies a weighting vector w with length f, which contains the real, positive-valued weights to be used during the least-squares minimization.

[b,g] = sgolay(...) [returns the matrix g of differentiation filters. Each column of g is a differentiation filter for derivatives of order p-1 where p is the column index. Given a signal x of length f, you can find an estimate of the p^th order derivative, xp, of its middle value from:

xp((f+1)/2) = (factorial(p)) * g(:,p+1)' * x

Remarks

Savitzky-Golay smoothing filters (also called digital smoothing polynomial filters or least squares smoothing filters) are typically used to "smooth out" a noisy signal whose frequency span (without noise) is large. In this type of application, Savitzky-Golay smoothing filters perform much better than standard averaging FIR filters, which tend to filter out a significant portion of the signal's high frequency content along with the noise. Although Savitzky-Golay filters are more effective at preserving the pertinent high frequency components of the signal, they are less successful than standard averaging FIR filters at rejecting noise.

Savitzky-Golay filters are optimal in the sense that they minimize the least-squares error in fitting a polynomial to each frame of noisy data.

Examples

Use sgolay to smooth a noisy sinusoid and display the result and the first and second derivatives:

N = 4;
F = 21;
[b,g]=sgolay(N,F);
x=5*sin(.4*pi*0:.2:199);
y=x+randn(1,996); % Noisy sinusoid

for n = (F+1)/2:996-(F+1)/2,
% Zero-th order derivative (equivalent to sgolayfilt except
% that it doesn't compute transients)
    z0(n)=g(:,1)'*y(n - (F+1)/2 + 1: n + (F+1)/2 - 1)';
% 1st order derivative
    z1(n)=g(:,2)'*y(n - (F+1)/2 + 1: n + (F+1)/2 - 1)';
% 2nd order derivative
    z2(n)=2*g(:,3)'*y(n - (F+1)/2 + 1: n + (F+1)/2 - 1)';
end

plot([x(1:length(z0))',y(1:length(z0))',z0'])
legend('Noiseless sinusoid','Noisy sinusoid',...
'Smoothed sinusoid')
figure
plot([diff(x(1:length(z0)+1))',z1'])
legend('Noiseless first-order derivative',...
'Smoothed first-order derivative')
figure
plot([diff(diff(x(1:length(z0)+2)))',z2'])
legend('Noiseless second-order derivative',...
'Smoothed second-order derivative')

Note The figures below are zoomed to show more detail.

See Also

fir1, firls, filter, sgolayfilt

References

[1] Orfanidis, S.J., Introduction to Signal Processing, Prentice-Hall, Englewood Cliffs, NJ, 1996.

seqperiod sgolayfilt