%% Cross-Spectrogram of Two Sequences
% Compute and plot the cross-spectrogram of two sequences.
% Specify each sequence to be 4096 samples long.

N = 4096;

%
% To create the first sequence, generate a convex quadratic chirp embedded in white Gaussian noise and bandpass filter it.
% The chirp has an initial normalized frequency of 0.1π that increases to 0.8π by the end of the measurement.
% The 16th-order filter passes normalized frequencies between 0.2π and 0.4π rad/sample and has a stopband attenuation of 60 dB.

rx = chirp(0:N-1,0.1/2,N,0.8/2,'quadratic',[],'convex')'+randn(N,1)/100;
dx = designfilt('bandpassiir','FilterOrder',16, ...
    'StopbandFrequency1',0.2,'StopbandFrequency2',0.4, ...
    'StopbandAttenuation',60);
x = filter(dx,rx);

%
% To create the second sequence, generate a linear chirp embedded in white Gaussian noise and bandstop filter it.
% The chirp has an initial normalized frequency of 0.9π that decreases to 0.1π by the end of the measurement.
% The 16th-order filter stops normalized frequencies between 0.6π and 0.8π rad/sample and has a passband ripple of 1 dB.

ry = chirp(0:N-1,0.9/2,N,0.1/2)'+randn(N,1)/100;
dy = designfilt('bandstopiir','FilterOrder',16, ...
    'PassbandFrequency1',0.6,'PassbandFrequency2',0.8, ...
    'PassbandRipple',1);
y = filter(dy,ry);

%
% Plot the two sequences. Offset the second sequence vertically so both are visible.

plot([x y+2])

%
% Compute and plot the cross-spectrogram of x and y. Use a 512-sample Hamming window. Specify 500 samples of overlap between adjoining segments and 2048 DFT points.

xspectrogram(x,y,hamming(512),500,2048,'yaxis')

%
% Set to zero the cross-spectrogram values smaller than –50 dB.

xspectrogram(x,y,hamming(512),500,2048,'MinThreshold',-50,'yaxis')

% The spectrogram shows the frequency regions that are enhanced or suppressed by the filters.
% Copyright 2012 The MathWorks, Inc.