function spit(command) 

%SPIT System Plotting Interactive Tool
%     SPIT is a graphical user interface to plot the linear transfer
%     function of a Simulink system.  The interface is initialized and
%     a new window is opened for issuing commands upon calling SPIT.
%
%     The SPIT command window contains a number of text entries which may
%     be edited and buttons to generate a plot. With a Simulink model
%     open, click on the "Bode", "Nyquist", or "Nichols" button.
%
%     The controls let you customize the appearance of the plots.
%     Frequencies should be typed as logs (use -2 instead of 0.01).%
%
%     For Bode plots, you can constrain the plot ranges by setting the
%     magnitude and frequency ranges at the bottom of the window. You can
%     measure the transfer function between two output ports by setting 
%     the TF parameter.
%
%     SPIT can also create a Bode plot representing the transfer function
%     of a single element in the model.
%
%     SPIT will send its plots to a seperate plotting window and
%     subsequent plots will replace one another. Preserve a plot by
%     clicking on the "New Figure" button.
%
% Version 2.5 by Shanti Rao, 1999
% Version 2.0 by Jay Heefner,  1998
% Version 1.0 by Torrey Lyons, finished 10/24/94.

%  Interpret commands.

if nargin == 0 
	command = 'initialize'; 
end 

if isstr(command) 
	if strcmp(lower(command),'initialize') 
		command = 0; 
	elseif strcmp(lower(command),'newfig') 
		command = 1; 
	elseif strcmp(lower(command),'bode') 
		command = 2; 
	elseif strcmp(lower(command),'bodetf') 
		command = 3; 
	elseif strcmp(lower(command),'nyquist') 
		command = 4; 
	elseif strcmp(lower(command),'nichols') 
		command = 5; 
	end 
end 

%  Load stored data containing handles to GUI controls.

if command ~= 0 
	ui_handles = get(gcf,'UserData'); 
	in_nm = ui_handles(1);
	in_str = ui_handles(2);
	in_stp = ui_handles(3);
	in_num = ui_handles(4);
	in_inp = ui_handles(5);
	in_mch = ui_handles(6);
	in_dch = ui_handles(7);
	in_fg  = ui_handles(8);	
    plot_start = ui_handles(9);
	plot_stop = ui_handles(10);
	plot_min = ui_handles(11);
	plot_max = ui_handles(12);
    in_fixed = ui_handles(13);

end 

%---Command to initialize spit.

if command == 0 
%		figure screen position [ left bottom width height ]
	figure('Position',[ 358 56 300 300 ], 'NumberTitle','off', ...
	    'Name', 'System Plotting Interactive Tool'); 
	set(gca,'visible','off'); 

%****  GUI Object Creation  ****

% Create the frame, label, and control to specify the system.
	fr_nm = uicontrol('Units','normalized',...
                'Position',[ 0.02 0.60 0.96 0.39 ],... 
		'Style','frame');
% System name
	tx_nm = uicontrol('Units','normalized',...
                'Position',[ 0.05 0.90 0.15 0.06 ],...
		'String','System:',... 
		'Style','text','HorizontalAlignment','Left');
	in_nm = uicontrol('Units','normalized',...
                'Position',[ 0.30 0.90 0.40 0.07 ],...
		'Style','edit','HorizontalAlignment','Left');
	in_fixed = uicontrol('Units','normalized',...
                'Position',[ 0.71 0.90 0.10 0.07 ],...
		'Style','edit','HorizontalAlignment','Right');
	tx_fixed = uicontrol('Units','normalized',...
                'Position',[ 0.81 0.90 0.14 0.06 ],...
		'String','[s] step',... 
		'Style','text','HorizontalAlignment','Left');
% start, stop
	tx_str = uicontrol('Units','normalized',...
		'Position',[ 0.05 0.82 0.25 0.06 ],... 
		'String','Start Freq:',... 
		'Style','text','HorizontalAlignment','Left');
	in_str = uicontrol('Units','normalized',...
                'Position',[ 0.30 0.82 0.15 0.07 ],...
		'Style','edit','String','0','HorizontalAlignment','Left');
	tx_stp = uicontrol('Units','normalized',...
                'Position',[ 0.05 0.72 0.25 0.06 ],... 
		'String','Stop Freq:',...
		'Style','text','HorizontalAlignment','Left');
	in_stp = uicontrol('Units','normalized',...
                'Position',[ 0.30 0.72 0.15 0.07 ],...
		'Style','edit','String','3','HorizontalAlignment','Left');
% points, port
	tx_num = uicontrol('Units','normalized',...
                'Position',[ 0.55 0.82 0.25 0.06 ],... 
		'String','Points:',... 
		'Style','text','HorizontalAlignment','Left',...
		'Units','normalized');
	in_num = uicontrol('Units','normalized',...
                'Position',[ 0.80 0.82 0.15 0.07 ],...
		'Style','edit', 'String', 50,'HorizontalAlignment','Left');
	tx_inp = uicontrol('Units','normalized',...
                'Position',[ 0.55 0.72 0.25 0.06 ],... 
		'String','Test in port:',... 
		'Style','text','HorizontalAlignment','Left');
	in_inp = uicontrol('Units','normalized',...
                'Position',[ 0.80 0.72 0.15 0.07 ],...
		'String',1, 'Style','edit','HorizontalAlignment','Left');

% output divisor port
	tx_out = uicontrol('Units','normalized',...
                'Position',[ 0.50 0.62 0.25 0.06 ],... 
		'String','TF: <',... 
		'Style','text','HorizontalAlignment','Left');
	in_mch = uicontrol('Units','normalized',...
                'Position',[ 0.63 0.62 0.12 0.07 ],...
		'String','', 'Style','edit','HorizontalAlignment','Left');
	tx_out = uicontrol('Units','normalized',...
                'Position',[ 0.75 0.62 0.05 0.06 ],... 
		'String','/',... 
		'Style','text','HorizontalAlignment','Left');
	in_dch = uicontrol('Units','normalized',...
                'Position',[ 0.80 0.62 0.12 0.07 ],...
		'String','', 'Style','edit','HorizontalAlignment','Left');
	tx_out = uicontrol('Units','normalized',...
                'Position',[ 0.93 0.62 0.05 0.06 ],... 
		'String','>',... 
		'Style','text','HorizontalAlignment','Left');
% new figure
	bt_fg = uicontrol('Units','normalized',...
                'Position',[ 0.05 0.62 0.25 0.07 ],... 
		'String','New Figure', 'Style','pushbutton',... 
		'Callback', 'spit(''newfig'')');
	in_fg = uicontrol('Units','normalized',...
                'Position',[ 0.30 0.62 0.15 0.07 ],...
		'Style','edit');
% frame
	fr_nm = uicontrol('Units','normalized',...
                'Position',[ 0.02 0.50 0.96 0.09 ],... 
		'Style','frame');

%  Create plotting push buttons.
	tx_inp = uicontrol('Units','normalized',...
                'Position',[ 0.05 0.51 0.10 0.06 ],... 
		'String','Plot:',... 
		'Style','text','HorizontalAlignment','Left');
	bt_bod = uicontrol('Units','normalized',...
                'Position',[ 0.16 0.51 0.18 0.07 ],... 
		'String','Bode', 'Style','pushbutton',... 
		'Callback', 'spit(''bode'')');
	bt_bodtf = uicontrol('Units','normalized',...
                'Position',[ 0.36 0.51 0.18 0.07 ],... 
		'String','TF', 'Style','pushbutton',... 
		'Callback', 'spit(''bodetf'')');
	bt_nyq = uicontrol('Units','normalized',...
                'Position',[ 0.56 0.51 0.18 0.07 ],... 
		'String','Nyquist', 'Style','pushbutton',... 
		'Callback', 'spit(''nyquist'')');
	bt_nch = uicontrol('Units','normalized',...
                'Position',[ 0.76 0.51 0.18 0.07 ],... 
		'String','Nichols', 'Style','pushbutton',... 
		'Callback', 'spit(''nichols'')');

% Test of frame for bode plot inputs
% [ left bottom width height ]
	fr_bode= uicontrol('Units','normalized',...
                'Position',[ 0.02 0.20 0.96 0.29 ],... 
		'Style','frame');

% plot start/stop
	tx_str_frq= uicontrol('Units','normalized',...
                'Position',[ 0.30 0.22 0.20 0.06 ],... 
		'String','Frequency:',... 
		'Style','text','HorizontalAlignment','Left');
	plot_start = uicontrol('Units','normalized',...
                'Position',[ 0.50 0.22 0.15 0.07 ],...
		'Style','edit');
	plot_stop = uicontrol('Units','normalized',...
                'Position',[ 0.70 0.22 0.15 0.07 ],...
		'Style','edit');

% magnitude min/max
	tx_min_mag= uicontrol('Units','normalized',...
                'Position',[ 0.05 0.32 0.25 0.12 ],... 
		'String','Magnitude (dB):',... 
		'Style','text','HorizontalAlignment','Right');
	plot_min = uicontrol('Units','normalized',...
                'Position',[ 0.30 0.32 0.15 0.07 ],...
		'Style','edit');
	plot_max = uicontrol('Units','normalized',...
                'Position',[ 0.30 0.39 0.15 0.07 ],...
		'Style','edit');

% Help text at the bottom
	fr_bode= uicontrol('Units','normalized',...
                'Position',[ 0.02 0.00 0.96 0.19 ],... 
		'Style','frame');

helptext = [ ...
'Type in the log of the frequency range at the top and then ' ...
'click on a plot button. For Bode plots, you can enter the plot ' ... 
'ranges at the bottom. A transfer function (TF) measures the ratio ' ...
'of two output ports.' ];,

	tx_help= uicontrol('Units','normalized',...
                'Position',[ 0.05 0.02 0.90 0.16 ],... 
		'String', helptext, ...
		'Style','text','HorizontalAlignment','Left');

%  Save handles to test input data in UserData matrix.

set(gcf,'UserData',[in_nm in_str in_stp in_num in_inp in_mch in_dch ...  
                    in_fg plot_start plot_stop plot_min plot_max ...
                    in_fixed]);

%---Command to open new figure.

elseif command == 1 
	fignum = figure('Name','System Plotting Output'); 
	set(in_fg, 'String', fignum);
else

%---Commands to plot transfer functions in various ways.

%  Load user entered data.
	name   = get(in_nm, 'String');
	start  = str2num(get(in_str, 'String'));
	stop   = str2num(get(in_stp, 'String'));
	numpt  = str2num(get(in_num, 'String'));
	inport = str2num(get(in_inp, 'String'));
	fignum = str2num(get(in_fg,  'String'));
	mch    = str2num(get(in_mch, 'String'));
	dch    = str2num(get(in_dch, 'String'));
    Ts = str2num( get(in_fixed,'String'));

	if (isempty(get(in_dch,'String')))
		if (isempty(get(in_mch,'String')))
			istransfer = 0;
		else
			istransfer = 3;
		end
	else
		if (isempty(get(in_mch,'String')))
			istransfer = 1;
		else 
			istransfer = 2;
		end
    end

    if (command == 2)
        istransfer = 0;
    end
% plot user data
	start_plot = str2num(get(plot_start, 'String'));
	stop_plot = str2num(get(plot_stop, 'String'));
	min_plot = str2num(get(plot_min, 'String'));
	max_plot = str2num(get(plot_max, 'String'));

%  Prepare for plotting commands.
	if (isempty(name)) 
	    name = gcs; 
	    set(in_nm,'String',name);
	end

	if (isempty(name))
	    disp(' ')
	    error('There is no active Simulink system.');
	    return
	end

	if (isempty(get(in_fg,  'String')))
	    fignum = figure('Name', 'System Plotting Output'); 
	    set(in_fg, 'String', fignum);
	end


    if isempty(Ts) || ~Ts
        [A,B,C,D] = linmod(name);
        isDiscrete = 0;
        sys = ss(A,B,C,D);
    else
        [A,B,C,D] = dlinmod(name,Ts);
        isDiscrete = 1;
        sys = ss(A,B,C,D,Ts);
        stop = min(stop,log10(1/Ts/2));
    end
    
	simple = 0;
	if (isempty(start) || isempty(stop))
	    simple = 1;
	else
 	    w = logspace(start + log10(2*pi), stop + log10(2*pi), numpt);
    end
%---Bode Plot
	if (command == 2 || command == 3)

        if (command == 2)
            set(fignum, 'Name', [name ': Bode Plot']);
        else
            set(fignum, 'Name', [name ': Transfer Function']);
        end
	figure(fignum);

	% set axis limits to evaluation limits if requests are out of range
	if (isempty(start_plot) | (start_plot < start) | (start_plot > stop))
             start_plot = start;
	end
	if ((isempty(stop_plot))|(stop_plot > stop) | (stop_plot < start))
            stop_plot = stop;
    end

	  if (simple == 1)                        
 		bode(sys);
	  else
 		[mag, phase, w] = bode(sys,w);
        
%         [H, w] = freqresp(tf,inport,w);
%         mag = abs(H(:,:));
%         phase = 180/pi * arg(H(:,:));

		if (istransfer == 0)
		   magdb = 20*log10(mag(inport,:));
           phase = phase(inport,:);
		elseif (istransfer == 1)
   error(['Sorry, you must provide a numerator in order to measure ' ...
          'a single channel transfer function.']);
		elseif (istransfer == 2)
		   magdb = 20*log10(mag(:,mch) ./ mag(:,dch));
		   phase = phase(:,mch) - phase(:,dch);
		elseif (istransfer == 3)
		   magdb = 20*log10(mag(:,mch));
		   phase = phase(:,mch);
		end

		v = w / (2*pi);

       rows = size(magdb);

       if (rows(1,1) == 1)
	    if (isempty(min_plot))  min_plot = min(magdb); end
		 if (isempty(max_plot))  max_plot = max(magdb); end
		 min_phase = min(phase);
         max_phase = max(phase);
        else
         if (isempty(min_plot))  min_plot = min(min(magdb)); end
         if (isempty(max_plot))  max_plot = max(max(magdb)); end
         min_phase = min(min(phase));
         max_phase = max(max(phase));
        end

		sp = 10^start_plot;
        ep = 10^stop_plot; 

		subplot(2,1,1); %plot magnitude on top plot
		semilogx(v, magdb);
        axis([sp ep min_plot max_plot]);
         
        gridcolor([.6 .6 .7],[.6 .6 .7],[.8 .8 .9],[.8 .8 .9],'-',':');
		       
		ylabel('Magnitude [dB]');
		xlabel('Frequency [Hz]');
        
		subplot(2,1,2); %plot phase on bottom plot
		semilogx(v, phase);
		axis([sp ep min_phase max_phase]);
        
        gridcolor([.6 .6 .7],[.6 .6 .7],[.8 .8 .9],[.8 .8 .9],'-',':');
        
		ylabel('Phase [deg]');
		xlabel('Frequency [Hz]');

        end

%---Nyquist Plot
	elseif command == 4 
	    if (isempty(fignum))
		fignum = figure('Name', [name ': Nyquist Plot']); 
		set(in_fg, 'String', fignum);
	    else
		set(fignum, 'Name', 'Nyquist Plot');
		figure(fignum);
	    end
	    if (simple)
		nyquist(A,B,C,D,inport);
	    else
		nyquist(A,B,C,D,inport,w);
	    end

%---Nichols Plot
	elseif command == 5 
	    if (isempty(fignum))
		fignum = figure('Name', [name ':Nichols Plot']); 
		set(in_fg, 'String', fignum);
	    else
		set(fignum, 'Name', 'Nichols Plot');
		figure(fignum);
	    end
	    if (simple)
%		nichols(A,B,C,D,inport);% no nichols commands in student edition
		[mag, phase] = bode(A,B,C,D,inport);
		magdb = 20*log10(mag);
 		plot(phase,magdb);
		xlabel('degrees');
		ylabel('dB');
		grid on;
	    else
%		nichols(A,B,C,D,inport,w); % no nichols commands in student edition
		[mag, phase, w] = bode(A,B,C,D,inport,w);
		magdb = 20*log10(mag);
 		plot(phase,magdb);
		xlabel('degrees');
		ylabel('dB');
		grid on;
	    end

%---Unrecognized command
	else 
	    error('Error: spit.m called with incorrect command.') 
	end
end 
end

function gridcolor(majorX, majorY, minorX, minorY, sMajor, sMinor)
if nargin < 2, majorY = majorX; end
if nargin < 3, minorX = majorX; end
if nargin < 4, minorY = minorX; end
if nargin < 5, sMajor = '-'; end
if nargin < 6, sMinor = sMajor; end

ax1 = gca;   %# get a handle to first axis

%# create a second transparent axis, same position/extents, same ticks and labels
ax2 = copyobj(ax1,gcf);
ax3 = copyobj(ax1,gcf);

delete(get(ax2,'Children'));
delete(get(ax3,'Children'));

set(ax2, 'Color','none', 'Box','off','YTickLabel',[],'YTickLabel',[],...
    'GridLineStyle', sMajor,...
    'XGrid','on','YGrid','on',...
    'XMinorGrid','off','YMinorGrid','off',...
    'XColor',majorX,'YColor',majorY);
set(ax3,'Box','off','YTickLabel',[],'YTickLabel',[],...
    'MinorGridLineStyle',sMinor,...
    'XGrid','off','YGrid','off',...
    'XMinorGrid','on','YMinorGrid','on',...
    'XColor',minorX,'YColor',minorY);

set(ax1, 'Color','none', 'Box','on')

handles = [ax3; ax2; ax1];
c = get(gcf,'Children');
for i=1:length(handles)
    c = c(find(c ~= handles(i)));
end
set(gcf,'Children',[c; flipud(handles)]);

linkaxes([ax1 ax2 ax3]);
end







