No BSD License  

Highlights from
Adventures in Bifurcations in Chaos++

Adventures in Bifurcations in Chaos++

by

 

24 Mar 2004 (Updated )

Simulation of Chua's Oscillator in MATLAB with Poincare Section analysis.

poincare_popup.m
function f = poincare_popup(decide)
% Run Poincar Application
%Peadar Forbes
global h a Timeseries approx_1 poincare_points V2vector plane H_ABC_pp vertex_matrix save1 line_drawn visible_line af2 ah aponc exp_data ab abponc

switch decide
    case 0 % Set up Poincar GUI 
        %Background Color
        back_clr = [0.7529 0.7929 0.97];
        % Figure
        f(1) = figure('Color',[0.83137 0.81569 0.78431], ...
            'Color',back_clr, ...
            'Name','Attractor Poincar Section', ...
            'NumberTitle','off', ...
            'Units','normalized','Position',[.01 .044 .974 .897], ...
            'Resize','off', ...
            'Menubar','none',...
            'Visible','off', ...
            'WindowStyle','normal');
        
        h(1) = axes('Parent',f(1), ...
            'Color',[1 1 1],'CameraUpVector',[0 1 0],...
            'XColor',[0 0 0],'YColor',[0 0 0],'ZColor',[0 0 0],...
            'Units','normalized','Position',[0.040 0.07 0.5 0.68]);
        
        h(2) = axes('Parent',f(1), ...
            'Color',[1 1 1],'CameraUpVector',[0 1 0],...
            'XColor',[0 0 0],'YColor',[0 0 0],'ZColor',[0 0 0],...
            'Units','normalized','Position',[0.64 0.095 0.35 0.40]);
        
        h(3) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame plane big
            'Units','normalized','Position',[0.01 0.919 0.21 0.04],'Style','frame');
        
        h(4) = axes('Parent',f(1), ...
            'Color',[1 1 1],'CameraUpVector',[0 1 0],...
            'XColor',[0 0 0],'YColor',[0 0 0],'ZColor',[0 0 0],...
            'Units','normalized','Position',[0.64 0.58 0.35 0.40]);
        
        h(5) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame Plane small
            'Units','normalized','Position',[0.01 0.959 0.13 0.03],'Style','frame');
        
        h(6) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame slope big
            'Units','normalized','Position',[0.22 0.759 0.238 0.20],'Style','frame');
        
        h(7) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame Slope
            'Units','normalized','Position',[0.22 0.96 0.08 0.03],'Style','frame');
        
        h(8) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame Options small
            'Units','normalized','Position',[0.01 0.888 0.13 0.03],'Style','frame');
        
        h(9) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame Options big
            'Units','normalized','Position',[0.01 0.759 0.21 0.13],'Style','frame');
        
        h(10) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame Pusbuttons big
            'Units','normalized','Position',[0.458 0.759 0.10 0.14],'Style','frame');
        
        
        % Recalculate due to changed step_size or changed plane
        H_ABC_pp(502) = uicontrol('Parent',f(1),'Style','pushbutton',...
            'Units','normalized','Position',[0.468 0.765 0.08 0.03125],...
            'ForegroundColor',[0 0 0],...
            'FontWeight','bold',...
            'String','Recalculate',...
            'Visible','on',...
            'CallBack','poincare_popup(1)');
        
        % Save Figures
        H_ABC_pp(504) = uicontrol('Parent',f(1),'Style','pushbutton',...
            'Units','normalized','Position',[0.468 0.815 0.08 0.03125],...
            'ForegroundColor',[0 0 0],...
            'FontWeight','bold',...
            'String','Screenshot',...
            'Visible','on',...
            'CallBack','poincare_popup(5)');
        
        % Calculate Slope - Pushbutton
        H_ABC_pp(511) = uicontrol('Parent',f(1),'Style','pushbutton',...
            'Units','normalized','Position',[0.24 0.765 0.1 0.03125],...
            'ForegroundColor',[0 0 0],...
            'FontWeight','bold',...
            'String','Calculate Slope',...
            'Visible','on',...
            'CallBack','poincare_popup(10)');
        
        
        
        % Clear Slope - Pushbutton
        H_ABC_pp(510) = uicontrol('Parent',f(1),'Style','pushbutton',...
            'Units','normalized','Position',[0.365 0.765 0.07 0.03125],...
            'ForegroundColor',[0 0 0],...
            'FontWeight','bold',...
            'String','Clear Slope',...
            'Visible','on',...
            'CallBack','poincare_popup(1)');
        
                
        % Load Experimental - Pushbutton
        H_ABC_pp(513) = uicontrol('Parent',f(1),'Style','pushbutton',...
            'Units','normalized','Position',[0.92 0.5 0.07 0.03125],...
            'ForegroundColor',[0 0 0],...
            'FontWeight','bold',...
            'String','Load Exp',...
            'Visible','on',...
            'CallBack','poincare_popup(11)');

        
        % Export Figures
        H_ABC_pp(515) = uicontrol('Parent',f(1),'Style','pushbutton',...
            'Units','normalized','Position',[0.468 0.86 0.08 0.03125],...
            'ForegroundColor',[0 0 0],...
            'FontWeight','bold',...
            'String','Export',...
            'Visible','on',...
            'CallBack','poincare_popup(15)');
       
 
        
        aponc(1) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Frame  (View angle)
            'Units','normalized','Position',[0.65 0.005 0.27 0.029],'Style','frame');
        
        aponc(2) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...%View Angle Text
            'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.653 0.009 0.08 0.02],'String','View Angle:','Style','text');
        
        
        
        
        aponc(4) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% I3-V2 Text
            'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.74 0.009 0.03 0.02],'String','I3-V2','Style','text');
        
        
        aponc(5) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% I3-V1 Text
            'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.8 0.009 0.03 0.02],'String','I3-V1','Style','text');
        
        aponc(6) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% V2-V1 Text
            'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.86 0.009 0.05 0.02],'String','V2-V1','Style','text');
        
        a(11) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Plane Equation Text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'FontWeight','bold','ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.02 0.965 0.11 0.020],'String','Plane Equation:','Style','text');
        
        a(23) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Options Text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'FontWeight','bold','ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.02 0.893 0.11 0.020],'String','Options:','Style','text');
        
        a(12) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% SlopeText
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'FontWeight','bold','ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.967 0.05 0.021],'String','Slope:','Style','text');
        
        
        a(13) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Step Size Reduction: Text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.02 0.827 0.13 0.06],'String','Step Size Reduction:','Style','text');
        
                
        a(22) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Ramp up text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.02 .825 0.1 0.035],'String','Ramp up:','Style','text');
        
        a(24) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Use smart Ramp Up Text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.039 0.79 0.13 0.035],'String','Use Smart Ramp-up','Style','text');
        
        a(25) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% View "Red" Points Text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.039 0.76 0.13 0.035],'String','View "Red" Points','Style','text');
        
        a(14) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Reduce Step Size input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.155 0.862 0.035 0.025],'Style','edit','String','10');
        
        a(15) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% x coefficient input box 0.70 0.78 0.025 0.025
            'HorizontalAlignment','left','Units','normalized','Position',[0.02 0.925 0.025 0.025],'Style','edit','String','1');
        
        a(24) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Remove Ramp-up input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.155 0.835 0.035 0.025],'Style','edit','String','1');
        
        a(19) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% x coefficient text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.045 0.925 0.021 0.025],'String','x +','Style','text');
    
        
        a(16) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% y coefficient input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.07 0.925 0.025 0.025],'Style','edit','String','0');
        
        a(20) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% y coefficient text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.095 0.925 0.021 0.025],'String','y +','Style','text');
        
        a(17) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% z coefficient input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.12 0.925 0.025 0.025],'Style','edit','String','0');
        
        a(23) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% z coefficient text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.145 0.925 0.021 0.025],'String','z +','Style','text');
        
        a(18) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% d coefficient input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.17 0.925 0.025 0.025],'Style','edit','String','-.00137');
        
        a(21) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% = 0 text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.195 0.925 0.021 0.025],'String','= 0','Style','text');


        
        
        %Slope Calculation GUI components:
        
        af2(2) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...%Linear Region of Graph text
            'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.92 0.15 0.03],'String','Linear Region of Graph: ','Style','text');
        
        af2(3) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Vxn min value text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.905 0.1 0.02],'String','V2(xn) min:','Style','text');
        
        af2(4) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Vxmin value input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.295 0.901 0.025 0.025],'Style','edit','String','0.2');
        
        
        af2(5) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Vxn max value text
            'FontName','Helvetica','HorizontalAlignment','left',...
            'FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.35 0.905 0.1 0.02],'String','V2(xn) max:','Style','text');
        
        af2(6) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...% Vxmax value input box
            'HorizontalAlignment','left','Units','normalized','Position',[0.42 0.901 0.025 0.025],'Style','edit','String','0.6');
        

        
        af2(9) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...%Results: text
            'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.87 0.1 0.02],'String','Results: ','Style','text');
        
        af2(10) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...%Least Squares Approximation: text
            'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.845 0.22 0.02],'String','Least Squares Approximation: ','Style','text');
        
        af2(11) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...%Calculated Slope: text
            'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.82 0.21 0.02],'String','Calculated Slope: ','Style','text');
        
        af2(12) = uicontrol('Parent',f(1),'BackgroundColor',[1 1 1],...%Percentage difference: text
            'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
            'Units','normalized','Position',[0.23 0.795 0.2 0.02],'String','% Diff: ','Style','text');
        
        
        
        

        
        %RADIOBUTTONS-----------------------------------------------
        
        %I3-V2 Radio Button
        H_ABC_pp(506) = uicontrol('Style','radiobutton','Min',0,'Max',1,... % Radio
            'Units','normalized','Position',[0.775 0.006 0.014 0.024],...
            'BackgroundColor',[1 1 1],'Visible','on',...
            'CallBack','poincare_popup(7)'); 
        
        %I3-V1 Radio Button
        H_ABC_pp(507) = uicontrol('Style','radiobutton','Min',0,'Max',1,... % Radio
            'Units','normalized','Position',[0.835 0.006 0.014 0.024],...
            'BackgroundColor',[1 1 1],'Visible','on',...
            'Value',1,...
            'CallBack','poincare_popup(8)'); 
        
        
        %V2-V1 Radio Button
        H_ABC_pp(508) = uicontrol('Style','radiobutton','Min',0,'Max',1,... % Radio
            'Units','normalized','Position',[0.90 0.006 0.014 0.024],...
            'BackgroundColor',[1 1 1],'Visible','on',...
            'CallBack','poincare_popup(9)'); 
        
        
        
        
        %CHECKBOXES--------------------------------------------------------
                
        %Smart Ramp-up checkbox
        H_ABC_pp(510) = uicontrol('Style','checkbox','Min',0,'Max',1,... % Checkbox
            'Units','normalized','Position',[0.02 0.80 0.014 0.02],...
            'BackgroundColor',[1 1 1],'Visible','on',...
            'CallBack','poincare_popup(1)');  
        
        
        
        %Radio button to select which direction through the plane you want
        %to plot in the poincare map(top right) window.
        H_ABC_pp(503) = uicontrol('Style','checkbox','Min',0,'Max',1,... % Checkbox
            'Units','normalized','Position',[0.02 0.77 0.014 0.02],...
            'BackgroundColor',[1 1 1],'Visible','on',...
            'CallBack','poincare_popup(1)');
        %CHECKBOXES END--------------------------------------------
        
        
        
        set(f(1),'Visible','on');
        set(af2(9:12),'Visible','off'); 
        % Read in the variables passed to poincar_popup that were temp saved to disk!!
        S = load('data\temp.mat');
        Timeseries = S.Timeseries;
        sys_variables = S.sys_variables;
        integ_variables = S.integ_variables;
        view_angle = S.view_angle;
        eigen_value = S.eigen_value;
 
        plane = [str2num(get(a(15),'String')),str2num(get(a(16),'String')),str2num(get(a(17),'String')),str2num(get(a(18),'String'))];
        integ_variables(6) = str2num(get(a(14),'String')); %Gets value for the step_size_divisor
        poincare_points_a = poincare(Timeseries,plane,sys_variables,integ_variables);
        
        poincare_points = poincare_points_a(:,1:4);
        approx_1 = poincare_points_a(:,5:8);
        
        
        if isreal(poincare_points)
            axes(h(2)); %Point to Poincare Section - bottom right
            for i=1:length(poincare_points(:,1))
                if poincare_points(i,4) == 1
                    plot3(poincare_points(i,1),poincare_points(i,2),poincare_points(i,3),'or');
                    hold on;
                else
                    plot3(poincare_points(i,1),poincare_points(i,2),poincare_points(i,3),'ob');
                    hold on;
                end
            end  
            hold off;
            xmax = max(poincare_points(:,1));
            xmin = min(poincare_points(:,1));
            zmax = max(poincare_points(:,3));
            zmin = min(poincare_points(:,3));
            ymax = max(poincare_points(:,2));
            ymin = min(poincare_points(:,2));
            
            if xmax > 0
                xmax=xmax*1.1;
            else
                xmax=xmax*.9;
            end
            if xmin > 0
                xmin=xmin*.9;
            else
                xmin=xmin*1.1;
            end
            if zmax > 0
                zmax=zmax*1.1;
            else
                zmax=zmax*.9;
            end
            if zmin > 0
                zmin=zmin*.9;
            else
                zmin=zmin*1.1;
            end
            if ymin > 0
                ymin=ymin*.9;
            else
                ymin=ymin*1.1;
            end
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),xmin,xmax,zmin,zmax,ymin,ymax);
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',.5,'EdgeAlpha',0.1);
            set(h(2),'CameraUpVector',[plane(1) plane(2) plane(3)]);
            xlabel('I3');
            ylabel('V2');
            zlabel('V1');
            view(0,0);axis tight;
            rotate3d on;
            
            axes(h(4)); % Point to Poincare Map - top right 
            
            %This separates the poincare map points into two bins - one for each 
            %direction that they cross the plane in
            
            poin = poincare_points;
            
            bindex=1;
            aindex=1;
            
            
            for i=1:length(poin(:,1))           
                if poin(i,4)==0  %this finds what direction they crossed in 
                    a_poin(aindex,:) = poin(i,:);
                    aindex = aindex+1;
                else
                    b_poin(bindex,:) = poin(i,:);
                    bindex = bindex+1;
                end
                
            end
            
            if(exist('a_poin') && length(a_poin(:,1)) >1) %if a_poin exists
              
                %poin_xnplusone is the value at time t+1
                %this will be plotted against poin_xn, the value at time t
                
                poin_xnplusone(1:length(a_poin(:,1))-1,:) = a_poin(2:length(a_poin(:,1)),:);
                poin_xn(1:length(a_poin(:,1))-1,:) = a_poin(1:length(a_poin(:,1))-1,:);
                if length(poin_xn)~=0
                    I3vector = [poin_xn(:,1),poin_xnplusone(:,1)]; %Sorts them into vectors that can be plotted
                    V2vector = [poin_xn(:,2),poin_xnplusone(:,2)];
                    V1vector = [poin_xn(:,3),poin_xnplusone(:,3)];
                    
                    plot(V2vector(:,1),V2vector(:,2),'.') %Plots the V2 vector at xn+1 against xn
                    xlabel('V2(xn)');
                    ylabel('V2(xn+1)');
                    clear poin;
                else
                    text(0.5,0.5,'Not enough points to plot');
                    axis tight
                end
            end
            axes(h(1)); % Point to main stage - i.e. TimeSeries Trajectory Plotting            
            plot3(Timeseries(1,:),Timeseries(2,:),Timeseries(3,:),'Color','k');
            xmin=min(Timeseries(1,:));
            xmax=max(Timeseries(1,:));
            zmin=min(Timeseries(3,:));
            zmax=max(Timeseries(3,:));
            
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),(xmin+xmin/10),(xmax+xmax/10),(zmin+zmin/10),(zmax+zmax/10),(ymin+ymin/10),(ymax+ymax/10));
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',0.7,'EdgeAlpha',0.1);
            %text(,'  E^{c}(0)','Color','m','FontWeight','bold');

            xlabel('I3');
            ylabel('V2');
            zlabel('V1');
            view(view_angle);%axis off;
            rotate3d(h(1));
            
            

        else
           
            % Read in warning image and place it Small Axes
            warning off all % Suppress Warning
            axes(h(2));
            imagen=imread('etc\warning.png');
            image(imagen);
            axis off;
            warning on all % Turn back on warnings
            axes(h(1)); % Point to main stage            
            plot3(Timeseries(1,:),Timeseries(2,:),Timeseries(3,:),'Color','k');

            xmin=min(Timeseries(1,:));
            xmax=max(Timeseries(1,:));
            zmin=min(Timeseries(3,:));
            zmax=max(Timeseries(3,:));
            ymin=min(Timeseries(2,:));
            ymax=max(Timeseries(2,:));
            
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),(xmin+xmin/10),(xmax+xmax/10),(zmin+zmin/10),(zmax+zmax/10),(ymin+ymin/10),(ymax+ymax/10));
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',1,'EdgeAlpha',0.1);
            view(view_angle);%axis off;
            rotate3d(h(1));
            
            axes(h(4));
            warning off all % Suppress Warning
            imagen=imread('etc\warning.png');
            image(imagen);
            axis off;
            warning on all % Turn back on warnings
            
        end
        set(H_ABC_pp(503),'Visible','on');
        set(H_ABC_pp(503),'Value',0);
        line_drawn=0;
        
    case 1
        %set(af2(8:12),'Visible','off'); %setting results in slope popup to invisible
        clear('poincare_points');
        S = load('data\temp.mat');
        Timeseries = S.Timeseries;
        sys_variables = S.sys_variables;
        plane = [str2num(get(a(15),'String')),str2num(get(a(16),'String')),str2num(get(a(17),'String')),str2num(get(a(18),'String'))];
        integ_variables = S.integ_variables;
        eigen_value = S.eigen_value;
        %set(a(22),'String',[get(a(15),'String'),'x+',get(a(16),'String'),'y+',get(a(17),'String'),'z+',get(a(18),'String'),'=0']);
        plane = [str2num(get(a(15),'String')),str2num(get(a(16),'String')),str2num(get(a(17),'String')),str2num(get(a(18),'String'))];
        
        integ_variables(6) = str2num(get(a(14),'String')); %Gets the step_size_divisor variable from user
        
        remove_ramp_up = str2num(get(a(24),'String'));
        Timeseries=Timeseries(:,remove_ramp_up:length(Timeseries));

        poincare_points_a = poincare(Timeseries,plane,sys_variables,integ_variables);

        poincare_points = poincare_points_a(:,1:4);
        approx_1 = poincare_points_a(:,5:8);
        
        if isreal(poincare_points)
            axes(h(2)); % Point to Poincare Section bottom right
            
            for i=1:length(poincare_points(:,1)) %plots the points in the window, red for one direction, blue for the other
                if poincare_points(i,4) == 1
                    plot3(poincare_points(i,1),poincare_points(i,2),poincare_points(i,3),'or');
                    hold on;
                else
                    plot3(poincare_points(i,1),poincare_points(i,2),poincare_points(i,3),'ob');
                    hold on;
                end
            end  
            hold off;   
            
            xmax = max(poincare_points(:,1));
            xmin = min(poincare_points(:,1));
            ymin = min(poincare_points(:,2));
            ymax = max(poincare_points(:,2));
            zmax = max(poincare_points(:,3));
            zmin = min(poincare_points(:,3));
            if xmax > 0
                xmax=xmax*1.1;
            else
                xmax=xmax*.9;
            end
            if xmin > 0
                xmin=xmin*.9;
            else
                xmin=xmin*1.1;
            end
            if zmax > 0
                zmax=zmax*1.1;
            else
                zmax=zmax*.9;
            end
            if zmin > 0
                zmin=zmin*.9;
            else
                zmin=zmin*1.1;
            end
            if ymin > 0
                ymin=ymin*.9;
            else
                ymin=ymin*1.1;
            end
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),xmin,xmax,zmin,zmax,ymin,ymax);
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',.5,'EdgeAlpha',0.1);
            set(h(2),'CameraUpVector',[plane(1) plane(2) plane(3)]);
            xlabel('I3');
            ylabel('V2');
            zlabel('V1');
            
            if get(H_ABC_pp(506),'Value') == 1
                view(0,90);
                axis tight;
            elseif get(H_ABC_pp(507),'Value') == 1
                view(0,0);  
                axis tight;
            elseif get(H_ABC_pp(508),'Value') == 1
                view(90,0);
                axis tight;
            end
            
            axes(h(4)); % Point to Poincare Map - top right 
            
            %This separates the poincare map points into two bins - one for each 
            %direction that they cross the plane in
            
            poin = poincare_points;
            
            bindex=1;
            aindex=1;
            for i=1:length(poin(:,1))           
                if poin(i,4)==0  %this finds what direction they crossed in 
                    a_poin(aindex,:) = poin(i,:);
                    aindex = aindex+1;
                else
                    b_poin(bindex,:) = poin(i,:);
                    bindex = bindex+1;
                end
                
            end
            
            if(exist('a_poin') && length(a_poin(:,1)) >1) %if a_poin exists
                %poin_xnplusone is the value at time t+1
                %this will be plotted against poin_xn, the value at time t
                
                if get(H_ABC_pp(503),'Value') == 0
                    
                    poin_xnplusone(1:length(a_poin(:,1))-1,:) = a_poin(2:length(a_poin(:,1)),:);
                    poin_xn(1:length(a_poin(:,1))-1,:) = a_poin(1:length(a_poin(:,1))-1,:);
                    

                        I3vector = [poin_xn(:,1),poin_xnplusone(:,1)]; %Sorts them into vectors that can be plotted
                        V2vector = [poin_xn(:,2),poin_xnplusone(:,2)];
                        V1vector = [poin_xn(:,3),poin_xnplusone(:,3)];
                        plot(V2vector(:,1),V2vector(:,2),'.') %Plots the V2 vector at xn+1 against xn


              
                    
                    
                else
                    poin_xnplusone(1:length(b_poin(:,1))-1,:) = b_poin(2:length(b_poin(:,1)),:);
                    poin_xn(1:length(b_poin(:,1))-1,:) = b_poin(1:length(b_poin(:,1))-1,:);

                        I3vector = [poin_xn(:,1),poin_xnplusone(:,1)]; %Sorts them into vectors that can be plotted
                        V2vector = [poin_xn(:,2),poin_xnplusone(:,2)];
                        V1vector = [poin_xn(:,3),poin_xnplusone(:,3)];
                        
                        plot(V2vector(:,1),V2vector(:,2),'.') %Plots the V2 vector at xn+1 against xn
                end
            end
            xlabel('V2(xn)');
            ylabel('V2(xn+1)');
           
            
            clear V2vector;
            clear poin_xnplusone;
            clear poin_xn;
            clear poin;
            
            
            axes(h(1)); % Point to main stage            
            plot3(Timeseries(1,:),Timeseries(2,:),Timeseries(3,:),'Color','k');
            xlabel('I3');
            ylabel('V2');
            zlabel('V1');
            
            xmin=min(Timeseries(1,:));
            xmax=max(Timeseries(1,:));
            zmin=min(Timeseries(3,:));
            zmax=max(Timeseries(3,:));
            ymin=min(Timeseries(2,:));
            ymax=max(Timeseries(2,:));
            
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),(xmin+xmin/10),(xmax+xmax/10),(zmin+zmin/10),(zmax+zmax/10),(ymin+ymin/10),(ymax+ymax/10));
            %            vertex_matrix = [xmin 0 zmin; xmin 0 zmax; xmax 0 zmin; xmax 0 zmax]

            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',.8,'EdgeAlpha',0.1);
            %view(view_angle);%axis off;
            rotate3d(h(1));
            
             
          
            

            
            
        else
            % Read in warning image and place it Small Axes

            warning off all % Suppress Warning
            axes(h(2));
            imagen=imread('etc\warning.png');
            image(imagen);
            axis off;
            warning on all % Turn back on warnings
            
            
            
            axes(h(4));
            warning off all % Suppress Warning
            imagen=imread('etc\warning.png');
            image(imagen);
            axis off;
            warning on all % Turn back on warnings
            axes(h(1)); % Point to main stage            
            plot3(Timeseries(1,:),Timeseries(2,:),Timeseries(3,:),'Color','k');
            xlabel('Current');
            xmin=min(Timeseries(1,:));
            xmax=max(Timeseries(1,:));
            zmin=min(Timeseries(3,:));
            zmax=max(Timeseries(3,:));
            ymin=min(Timeseries(2,:));
            ymax=max(Timeseries(2,:));
            
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),(xmin+xmin/10),(xmax+xmax/10),(zmin+zmin/10),(zmax+zmax/10),(ymin+ymin/10),(ymax+ymax/10));
            %vertex_matrix = [xmin 0 zmin; xmin 0 zmax; xmax 0 zmin; xmax 0 zmax]
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',1,'EdgeAlpha',0.1);
            %view(view_angle);%axis off;
            axis on;
            rotate3d(h(1));
            

        end
        line_drawn=0;
        
    
        
     %When toggling between which points, red or blue to show in poincare map   
    case 3
        
        %This separates the poincare map points into two bins - one for each 
        %direction that they cross the plane in
        
        poin = poincare_points;
        
        bindex=1;
        aindex=1;
        for i=1:length(poin)           
            if poin(i,4)==0  %this finds what direction they crossed in,0 means the blue points
                a_poin(aindex,:) = poin(i,:);
                aindex = aindex+1;
            else
                b_poin(bindex,:) = poin(i,:);
                bindex = bindex+1;
            end
            
        end
     

        if get(H_ABC_pp(503),'Value') == 1 % If (Plot Red points in Poincare map) radio selected
              
            %poin_xnplusone is the value at time t+1
            %this will be plotted against poin_xn, the value at time t
            
            poin_xnplusone(1:length(b_poin(:,1))-1,:) = b_poin(2:length(b_poin(:,1)),:);
            poin_xn(1:length(b_poin(:,1))-1,:) = b_poin(1:length(b_poin(:,1))-1,:);
            
            
        else
            
            poin_xnplusone(1:length(a_poin(:,1))-1,:) = a_poin(2:length(a_poin(:,1)),:);
            poin_xn(1:length(a_poin(:,1))-1,:) = a_poin(1:length(a_poin(:,1))-1,:);
            
        end
            I3vector = [poin_xn(:,1),poin_xnplusone(:,1)]; %Sorts them into vectors that can be plotted
            V2vector = [poin_xn(:,2),poin_xnplusone(:,2)];
            V1vector = [poin_xn(:,3),poin_xnplusone(:,3)];
            axes(h(4));
            plot(V2vector(:,1),V2vector(:,2),'.') %Plots the V2 vector at xn+1 against xn
            axes(h(1));
            axis on;
            rotate3d(h(1));
            clear poin_xnplusone;
            clear poin_xn;
            clear V2vector;
            
            
            %Save the figures - code in saving.m file 
        case 4
            savefile = 'data\temp2.mat';
            save(savefile,'poincare_points','plane')
          
            % Add etc directory to Matlab's Path: so as to call poincare_popup!
            saving(0,h,Timeseries,vertex_matrix,plane);
         
            %Screenshot
        case 5
            home = pwd;

            
            
            [filename, pathname, filterindex] = uiputfile( ...
                {
                '*.eps','Encapsulated Postscript (*.eps)'; ...
                    '*.jpg','JPEG image (*.jpg)'; ...
                    '*.emf','Enhanced metafile (*.emf)'}, ...
                'Save Figure as');
            
            if filename ~= 0 % Yes Save figure
                cd(pathname);
                if filterindex == 1
                    extension = 'eps';
                elseif filterindex == 2
                    extension = 'jpg';
                elseif filterindex == 3
                    extension = 'emf';
                end
                filename1 = [filename,'_colour'];
                saveas(h(1),filename,extension);
                saveas(h(1),filename1,'epsc');
                cd(home);
            else
                cd(home);
            end
            

            
            
            
            
            %Toggling the View angle on the Poincare Section
        case 7
            if get(H_ABC_pp(506),'Value') == 1
                
                set(H_ABC_pp(507),'Value',0);
                set(H_ABC_pp(508),'Value',0);
                
                axes(h(2)); %Point to poincare section-bottom right
                view(0,90);
                axis tight;
                
            end
            
            
            %Toggling the View angle on the Poincare Section
        case 8       
            if get(H_ABC_pp(507),'Value') == 1
                
                set(H_ABC_pp(506),'Value',0);
                set(H_ABC_pp(508),'Value',0);
                
                axes(h(2)); %Point to poincare section-bottom right
                view(0,0); 
                axis tight;
            end
            
            
            %Toggling the View angle on the Poincare Section  
        case 9
            
            if get(H_ABC_pp(508),'Value') == 1
                
                set(H_ABC_pp(506),'Value',0);
                set(H_ABC_pp(507),'Value',0);
                
                axes(h(2)); %Point to poincare section-bottom right
                view(90,0); 
                axis tight;
            end
            
        %Draw least squares approximation on poincare map - top right    
        case 10
                        
            S = load('data\temp.mat');
            eigen_value = S.eigen_value;
            axes(h(4)); %Point to Poincare map - top right
            
            V2_xmin = str2num(get(af2(4),'String'));
            V2_xmax = str2num(get(af2(6),'String'));
            
            %This code only relevant to slope line when drawn by hand:
            %if(line_drawn~=0)
                
            %    set(visible_line,'Visible','off');
            %end
            
            %Code to manually draw the slope by picking two points on the
            %line - another option
            %x=ginput(1);
            %y=ginput(1);
            %slope = (x(2) - y(2))/(x(1) - y(1))
            %line1 = [x(1) y(1)];
            %line2 = [x(2) y(2)];
            %visible_line= line(line1,line2,'Color',[1 0 0]);
            %set(visible_line,'Visible','on');
            %Calculating slope from eigenvalues:
            
            im = (imag(eigen_value));
            re = (real(eigen_value));
            
            tot = (re*2*pi);
            calculated_slope = exp(tot/im);
            line_drawn=1;
            
            %Calculating slope of the linear region using least squares
            %approximation - EXPLAIN!!!!!!!!!!!
            
            %this FOR statement seperates the x values according to the
            %user specs
            index=1;
            for i=1:length(V2vector(:,1))
                if ((V2vector(i,1) >=V2_xmin) && (V2vector(i,1) <=V2_xmax))
                    A(index,:) = [V2vector(i,1) 1];
                    Y(index) = V2vector(i,2);
                    index=index+1;
                end
            end
            
            Y = Y';
            x = (inv(A'*A))*A'*Y; %x is a 2x1 matrix containing the slope and the constant, i.e. mx+c
            xvector=min(V2vector(:,1)):0.05:max(V2vector(:,1));
            y = x(1)*xvector+x(2);
            
            hold on
            plot(xvector,y,'r');
            axis tight;
            hold off
            least_squares_slope = x(1);
            
            
            if(least_squares_slope >= calculated_slope)
                percentage_diff = 100 - (calculated_slope/least_squares_slope)*100;
            else
                percentage_diff = 100 - (least_squares_slope/calculated_slope)*100;
            end
            
            set(af2(10),'String',['Least Squares Approximation: ',num2str(least_squares_slope)]);
            set(af2(11),'String',['Calculated Slope: ',num2str(calculated_slope)]);
            set(af2(12),'String',['% Diff: ',num2str(percentage_diff),'%']);
            set(af2(9:12),'Visible','on');

            axes(h(1));
            rotate3d(h(1));
            
            %Draw the expermental results    
        case 11
            %|||||||||||||||||||||||||||||||||||||||||||||NB|||||||||||||||||||||||||||||||||||||||||||||||||||||
            %||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||||
            %This assumes that the data is in the form colunm1: I3,
            %column2: V2, and column3: V1
            
            % Home Directory
            home = pwd;
            %Open File Open Dialog
            [filename, pathname] = uigetfile('*.mat','Select Data File');   
            if filename ~= 0 % Yes file exists
                % Navigate to data directory - if different from default
                cd(pathname);
                % Load the .mat file and assign it to struct
                % exp_data.(vec)
                exp_data = load(filename);
                names = fieldnames(exp_data);
                vec = char(names);
                
                axes(h(4));
                hold on
                plot(exp_data.(vec)(1:length(exp_data.(vec))-1,2), exp_data.(vec)(2:length(exp_data.(vec)),2),'.r');
                axis tight
                hold off
                cd(home);
            end
            
            axes(h(1));
            rotate3d(h(1));
            
            %Prepare Experimental Results for program
        case 12
            % Display ChuaDiode DrivingPoint Characteristic in popup window
            
            %Background Color
            back_clr = [0.7529 0.7929 0.97];
            % Figure
            f(4) = figure('Color',[0.83137 0.81569 0.78431], ...
                'Color',back_clr, ...
                'Name','Prepare Experimental Results', ...
                'NumberTitle','off', ...
                'Units','normalized','Position',[.3 .1 .42 .84], ...
                'Resize','off', ...
                'ToolBar','none', ...
                'MenuBar','none',...
                'Visible','off', ...
                'WindowStyle','normal');
            
            h(5) = axes('Parent',f(4), ...
                'Color',[1 1 1],'CameraUpVector',[0 1 0],...
                'XColor',[0 0 0],'YColor',[0 0 0],'ZColor',[0 0 0],...
                'Units','normalized','Position',[0.1 0.49 0.8 0.35]);
          
            
            h(6) = axes('Parent',f(4), ...
                'Color',[1 1 1],'CameraUpVector',[0 1 0],...
                'XColor',[0 0 0],'YColor',[0 0 0],'ZColor',[0 0 0],...
                'Units','normalized','Position',[0.1 0.05 0.8 0.35]);
            
            abponc(1) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% Frame  (Load Data1)
                'Units','normalized','Position',[0.025 0.96 0.2 0.035],'Style','frame');
            
            abponc(2) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% Frame  (Load Data2)
                'Units','normalized','Position',[0.025 0.86 0.45 0.1],'Style','frame');
            
            abponc(6) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% Frame  (Remove Noise1)
                'Units','normalized','Position',[0.475 0.96 0.235 0.035],'Style','frame');
            
            abponc(7) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% Frame  (Remove Noise2)
                'Units','normalized','Position',[0.475 0.86 0.52 0.1],'Style','frame');
                        
            abponc(3) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%Filename Text
                'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.0265 0.92 0.15 0.03],'String','Filename:','Style','text');
            
            abponc(4) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%no. of files Text
                'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.0265 0.89 0.15 0.03],'String','No. of files:','Style','text');
           
            abponc(5) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%Load Data Text
                'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.035 0.962 0.15 0.03],'String','Load Data:','Style','text');
             
            abponc(8) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%Remove Noise Text
                'FontName','Helvetica','HorizontalAlignment','left','FontWeight','bold','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.485 0.962 0.22 0.03],'String','Remove Noise:','Style','text');
            
            abponc(9) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%Linear Region: Text
                'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.49 0.92 0.4 0.03],'String','Linear Region:','Style','text');
            
            abponc(10) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%V2(xn)min Text
                'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.49 0.90 0.15 0.03],'String','V2(xn)min','Style','text');
            
            abponc(11) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...%V2(xn)max Text
                'FontName','Helvetica','HorizontalAlignment','left','FontSize',10,'ForegroundColor',[0 0 0],...
                'Units','normalized','Position',[0.75 0.90 0.15 0.03],'String','V2(xn)max','Style','text');
           
            
            ab(1) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% Filename input box
                'HorizontalAlignment','left','Units','normalized','Position',[0.180 0.926 0.15 0.025],'Style','edit','String','data');
            
            ab(2) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% no. of files input box
                'HorizontalAlignment','left','Units','normalized','Position',[0.180 0.895 0.05 0.025],'Style','edit','String','8');
            
            ab(3) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% V2(xn)min input box
                'HorizontalAlignment','left','Units','normalized','Position',[0.63 0.906 0.065 0.025],'Style','edit','String','0.2');
            
            ab(4) = uicontrol('Parent',f(4),'BackgroundColor',[1 1 1],...% V2(xn)max input box
                'HorizontalAlignment','left','Units','normalized','Position',[0.89 0.906 0.065 0.025],'Style','edit','String','0.8');

            % Plot Data - Pushbutton
            H_ABC_pp(511) = uicontrol('Parent',f(4),'Style','pushbutton',...
                'Units','normalized','Position',[0.329 0.864 0.12 0.0325],...
                'ForegroundColor',[0 0 0],...
                'FontWeight','bold',...
                'String','Plot Data',...
                'Visible','on',...
                'CallBack','poincare_popup(13)');
            
            
            % Remove Noise - Pushbutton
            H_ABC_pp(512) = uicontrol('Parent',f(4),'Style','pushbutton',...
                'Units','normalized','Position',[0.5 0.864 0.2 0.0325],...
                'ForegroundColor',[0 0 0],...
                'FontWeight','bold',...
                'String','Remove Noise',...
                'Visible','on',...
                'CallBack','poincare_popup(14)');
            
            
            set(f(4),'Visible','on');
          
            %Plotting data in experimental results popup window    
        case 13
            home = pwd;
            cd results;
            filename = get(ab(1),'String');
            no_of_files= str2num(get(ab(2),'String'));
            
            %Read in files
            for i=1:no_of_files
                name = [filename,num2str(i)];
                name = [name,'.tab'];
                if exist(name)==2 %Check if it exists
                    data{i} = load(name);
                else
                    data{i}=[];
                end
            end
            cd(home);
            for i=no_of_files:10 %10 is the max number of files you can have
                data{i}=[];
            end
            
            results = [data{1};data{2};data{3};data{4};data{5};data{6};data{7};data{8};data{9};data{10}];
            if ~isempty(results)
                %plotting samples
                axes(h(5)) %Point to top graph in experimental results window
                plot(results(:,2),results(:,3),'b.');
            end
            
            
            
            %Removing Noise
        case 14
            home = pwd;
            cd results;
            filename = get(ab(1),'String');
            no_of_files= str2num(get(ab(2),'String'));
            minV2value = str2num(get(ab(3),'String'));
            maxV2value = str2num(get(ab(4),'String'));
            
            %Read in files
            for i=1:no_of_files
                name = [filename,num2str(i)];
                name = [name,'.tab'];
                data{i} = load(name);
            end
            
            for i=no_of_files:10 %10 is the max number of files you can have
                data{i}=[];
            end
            cd(home);
            results = [data{1};data{2};data{3};data{4};data{5};data{6};data{7};data{8};data{9};data{10}];
            results = sortrows(results,2);
            
            vector_index=1;
            section1 = minV2value-0.02;
            section2 = minV2value;
            i = 1;
            noise=0;
            noise_index=1;
            through=0;
           
            % PART 1 OF ALGORITHM - REMOVE NOISE FROM HORIZONTAL SECTION OF GRAPH
            while (i < (length(results(:,1)))) & (section2 < maxV2value)
                tempempty=0;
                passedthrough = 0;
                flag =0;
                index = 1;
                
                section1 = section1+.02;
                section2 = section1+.02;
                
                while (flag==0) & (i < (length(results(:,1)))) & (results(i,2)<section2)  %include - "& results(i,2) < section2"
                    %this should stop
                    %errors ocurring when
                    %there's no data in the
                    %interval selected
                    if (results(i,2) <= section2) & (results(i,2) >= section1)
                        passedthrough = 1;
                        temp(index,:) = [results(i,:)];
                        index=index+1;
                        i=i+1;
                        tempempty=1;
                    else
                        if passedthrough==1
                            flag=1;
                        else
                            i=i+1;
                        end
                    end
                end
                
                if tempempty~=0
                    average = mean(temp(:,3),1);
                    upperlim =average*1.05;
                    lowerlim =average*0.95;
                    for j=1:length(temp(:,1))
                        if (temp(j,3) <= upperlim) & (temp(j,3) >= lowerlim)
                            finalvector(vector_index,:) = temp(j,:);
                            vector_index=vector_index+1;
                            through = through+1;
                        else
                            noisevector(noise_index,:) = temp(j,:);
                            noise=noise+1;
                            noise_index = noise_index+1;
                        end
                    end
                    clear temp;
                    
                else
                    i=i+1;
                end
                
            end
            
            
            first_through = through;
            first_noise = noise;
 
            %PART ii of algorithm - remove noise from vertical part of graph
            ceilingvalue = max(results(:,2));
            while( i < length(results(:,1)))
                if(results(i,2) < 1.02) & (results(i,2) >maxV2value)
                    finalvector(vector_index,:) = results(i,:);
                    vector_index=vector_index+1;
                    through = through+1;
                    i=i+1;
                else
                    noisevector(noise_index,:)= results(i,:);
                    noise=noise+1;
                    noise_index=noise_index+1;
                    i=i+1;
                end
            end       
            
            second_through = through;
            second_noise = noise;
            
            axes(h(6)) %Point to bottom graph in experimental results window
            plot(finalvector(:,2),finalvector(:,3),'b.');
           
            
            finalvector = [finalvector(:,2),finalvector(:,3)];
           
            savefile = [filename,'_finalvector','.mat'];
            save(savefile,'finalvector'); %Save the finalvector as NAME_finalvector.mat file
            
            
       
    %Export each figure to another window for individual saving    
    case 15
        figure(11)
        
        clear('poincare_points');
        S = load('data\temp.mat');
        Timeseries = S.Timeseries;
        sys_variables = S.sys_variables;
        plane = [str2num(get(a(15),'String')),str2num(get(a(16),'String')),str2num(get(a(17),'String')),str2num(get(a(18),'String'))];
        integ_variables = S.integ_variables;
        eigen_value = S.eigen_value;
        %set(a(22),'String',[get(a(15),'String'),'x+',get(a(16),'String'),'y+',get(a(17),'String'),'z+',get(a(18),'String'),'=0']);
        plane = [str2num(get(a(15),'String')),str2num(get(a(16),'String')),str2num(get(a(17),'String')),str2num(get(a(18),'String'))];
        
        integ_variables(6) = str2num(get(a(14),'String')); %Gets the step_size_divisor variable from user
        
        remove_ramp_up = str2num(get(a(24),'String'));
        Timeseries=Timeseries(:,remove_ramp_up:length(Timeseries));

        poincare_points_a = poincare(Timeseries,plane,sys_variables,integ_variables);

        poincare_points = poincare_points_a(:,1:4);
        approx_1 = poincare_points_a(:,5:8);
        
        if isreal(poincare_points)
            for i=1:length(poincare_points(:,1)) %plots the points in the window, red for one direction, blue for the other
                if poincare_points(i,4) == 1
                    plot3(poincare_points(i,1),poincare_points(i,2),poincare_points(i,3),'or');
                    hold on;
                else
                    plot3(poincare_points(i,1),poincare_points(i,2),poincare_points(i,3),'ob');
                    hold on;
                end
            end  
            hold off;   
            
            xmax = max(poincare_points(:,1));
            xmin = min(poincare_points(:,1));
            ymin = min(poincare_points(:,2));
            ymax = max(poincare_points(:,2));
            zmax = max(poincare_points(:,3));
            zmin = min(poincare_points(:,3));
            if xmax > 0
                xmax=xmax*1.1;
            else
                xmax=xmax*.9;
            end
            if xmin > 0
                xmin=xmin*.9;
            else
                xmin=xmin*1.1;
            end
            if zmax > 0
                zmax=zmax*1.1;
            else
                zmax=zmax*.9;
            end
            if zmin > 0
                zmin=zmin*.9;
            else
                zmin=zmin*1.1;
            end
            if ymin > 0
                ymin=ymin*.9;
            else
                ymin=ymin*1.1;
            end
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),xmin,xmax,zmin,zmax,ymin,ymax);
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',.5,'EdgeAlpha',0.1);
           % set(gcf,'CameraUpVector',[plane(1) plane(2) plane(3)]);
            xlabel('I3');
            ylabel('V2');
            zlabel('V1');
            
            if get(H_ABC_pp(506),'Value') == 1
                view(0,90);
                axis tight;
            elseif get(H_ABC_pp(507),'Value') == 1
                view(0,0);  
                axis tight;
            elseif get(H_ABC_pp(508),'Value') == 1
                view(90,0);
                axis tight;
            end
            
            figure(12); % Point to Poincare Map - top right 
            
            %This separates the poincare map points into two bins - one for each 
            %direction that they cross the plane in
            
            poin = poincare_points;
            
            bindex=1;
            aindex=1;
            for i=1:length(poin(:,1))           
                if poin(i,4)==0  %this finds what direction they crossed in 
                    a_poin(aindex,:) = poin(i,:);
                    aindex = aindex+1;
                else
                    b_poin(bindex,:) = poin(i,:);
                    bindex = bindex+1;
                end
                
            end
            
            if(exist('a_poin')) %if a_poin exists
                %poin_xnplusone is the value at time t+1
                %this will be plotted against poin_xn, the value at time t
                
                if get(H_ABC_pp(503),'Value') == 0
                    poin_xnplusone(1:length(a_poin(:,1))-1,:) = a_poin(2:length(a_poin(:,1)),:);
                    poin_xn(1:length(a_poin(:,1))-1,:) = a_poin(1:length(a_poin(:,1))-1,:);
                    
                    I3vector = [poin_xn(:,1),poin_xnplusone(:,1)]; %Sorts them into vectors that can be plotted
                    V2vector = [poin_xn(:,2),poin_xnplusone(:,2)];
                    V1vector = [poin_xn(:,3),poin_xnplusone(:,3)];
                    
                    plot(V2vector(:,1),V2vector(:,2),'.') %Plots the V2 vector at xn+1 against xn
                    
                    
                else
                    poin_xnplusone(1:length(b_poin(:,1))-1,:) = b_poin(2:length(b_poin(:,1)),:);
                    poin_xn(1:length(b_poin(:,1))-1,:) = b_poin(1:length(b_poin(:,1))-1,:);
                    
                    I3vector = [poin_xn(:,1),poin_xnplusone(:,1)]; %Sorts them into vectors that can be plotted
                    V2vector = [poin_xn(:,2),poin_xnplusone(:,2)];
                    V1vector = [poin_xn(:,3),poin_xnplusone(:,3)];
                    
                    plot(V2vector(:,1),V2vector(:,2),'.') %Plots the V2 vector at xn+1 against xn
                    
                end
            end
            xlabel('V2(xn)');
            ylabel('V2(xn+1)');
           
            
            clear V2vector;
            clear poin_xnplusone;
            clear poin_xn;
            clear poin;
            
            
            figure(13); % Point to main stage            
            plot3(Timeseries(1,:),Timeseries(2,:),Timeseries(3,:),'Color','k');
            xlabel('I3');
            ylabel('V2');
            zlabel('V1');
            
            xmin=min(Timeseries(1,:));
            xmax=max(Timeseries(1,:));
            zmin=min(Timeseries(3,:));
            zmax=max(Timeseries(3,:));
            ymin=min(Timeseries(2,:));
            ymax=max(Timeseries(2,:));
            
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),(xmin+xmin/10),(xmax+xmax/10),(zmin+zmin/10),(zmax+zmax/10),(ymin+ymin/10),(ymax+ymax/10));
            %            vertex_matrix = [xmin 0 zmin; xmin 0 zmax; xmax 0 zmin; xmax 0 zmax]

            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',.8,'EdgeAlpha',0.1);
            %view(view_angle);%axis off;
            rotate3d;
            
                    else
            % Read in warning image and place it Small Axes

            warning off all % Suppress Warning
            figure(12);
            imagen=imread('etc\warning.png');
            image(imagen);
            axis off;
            warning on all % Turn back on warnings
            
            
            
            figure(11);
            warning off all % Suppress Warning
            imagen=imread('etc\warning.png');
            image(imagen);
            axis off;
            warning on all % Turn back on warnings
            figure(10); % Point to main stage            
            plot3(Timeseries(1,:),Timeseries(2,:),Timeseries(3,:),'Color','k');
            xlabel('Current');
            xmin=min(Timeseries(1,:));
            xmax=max(Timeseries(1,:));
            zmin=min(Timeseries(3,:));
            zmax=max(Timeseries(3,:));
            ymin=min(Timeseries(2,:));
            ymax=max(Timeseries(2,:));
            
            vertex_matrix = plotplane(plane(1),plane(2),plane(3),plane(4),(xmin+xmin/10),(xmax+xmax/10),(zmin+zmin/10),(zmax+zmax/10),(ymin+ymin/10),(ymax+ymax/10));
            %vertex_matrix = [xmin 0 zmin; xmin 0 zmax; xmax 0 zmin; xmax 0 zmax]
            patch('Vertices',vertex_matrix,'Faces',[1 2 4 3],'EdgeColor','r','FaceColor',[ 0.000 0.502 0.000 ],'FaceAlpha',1,'EdgeAlpha',0.1);
            %view(view_angle);%axis off;
            axis on;
            rotate3d(h(1));
            

        end
            
        end
        
        
        
        
        
    
    
function regional = intersection_finder(time_series,plane)

% 'regional' each row of this returned matrix details the intersection of the trajectory 'time_Series' 
% with the input 'plane'
% The first three columns (n,1:3) of each row are the xyz co-ordinates of the point one side of the plane.
% The next three columns (n,4:6) detail the first point in the data set after the intersection with the plane.
% The seventh column (n,10)is updated every time an intersection is found and details the total number of 
% intersections in that data set.

% Consider two point vectors P1(x1,y1,z1) and P2(x2,y2,z2) taken from the time_series data set
% Now consider a plane of the form (aX + bY + cZ - d = 0)
% If the vector P1P2 intersects the plane then one of the following
% conditions will be true
% condition1:
%  (a*x1 + b*y1 + c*z1 -d > 0) and (a*x2 + b*y2 + c*z2 -d < 0)
% condition2:
%  (a*x1 + b*y1 + c*z1 -d > 0) and (a*x2 + b*y2 + c*z2 -d < 0)
% To complete the algorithm the condition where one of the points lies upon the
% plane is accounted for
% condition3:
% (a*x1 + b*y1 + c*z1 -d == 0) or (a*x2 + b*y2 + c*z2 -d == 0)

a=plane(1);
b=plane(2);
c=plane(3);
d=plane(4);

check='voido';
time_series=time_series';
index=1;
for n=2:length(time_series)
    
    cond11=(a*time_series(n-1,1) + b*time_series(n-1,2) + c*time_series(n-1,3) +d > 0);
    cond12=(a*time_series(n,1) + b*time_series(n,2) + c*time_series(n,3) +d < 0);
    cond1=cond11 && cond12;
    % Checks on direction of the vector for intersection
    
    cond21=(a*time_series(n-1,1) + b*time_series(n-1,2) + c*time_series(n-1,3) +d < 0);
    cond22=(a*time_series(n,1) + b*time_series(n,2) + c*time_series(n,3) +d > 0);
    cond2=cond21 && cond22;
    % Checks the opposite direction of the vector for intersection
    
    cond3=(a*time_series(n,1) + b*time_series(n,2) + c*time_series(n,3) +d == 0);
    % Checks if either point lies upon the plane
    
    if(cond1)
        regional(index,1:3)=time_series(n-1,1:3);
        regional(index,4:6)=time_series(n,1:3);
        regional(index,7)=1;
        regional(index,8)=0;
        index=index+1;
        check='valid';
    end;
    
    if(cond2)
        regional(index,1:3)=time_series(n-1,1:3);
        regional(index,4:6)=time_series(n,1:3); 
        regional(index,7)=0;
        regional(index,8)=0;
        index=index+1;
        check='valid';
    end;    
    
    if(cond3)
        
        regional(index,1:3)=time_series(n,1:3);
        regional(index,4:6)=time_series(n,1:3);
        
        cond23=(a*time_series(n+1,1) + b*time_series(n+1,2) + c*time_series(n+1,3) +d > 0);
        cond25=cond21&&cond23;       
        cond13=(a*time_series(n+1,1) + b*time_series(n+1,2) + c*time_series(n+1,3) +d < 0);
        cond15=cond11 && cond13;
        
        if(cond15)
            regional(index,7)=1;
        end;
        if(cond25)
            regional(index,7)=0;
        end;
        regional(index,8) = 1; %If the point lies on the plane, column 8 indicates this with a one
        index=index+1;
        check='valid';
    end;
   
end;


if(check=='voido')
    clear regional;
    regional(1,1:8)=i; %Sets regional to all imaginary numbers - this occurs when the plane does not
                       %intersect the trajectory

end

%%%%%%%%%%%%%%%%%%%&

%%%%%%%%%%%%%%%%%%%%

function poin = poincare(time_series_1,plane,sys_variables,integ_variables)

% 'time_series' is Nx1 matrix where N is the length of the timeseries
% 'plane' is the plane equation coefficients in matrix form
%  if a plane is desccribed by ax+by+cz-d=0 then the matrix is of the form [a b c d]
%  sys_variables = [L,R0,C2,G,Ga,Gb,C1,E];
%  integ_variables = [x0,y0,z0,dataset_size,step_size];

a=plane(1,1);
b=plane(1,2);
c=plane(1,3);
d=plane(1,4);

execute_flag=1;


approx_1=intersection_finder(time_series_1,plane);

if(~isreal(approx_1)) %if approx_1 is not real, the plane doesnt intersect the trajectory
    poin = approx_1;
    
elseif (length(approx_1(:,1)) ~= 0) 
    
    %this will execute for every set of intersection points found except
    %for the ones which are found to have been lying on the plane already
    %i.e. those with a 1 in the 8th column of approx_1
    %It uses the first point as an initial condition, and then integrates
    %from this initial condition until it crosses the plane
    for i=1:length(approx_1(:,1))
        if(approx_1(i,8) ~= 1)
            approx_2(i,:) = chua(approx_1(i,:),sys_variables,integ_variables,plane);
            
            
        elseif(approx_1(i,8)==1)
            approx_2(i,:) = [approx_1(i,1:3), approx_1(i,8)];
            
        end
        
        
        i=i+1;
    end
        
    poin=[approx_2,approx_1(:,1:4)];
    
end

        
        





function point = chua(approx_1,sys_variables,integ_variables,plane)
%This will integrate from the initial conditions i.e. the point on one side
%of the plane, to the plane until it either lands on the plane or crosses
%to the other side of it
%it returns an n by 4 matrix "point" which contains the x y z co-ordinates
%for the point and also what direction it was going in indicated by a 1 or
%a 0

% Syntax: TimeSeries=chua(sys_variables,integ_variables)
% sys_variables = [L,R0,C2,G,Ga,Gb,E,C1];
% integ_variables = [x0,y0,z0,dataset_size,step_size];


% Models Initial Variables
%-------------------------
%NB C2 and C1 were passed incorrectly to this function(as usual) so be
%careful of this when editing - peadar
L  =           sys_variables(1);
R0 =           sys_variables(2);
C1 =           sys_variables(3);
G  =           sys_variables(4);
Ga =           sys_variables(5);
Gb =           sys_variables(6);
E  =           sys_variables(7);
C2 =           sys_variables(8);
x0 =           integ_variables(1);
y0 =           integ_variables(2);
z0 =           integ_variables(3);
dataset_size = integ_variables(4);
step_size =    integ_variables(5);
step_size_divisor = integ_variables(6);

a=plane(1,1);
b=plane(1,2);
c=plane(1,3);
d=plane(1,4);


x0 = approx_1(1);
y0 = approx_1(2);
z0 = approx_1(3);

TimeSeries = [x0, y0, z0]'; % these are the initial conditions on one side of the plane
                            

% Optimized Runge-Kutta Variables
%--------------------------------

h = step_size*G/C2;

h = h/step_size_divisor;
h2 = (h)*(.5);
h6 = (h)/(6);

anor = Ga/G;
bnor = Gb/G;
bnorplus1 = bnor + 1;
alpha = C2/C1;
beta = C2/(L*G*G);
gammaloc = (R0*C2)/(L*G);

bh = beta*h;
bh2 = beta*h2;
ch = gammaloc*h;
ch2 = gammaloc*h2;
omch2 = 1 - ch2;

k1 = [0 0 0]';
k2 = [0 0 0]';
k3 = [0 0 0]';
k4 = [0 0 0]';
M = [0 0 0]';

% Calculate Time Series
%----------------------

M(1) = TimeSeries(3)/E;
M(2) = TimeSeries(2)/E;
M(3) = TimeSeries(1)/(E*G);

check = 'voido';
i=1;
index = 1;
while check == 'voido'
    % Runge Kutta
    % Round One
    k1(1) = alpha*(M(2) - bnorplus1*M(1) - (.5)*(anor - bnor)*(abs(M(1) + 1) - abs(M(1) - 1)));
    k1(2) = M(1) - M(2) + M(3);
    k1(3) = -beta*M(2) - gammaloc*M(3);
    % Round Two
    temp = M(1) + h2*k1(1);
    k2(1) = alpha*(M(2) + h2*k1(2) - bnorplus1*temp - (.5)*(anor - bnor)*(abs(temp + 1) - abs(temp - 1)));
    k2(2) = k1(2) + h2*(k1(1) - k1(2) + k1(3));
    k2(3) = omch2*k1(3) - bh2*k1(2);
    % Round Three
    temp = M(1) + h2*k2(1);
    k3(1) = alpha*(M(2) + h2*k2(2) - bnorplus1*temp - (.5)*(anor - bnor)*(abs(temp + 1) - abs(temp - 1)));
    k3(2) = k1(2) + h2*(k2(1) - k2(2) + k2(3));
    k3(3) = k1(3) - bh2*k2(2) - ch2*k2(3);
    % Round Four
    temp = M(1) + h*k3(1);
    k4(1) = alpha*(M(2) + h*k3(2) - bnorplus1*temp - (.5)*(anor - bnor)*(abs(temp + 1) - abs(temp - 1)));
    k4(2) = k1(2) + h*(k3(1) - k3(2) + k3(3));
    k4(3) = k1(3) - bh*k3(2) - ch*k3(3);
    
    M = M + (k1 + 2*k2 + 2*k3 + k4)*(h6);
    
    TimeSeries(3,i+1) = E*M(1);
    TimeSeries(2,i+1) = E*M(2); 
    TimeSeries(1,i+1) = (E*G)*M(3);
  

    %Code for checking whether or not the point has passed the plane

    cond11=(a*TimeSeries(1,i) + b*TimeSeries(2,i) + c*TimeSeries(3,i) +d > 0);
    cond12=(a*TimeSeries(1,i+1) + b*TimeSeries(2,i+1) + c*TimeSeries(3,i+1) +d < 0);
    cond1=cond11 && cond12;
    % Checks on direction of the vector for intersection
    
    cond21=(a*TimeSeries(1,i) + b*TimeSeries(2,i) + c*TimeSeries(3,i) +d < 0);
    cond22=(a*TimeSeries(1,i+1) + b*TimeSeries(2,i+1) + c*TimeSeries(3,i+1) +d > 0);
    cond2=cond21 && cond22;
    % Checks the opposite direction of the vector for intersection
    
    cond3=(a*TimeSeries(1,i+1) + b*TimeSeries(2,i+1) + c*TimeSeries(3,i+1) +d == 0);
    % Checks if either point lies upon the plane
    
    if(cond1)
        point(1,1:3)=TimeSeries(1:3,i)';
        point(1,4)=1;
        check='valid';
    end;
    
    if(cond2)
        point(1,1:3)=TimeSeries(1:3,i)'; 
        point(1,4)=0;
        check='valid';
    end;
    
    if(cond3)
        
        point(1,1:3)=TimeSeries(1:3,i)';
        point(1,4) = approx(4); % Sign in this case is the same as the last time
        check='valid';
    end;

    i=i+1;
 
end



function vertex_matrix = plotplane(a,b,c,d,xmin,xmax,zmin,zmax,ymin,ymax)
% Uses simple geometry:'plane equation' manipulation to return a vertex
% matrix using plane equation ax+by+cz+d=0 coefficients to calculate
% unknown vertices - ymin ymax



if b~=0
vertex_matrix = [xmin ((-d-a*xmin-c*zmin)/b) zmin; xmin ((-d-a*xmin-c*zmax)/b) zmax; xmax ((-d-a*xmax-c*zmin)/b) zmin; xmax ((-d-a*xmax-c*zmax)/b) zmax];
end

if (b==0 & a~=0)
vertex_matrix = [((-d-b*ymin-c*zmin)/a) ymin zmin;((-d-b*ymin-c*zmax)/a) ymin zmax;((-d-b*ymax-c*zmin)/a) ymax zmin;((-d-b*ymax-c*zmax)/a) ymax zmax];
end

if a==0
vertex_matrix = [xmin ymin ((-d-b*ymin-a*xmin)/c);xmax ymin ((-d-b*ymin-a*xmax)/c);xmin ymax ((-d-b*ymax-a*xmin)/c);xmax ymax ((-d-b*ymax-a*xmax)/c)];
end

Contact us