Blackman Window Example

Below is the Matlab script for creating Figures 4.4 and 4.5 in §4.4.

% Illustrate zero-phase zero-padding around a Blackman window

% Analysis parameters:
M = 31;         % Window length
N = 64;         % FFT length (zero padding factor = N/M)
Mo2 = (M-1)/2;  % Shorthand
dBCut = -100;   % Clip dB at this level

% Signal parameters (synthetic sinusoid):
wxT = pi/exp(1);% Sinusoid frequency in rad/sample
A = 1;          % Sinusoid amplitude
phix = pi/3;    % Sinusoid phase

% Compute the signal x:
n = [0:N-1];    % time indices for sinusoid and FFT
x = A * cos(wxT*n+phix); % a sinusoid

% Compute a causal Blackman window:
% w = blackman(M); % Matlab signal processing toolbox, or
w = .42-.5*cos(2*pi*(0:M-1)/(M-1))+.08*cos(4*pi*(0:M-1)/(M-1));

% Apply the window:
xw = w .* x(1:M);
% Note that we have "skipped" the first Mo2 samples of x
% since the center of the window is its proper time origin.

% Form the zero-padded FFT input buffer.
% Note how the negative-time portion goes on the right:
xwzp = [xw(Mo2+1:M),zeros(1,N-M),xw(1:Mo2)];
wzp = [w(Mo2+1:M),zeros(1,N-M),w(1:Mo2)];

figure(1);
subplot(1,1,1);  % force a clear
subplot(2,1,1); 
n = -Mo2:Mo2;    % time axis for plot
stem(n,xw,'ok'); % windowed data
hold on;
plot(n,w,':k');  % window
plot(n,zeros(1,M),'-k');  % zero line
title('Blackman Windowed Sinusoid');
xlabel('Time (samples)'); 
ylabel('Amplitude');
text(-8,1,'(a)'); 
hold off;

subplot(2,1,2); 
n = 0:N-1;     % FFT buffer time axis 
stem(n,xwzp,'ok'); 
hold on;
plot(n,zeros(1,N),'-k');  % zero line
plot(n,wzp,':k');         % window
plot([N/2,N/2],[-1,1],'--k');  % dividing line
text(N/4,0.2,'positive time'); 
text(N/2+N/8,0.2,'negative time'); 
xlabel('Time (samples)'); 
ylabel('Amplitude');
text(-8,1,'(b)'); 

saveplot('../eps/zpblackmanT.eps');
hold off;

figure(2);
subplot(1,1,1); % force a clear

% Now show the window transform:
Xwzp = fft(xwzp);             % Blackman window transform
spec = 20*log10(abs(Xwzp));   % Spectral magnitude in dB
spec = spec - max(spec);      % Normalize to 0 db max
spec = max(spec,dBCut*ones(1,N)); % clip to dBCut dB

subplot(2,1,1); 
fni = [0:1.0/N:1-1.0/N];   % Normalized frequency axis
bins = [0:N-1];            % Bin axis
Xwzpa = abs(Xwzp)
stem(bins,Xwzpa,'ok'); 
hold on;
plot([N/2,N/2],YLim,'--k'); % dividing line
text(5,7,'positive frequencies'); 
text(45,7,'negative frequencies'); 
%zpblackmanplot3d
axis([0,N-1,YLim]); 
xlabel('Frequency (bins))'); 
ylabel('Magnitude (linear)');
text(-8,8,'(a)'); 

% Replot interpreting upper bin numbers as (-) frequencies:
subplot(2,1,2);
sbins = bins - N/2;
nh = N/2;
sXwzpa = [Xwzpa(nh+1:N),Xwzpa(1:nh)]; % see also fftshift()
%fninf = fni - 0.5;               % Plot (-) freqs on the left
%plot(fninf,sXwzpa,'-k'); axis([-0.5,0.5,dBCut,10]); grid;
%plot(fni,Xwzpa,'-k'); grid;
stem(sbins,abs(sXwzpa),'ok');
axis([-N/2,N/2-1,YLim]); 
hold on;
plot([0,0],YLim,'--k');  % dividing line
text(-20,7,'negative frequencies'); 
text(5,7,'positive frequencies'); 
xlabel('Frequency (bins))'); 
ylabel('Magnitude (linear)');
text(-48,8,'(b)'); 
saveplot('../eps/zpblackmanF.eps');
hold off;

---

Order a Hardcopy of Spectral Audio Signal Processing

Previous: Matlab for Spectrum Analysis Windows
Next: Matlab listing: dpssw.m

---

written by Julius Orion Smith III
Julius Smith's background is in electrical engineering (BS Rice 1975, PhD Stanford 1983). He is presently Professor of Music and Associate Professor (by courtesy) of Electrical Engineering at Stanford's Center for Computer Research in Music and Acoustics (CCRMA), teaching courses and pursuing research related to signal processing applied to music and audio systems. See http://ccrma.stanford.edu/~jos/ for details.