% This code has been adapted from 
%   http://www.seas.upenn.edu/~ese216/handouts/Chpt14_3DPlotTransferFunction.pdf
% This code is meant to facilitate the 3D visualization of the magnitude 
% plot of transfer function in s-domain. Also the its projection along 
% the frequency axis, or the Bode plot can also be plotted. 
%
%
%%%%  Transfer function  is  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
%               
%           s^2
% sys = --------------
%         s^2+2*s+5
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
clear all;close all;clc
  num = [1 0 0];
  den = [1 2 5];
  sys=tf(num,den);
  printsys(num,den)
  [z,p,k] = tf2zp(num,den)
  ltiview('pzmap',sys);
  
figure(1)
 opts = bodeoptions;
 opts.PhaseMatching = 'on';
 opts.Grid = 'on';
 [magni,phase,w] = bode(sys);
  bode(sys,opts)

mag=zeros(length(w));
 for i=1:length(w)
     mag(i) = 20*log10(magni(:,:,i)); %magnitude in dB
 end
 %
%%%%    FOR 3D VISUALIZATION     %%%%%

 [X,Y]=meshgrid(-20:0.2:20,-10:0.2:10); 
 F2=20*log10((X.^2+Y.^2)./(sqrt((X.^2-Y.^2+2*X+5).^2+(2*X.*Y+2*Y).^2)+eps));  
 
 %%%%    3D SURFACE PLOT OF THE SAME  TRANSFER FUNCTION    %%%%%%%%%%%%
  
 figure(2)
  s = surf(X,Y,F2);  %,'FaceAlpha',0.9 (transparencia)
  %s.EdgeColor = 'none';
  lim = caxis
  caxis([-10 10])
  
  %set(gca,'yscale','log')
  %set(gca,'zscale','log')
  axis([-10 10 -10 10 -10 20])
  xlabel('Parte Real')
  ylabel('Parte Imaginária')
  zlabel('Magnitude, em dB');
  title('3D Plot of Magnitude')
  
 figure(3)
  [X,Y]=meshgrid(-20:0.05:20,-100:0.05:100); 
  F2=20*log10((X.^2+Y.^2)./(sqrt((X.^2-Y.^2+2*X+5).^2+(2*X.*Y+2*Y).^2)+eps));  
  s = surf(X,Y,F2);  %,'FaceAlpha',0.9 (transparencia)
  %s.EdgeColor = 'none';
  lim = caxis
  caxis([-60 20])
  
  set(gca,'yscale','log')
  %set(gca,'zscale','log')
  axis([0 0.01 0 100 -60 20])
  view(90,0)
  xlabel('Parte Real')
  ylabel('Parte Imaginária')
  zlabel('Magnitude, em dB');
  title('3D Plot of Magnitude')