%% ------------------------------------------------------------------------

%%            Output Info about this m-file

fprintf('\n***********************************************************\n');

fprintf('        <DSP using MATLAB> Exameple 9.10 \n\n');

time_stamp = datestr(now, 31);

[wkd1, wkd2] = weekday(today, 'long');

fprintf('      Now is %20s, and it is %7s  \n\n', time_stamp, wkd2);

%% ------------------------------------------------------------------------

% Given parameters:

n = [0:50]; wxp = 0.5*pi; x = cos(wxp*n);

n1 = n(1:9); x1 = x(9:17);             % for plotting purposes

I = 5; Rp = 0.01; As = 50; wp = wxp/I; ws = (2*pi-wxp)/I;

[delta1, delta2] = db2delta(Rp, As); weights = [delta2/delta1, 1];

[N, Fo, Ao, weights] = firpmord([wp, ws]/pi, [1, 0], [delta1, delta2], 2); N = N + 2;

N

%% -----------------------------------------------------------------

%%                             Plot

%% -----------------------------------------------------------------  

% Input signal

Hf1 = figure('units', 'inches', 'position', [1, 1, 8, 6], ...

	'paperunits', 'inches', 'paperposition', [0, 0, 6, 4], ...

	'NumberTitle', 'off', 'Name', 'Exameple 9.10');

set(gcf,'Color','white'); 

TF = 10;

subplot(2, 2, 1);

Hs1 = stem(n1, x1, 'filled'); set(Hs1, 'markersize', 2, 'color', 'g');

axis([-1, 9, -1.2, 1.2]); grid on;

xlabel('n', 'vertical', 'middle'); ylabel('Amplitude');

title('Input Signal x(n)', 'fontsize', TF);

set(gca, 'xtick', [0:4:16]);

set(gca, 'ytick', [-1, 0, 1]);

% Interpolation with Filter Design: Length M=31

%M = 31; F = [0, wp, ws, pi]/pi; A = [I, I, 0, 0];

h = firpm(N, Fo, I*Ao, weights); y = upfirdn(x, h, I);

delay = (N)/2;                                                   % Delay imparted by the filter

m = delay+1:1:50*I+delay+1; y = y(m); m = 0:40; y = y(81:121);   % for plotting

subplot(2, 2, 3);

Hs2 = stem(m, y, 'filled'); axis([-5, 45, -1.2, 1.2]); grid on;

set(Hs2, 'markersize', 2, 'color', 'm');

xlabel('m', 'vertical', 'middle'); ylabel('Amplitude');

title(' Output Singal y(n): I = 5', 'fontsize', TF);

set(gca, 'xtick', [0:4:16]*I);

set(gca, 'ytick', [-1, 0, 1]);

% Filter Desing Plots

[Hr, w, a, L] = Hr_Type1(h); Hr_min = min(Hr); w_min = find(Hr == Hr_min);

H = abs(freqz(h, 1, w)); Hdb = 20*log10(H/max(H)); min_attn = Hdb(w_min);

subplot(2, 2, 2);

plot(w/pi, Hr, 'm', 'linewidth', 1.0); axis([0, 1, -1, 6]); grid on;

xlabel('Frequency in \pi units', 'vertical', 'middle'); ylabel('Amplitude');

title('Amplitude Response ', 'fontsize', TF);

set(gca, 'xtick', [0, wp/pi, ws/pi, 1]);

set(gca, 'ytick', [0, I]);

subplot(2, 2, 4);

plot(w/pi, Hdb, 'm', 'linewidth', 1.0); axis([0, 1, -70, 10]); grid on;

xlabel('Frequency in \pi units', 'vertical', 'middle'); ylabel('Decibels');

title('Log-magnitude Response ', 'fontsize', TF);

set(gca, 'xtick', [0, wp/pi, ws/pi, 1]);

set(gca, 'ytick', [-70, round(min_attn), 0]);