% ----------------------------------------------------------------------------
% Example output from hpeq.m
%
% Butterworth equalizer - boost & cut
% ----------------------------------------------------------------------------
% Sophocles J. Orfanidis
% ECE Department, Rutgers University
% 94 Brett Road, Piscataway, NJ 08854-8058
%
% Email:   orfanidi@ece.rutgers.edu
% Date:    June 15, 2005
%
% Reference: Sophocles J. Orfanidis, "High-Order Digital Parametric Equalizer
%            Design," J. Audio Eng. Soc., vol.53, pp. 1026-1046, November 2005.
%
% Web Page: http://www.ece.rutgers.edu/~orfanidi/hpeq
%
% tested with MATLAB R11.1 and R14
% ----------------------------------------------------------------------------


clear all;

disp('% ------------- Boost case -------------------------------------------------------------');

N = 5; fs = 40; f0 = 5; Df = 4; w0 = 2*pi*f0/fs; Dw = 2*pi*Df/fs;
G0 = 0; G = 12; GB = 9; type = 0;

[B,A,Bh,Ah] = hpeq(N, G0, G, GB, w0, Dw, type)

f = linspace(0,20,1001); w = 2*pi*f/fs;

H = 20*log10(abs(fresp(B,A,w)));

% H = 20*log10(abs(fresp(Bh,Ah,w,w0)));  % alternative computation of frequency response

[w1,w2] = bandedge(w0,Dw); f1 = fs * w1/2/pi; f2 = fs * w2/2/pi;

figure;
plot(f,H,'r-', [f1,f2],([GB,GB]),'b.');
ylim([-8 14]); ytick(-6:3:12);
xlim([0,20]); xtick(0:2:20);
title('Butterworth Boost {\it N} = 5');
xlabel('{\it f}  (kHz)'); ylabel('dB');
grid;

disp('% ------------- Cut case ---------------------------------------------------------------');

N = 5; fs=40; f0 = 5; Df = 4; w0 = 2*pi*f0/fs; Dw = 2*pi*Df/fs;
G0 = 0; G = -12; GB = -9; type = 0;

[B,A,Bh,Ah] = hpeq(N, G0, G, GB, w0, Dw, type)

f = linspace(0,20,1001); w = 2*pi*f/fs;

H = 20*log10(abs(fresp(B,A,w)));

% H = 20*log10(abs(fresp(Bh,Ah,w,w0)));  % alternative computation of frequency response

[w1,w2] = bandedge(w0,Dw); f1 = fs * w1/2/pi; f2 = fs * w2/2/pi;

figure;
plot(f,H,'r-', [f1,f2],[GB,GB],'b.');
ylim([-14 14]); ytick(-12:3:12);
xlim([0,20]); xtick(0:2:20);
title('Butterworth Cut {\it N} = 5');
xlabel('{\it f}  (kHz)'); ylabel('dB');
grid;
% ------------- Boost case -------------------------------------------------------------

B =

    1.0773   -1.0708    0.4370         0         0
    1.1068   -2.4579    2.7996   -1.8899    0.7051
    1.1486   -2.3405    2.2459   -1.1474    0.3050


A =

    1.0000   -1.0708    0.5143         0         0
    1.0000   -2.3893    2.9162   -1.9584    0.6953
    1.0000   -2.1965    2.3395   -1.2915    0.3600


Bh =

    1.0773   -0.4370         0
    1.1068   -1.2624    0.7051
    1.1486   -1.0127    0.3050


Ah =

    1.0000   -0.5143         0
    1.0000   -1.3790    0.6953
    1.0000   -1.1063    0.3600

% ------------- Cut case ---------------------------------------------------------------

B =

    0.9283   -0.9940    0.4774         0         0
    0.9035   -2.1588    2.6348   -1.7694    0.6282
    0.8706   -1.9123    2.0368   -1.1243    0.3134


A =

    1.0000   -0.9940    0.4057         0         0
    1.0000   -2.2207    2.5294   -1.7075    0.6371
    1.0000   -2.0377    1.9553   -0.9990    0.2655


Bh =

    0.9283   -0.4774         0
    0.9035   -1.2460    0.6282
    0.8706   -0.9632    0.3134


Ah =

    1.0000   -0.4057         0
    1.0000   -1.1406    0.6371
    1.0000   -0.8817    0.2655