%Chapter 10, Example 10.4, Fig 10.17
%Eric Dubois, updated 2018-12-21
% Least pth design for polygonal spec 
% Results are generated for L=6 and L=10.
% Change file name accordingly
% Impose quadrantal symmetry, but no constraints on the response
%Required functions: ideal_response_qs, filter_design_qs_unc
close all, clear all;
%Generate ideal filter response
Npt=128;
Hd = ideal_response_qs('polygonal_specs_example.txt',Npt);
L = 6; %Filter size will be 2L+1 x 2L+1
%window design
h1LP = fwind1(Hd,hamming(2*L+1));
%Compute and plot the frequency response
[H1LP,fx,fy] = freqz2(h1LP,64,64,[1 1]);
H1LPmag = abs(H1LP);
[Fx,Fy] = meshgrid(fx,fy);
%perspective plot
figure;
mesh(Fx(1:2:64,1:2:64),Fy(1:2:64,1:2:64),abs(H1LP(1:2:64,1:2:64)));            % generate the perspective plot
colormap([0 0 0]);     % use black only
axis([-.5 .5 -.5 .5 0 1.2]);
xlabel('\it u \rm (c/px)'), ylabel('\it v \rm (c/px)');
set(gca,'ydir','reverse'); 
set(gca,'fontname','times')          %set the plot font to Times
set(gcf,'Color',[1 1 1]);

h = filter_design_qs_unc(Hd,2*L+1,2*L+1);
%save('least_pth_filter_cof','h');
%Compute and plot the frequency response
[H,fx,fy] = freqz2(h,64,64,[1 1]);
Hmag = abs(H);
[Fx,Fy] = meshgrid(fx,fy);
%perspective plot
figure;
mesh(Fx(1:2:64,1:2:64),Fy(1:2:64,1:2:64),abs(H(1:2:64,1:2:64)));            % generate the perspective plot
colormap([0 0 0]);     % use black only
axis([-.5 .5 -.5 .5 0 1.2]);
xlabel('\it u \rm (c/px)'), ylabel('\it v \rm (c/px)');
set(gca,'ydir','reverse'); 
set(gca,'fontname','times')          %set the plot font to Times
set(gcf,'Color',[1 1 1]);