function varargout = DELAY_ODE_GUI(varargin)
% DELAY_ODE_GUI M-file for DELAY_ODE_GUI.fig
%      DELAY_ODE_GUI, by itself, creates a new DELAY_ODE_GUI or raises the existing
%      singleton*.
%
%      H = DELAY_ODE_GUI returns the handle to a new DELAY_ODE_GUI or the handle to
%      the existing singleton*.
%
%      DELAY_ODE_GUI('CALLBACK',hObject,eventData,handles,...) calls the local
%      function named CALLBACK in DELAY_ODE_GUI.M with the given input arguments.
%
%      DELAY_ODE_GUI('Property','Value',...) creates a new DELAY_ODE_GUI or raises the
%      existing singleton*.  Starting from the left, property value pairs are
%      applied to the GUI before DELAY_ODE_GUI_OpeningFcn gets called.  An
%      unrecognized property name or invalid value makes property application
%      stop.  All inputs are passed to DELAY_ODE_GUI_OpeningFcn via varargin.
%
%      *See GUI Options on GUIDE's Tools menu.  Choose "GUI allows only one
%      instance to run (singleton)".
%
% See also: GUIDE, GUIDATA, GUIHANDLES

% Edit the above text to modify the response to help DELAY_ODE_GUI

% Last Modified by GUIDE v2.5 21-Apr-2011 14:48:59

% Begin initialization code - DO NOT EDIT
gui_Singleton = 0;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @DELAY_ODE_GUI_OpeningFcn, ...
                   'gui_OutputFcn',  @DELAY_ODE_GUI_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin && ischar(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT


% --- Executes just before DELAY_ODE_GUI is made visible.
function DELAY_ODE_GUI_OpeningFcn(hObject, eventdata, handles, varargin)
% This function has no output args, see OutputFcn.
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
% varargin   command line arguments to DELAY_ODE_GUI (see VARARGIN)

% Choose default command line output for DELAY_ODE_GUI
handles.output = hObject;
 set(handles.no_MRS,'Value',1,'Callback',{@FIX_MRS_Callback, handles,1})
  set(handles.yes_MRS,'Value',0,'Callback',{@FIX_MRS_Callback, handles,2})
 set(handles.no_rotated_frame,'Value',1,'Callback',{@ROTATE_FRAME_Callback,handles,1})
  set(handles.yes_rotated_frame,'Value',0,'Callback',{@ROTATE_FRAME_Callback,handles,2})
 handles.fix_mrs = 0;
 handles.rotate_frame = 0;
% Update handles structure
guidata(hObject, handles);



% --- Outputs from this function are returned to the command line.
function varargout = DELAY_ODE_GUI_OutputFcn(hObject, eventdata, handles) 
% varargout  cell array for returning output args (see VARARGOUT);
% hObject    handle to figure
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Get default command line output from handles structure
varargout{1} = handles.output;
% UIWAIT makes DELAY_ODE_GUI wait for user response (see UIRESUME)
% uiwait(handles.figure1);



% function FIX_MRS_Callback(hObject,eventdata, handles,radio_value)
    function FIX_MRS_Callback( hObject,eventdata,handles,radio_value)
% hObject    handle to FIX_MRS (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% initialize plot without mrs fixed

 
switch radio_value
    
    case 1 % no, don't fix mrs
        set(handles.yes_MRS,'Value',0)
        set(handles.no_MRS,'Value',1)
        
%         if length(handles.ymax)==1
%         set(handles.ymax,max(get(gca,'YLim')))
%         else
%         set(handles.ymax(1:2),get(gca,'XLim'))
%         set(handles.ymax(3:4),get(gca,'YLim'))
%         end
        
        handles.fix_mrs = 0; 
        handles
        guidata(hObject, handles);
        
        
    case 2 % yes, fix mrs
        set(handles.yes_MRS,'Value',1)
        set(handles.no_MRS,'Value',0)
        
%         if length(handles.ymax)==1
%         set(handles.ymax,max(get(gca,'YLim')))
%         else
%         set(handles.ymax(1:2),get(gca,'XLim'))
%         set(handles.ymax(3:4),get(gca,'YLim'))
%         end
        
        handles.fix_mrs = 1;
        handles
        guidata(hObject, handles);

end

% --- Executes on button press in ROTATED_FRAME_NO.
function ROTATE_FRAME_Callback(hObject,eventdata,handles,radio_value)
% hObject    handle to ROTATED_FRAME_NO (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hint: get(hObject,'Value') returns toggle state of ROTATED_FRAME_NO
switch radio_value
    
    case 1 % no, don't rotated_frame
        set(handles.yes_rotated_frame,'Value',0)
        set(handles.no_rotated_frame,'Value',1)
        
%         if length(handles.ymax)==1
%         set(handles.ymax,max(get(gca,'YLim')))
%         else
%         set(handles.ymax(1:2),get(gca,'XLim'))
%         set(handles.ymax(3:4),get(gca,'YLim'))
%         end

        handles.rotate_frame = 0;
         handles
        guidata(hObject, handles); 
        
    case 2 % yes, rotated_frame
        set(handles.yes_rotated_frame,'Value',1)
        set(handles.no_rotated_frame,'Value',0)
        
%         if length(handles.ymax)==1
%         set(handles.ymax,max(get(gca,'YLim')))
%         else
%         set(handles.ymax(1:2),get(gca,'XLim'))
%         set(handles.ymax(3:4),get(gca,'YLim'))
%         end
%         
        handles.rotate_frame = 1;
         handles
        guidata(hObject, handles);     
end




function AMIN_Callback(hObject, eventdata, handles)
% hObject    handle to AMIN (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of AMIN as text
%        str2double(get(hObject,'String')) returns contents of AMIN as a double
handles.Amin = str2double(get(hObject,'String'));
% Update handles structure
guidata(hObject, handles);


function AMAX_Callback(hObject, eventdata, handles)
% hObject    handle to AMAX (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of AMAX as text
%        str2double(get(hObject,'String')) returns contents of AMAX as a double
handles.Amax = str2double(get(hObject,'String'));
% Update handles structure
guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function AMAX_CreateFcn(hObject, eventdata, handles)
% hObject    handle to AMAX (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end



function RMIN_Callback(hObject, eventdata, handles)
% hObject    handle to RMIN (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of RMIN as text
%        str2double(get(hObject,'String')) returns contents of RMIN as a double
handles.Rmin = str2double(get(hObject,'String'));
% Update handles structure
guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function RMIN_CreateFcn(hObject, eventdata, handles)
% hObject    handle to RMIN (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end



function RMAX_Callback(hObject, eventdata, handles)
% hObject    handle to RMAX (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of RMAX as text
%        str2double(get(hObject,'String')) returns contents of RMAX as a double
handles.Rmax = str2double(get(hObject,'String'));
% Update handles structure
guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function RMAX_CreateFcn(hObject, eventdata, handles)
% hObject    handle to RMAX (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end


% --- Executes on button press in pushbutton1.
function pushbutton1_Callback(hObject, eventdata, handles)
% hObject    handle to pushbutton1 (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)
handles.ymax
if(length(handles.ymax)==1)
dynamicBif([handles.Amin, handles.Amax],[handles.Rmin, handles.Rmax],...
    [-handles.ymax handles.ymax -handles.ymax handles.ymax],...
handles.curves,handles.fix_mrs,handles.rotate_frame);
else
   dynamicBif([handles.Amin, handles.Amax],[handles.Rmin, handles.Rmax],...
    handles.ymax,handles.curves,handles.fix_mrs,handles.rotate_frame);
end 
return


function AXIS_WINDOW_EDITABLE_TEXT_Callback(hObject, eventdata, handles)
% hObject    handle to AXIS_WINDOW_EDITABLE_TEXT (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of AXIS_WINDOW_EDITABLE_TEXT as text
%        str2double(get(hObject,'String')) returns contents of AXIS_WINDOW_EDITABLE_TEXT as a double
handles.ymax = eval(get(hObject,'String'))
% Update handles structure
guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function AXIS_WINDOW_EDITABLE_TEXT_CreateFcn(hObject, eventdata, handles)
% hObject    handle to AXIS_WINDOW_EDITABLE_TEXT (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end


% % --- Executes on button press in yes_MRS.
% function yes_MRS_Callback(hObject, eventdata, handles)
% % hObject    handle to yes_MRS (see GCBO)
% % eventdata  reserved - to be defined in a future version of MATLAB
% % handles    structure with handles and user data (see GUIDATA)
% 
% % Hint: get(hObject,'Value') returns toggle state of yes_MRS
% 
% 
% % --- Executes on button press in no_MRS.
% function no_MRS_Callback(hObject, eventdata, handles)
% % hObject    handle to no_MRS (see GCBO)
% % eventdata  reserved - to be defined in a future version of MATLAB
% % handles    structure with handles and user data (see GUIDATA)
% 
% % Hint: get(hObject,'Value') returns toggle state of no_MRS
% 
% 
% % --- Executes on button press in yes_rotated_frame.
% function yes_rotated_frame_Callback(hObject, eventdata, handles)
% % hObject    handle to yes_rotated_frame (see GCBO)
% % eventdata  reserved - to be defined in a future version of MATLAB
% % handles    structure with handles and user data (see GUIDATA)
% 
% % Hint: get(hObject,'Value') returns toggle state of yes_rotated_frame
% 
% 
% % --- Executes on button press in no_rotated_frame.
% function no_rotated_frame_Callback(hObject, eventdata, handles)
% % hObject    handle to no_rotated_frame (see GCBO)
% % eventdata  reserved - to be defined in a future version of MATLAB
% % handles    structure with handles and user data (see GUIDATA)
% 
% % Hint: get(hObject,'Value') returns toggle state of no_rotated_frame



function NUMB_CURVES_EDITABLE_TEXT_Callback(hObject, eventdata, handles)
% hObject    handle to NUMB_CURVES_EDITABLE_TEXT (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    structure with handles and user data (see GUIDATA)

% Hints: get(hObject,'String') returns contents of NUMB_CURVES_EDITABLE_TEXT as text
%        str2double(get(hObject,'String')) returns contents of NUMB_CURVES_EDITABLE_TEXT as a double
handles.curves = str2double(get(hObject,'String'));
% Update handles structure
guidata(hObject, handles);

% --- Executes during object creation, after setting all properties.
function NUMB_CURVES_EDITABLE_TEXT_CreateFcn(hObject, eventdata, handles)
% hObject    handle to NUMB_CURVES_EDITABLE_TEXT (see GCBO)
% eventdata  reserved - to be defined in a future version of MATLAB
% handles    empty - handles not created until after all CreateFcns called

% Hint: edit controls usually have a white background on Windows.
%       See ISPC and COMPUTER.
if ispc && isequal(get(hObject,'BackgroundColor'), get(0,'defaultUicontrolBackgroundColor'))
    set(hObject,'BackgroundColor','white');
end


function  dynamicBif(A,R,window,j,fix_mrs,rotate_frame)
%
% dynamicBif.m produces a dynamic plot of the first j bifurcation 
% curves of the two delay system in the BC plane, while A varies.
% It is also possible to vary R. Or if one wants a static plot,
% make A and R scalars.
% 
%
% INPUT VARIABLES
%--------------------------------------------------
% A - can be a single value or an interval [Ao Af]
% R - can be a single value or an interval [Ro Rf]
% j - the number of bifurcation curves
%
% window  - a positive scalar or a 1 by 4 vector that sets axis window
%           if window is a vector, then that is passed to the
%           'axis' command, i.e. set axis(window)
% nframes - determines the interval spacing for A or R
%
% speed_control - frames per second (fps) a floating point 
%                 number which controls the speed of the movie
%                 (not currently implemented)  
%-------------------------------------------------
% COLOR SCHEME - the color of the jth curve is the same as the color
%                   of the mod(j,12)'th curve, eg curve 13 has the same 
%                   color as curve 1, curve 15 has the same color of 
%                   curve 3, etc. 
%
% CURVE  COLOR
% ----------------
% A+B+C=0 violet 
%  1    blue   same color as the paper
%  2    green   same color as the paper
%  3    black   
%  4    red      same color as the paper
%  5    yellow
%  6    orange    same color as the paper
%  7    cyan
%  8    blue (changed to dark green 4/6/11  
%  9    grey
% 10    green   (changed to purple 4/6/11
% 11    magenta
% 12    red   (changed to grey blue 4/6/11 
%---------------------------------------------------
%
% call UISETCOLOR in the command window to determine
%                  an rbg triplet for a selected color
%
nframes = 10;

% Define A = [Ao Af], or as a single value
if(isequal(A(1),A(2)) && ~isequal(R(1),R(2)) )
    %then A is fixed and R is bifurcating  
    a = A(1)*ones(1,nframes);
    r = linspace(R(1),R(2),nframes);
elseif(isequal(R(1),R(2)) && ~isequal(A(1),A(2)) )
    %then R is fixed and A is bifurcating  
    a = linspace(A(1),A(2),nframes);
    r = R(1)*ones(1,nframes);
elseif(isequal(A(1),A(2)) && isequal(R(1),R(2)) )
   %then both A and R are fixed-1 still frame
   nframes=1;
   a = A(1)*ones(1,nframes);
   r = R(1)*ones(1,nframes);
elseif(~isequal(A(1),A(2)) && ~isequal(R(1),R(2)) )
%then both A and R are bifurcating
    a = linspace(A(1),A(2),nframes);
    r = linspace(R(1),R(2),nframes);
end
 

% the following is the rotated coordsys (B,C)->(x,y) with area preservation
B  = @(w,A,R)(  A*sin(w*R)+w.*cos(w*R))./sin(w*(1-R));
C  = @(w,A,R)(-(A*sin(w)  +w.*cos(w)  )./sin(w*(1-R)));

XX  = @(w,A,R)(  (B(w,A,R)+C(w,A,R))/sqrt(2)  );
YY  = @(w,A,R)(  (B(w,A,R)-C(w,A,R))/sqrt(2)  );


Btrans =  @(p,R)((-1)^p*( (1-R).*cos(p*pi*R./(1-R)) -p*pi*R.*csc(p*R*pi./(1-R)) )./ ((1-R).^2));
Ctrans =  @(p,R)((-1)^(p+1)*( (1-R).*cos(p*pi./(1-R)) -p*pi*csc(p*pi./(1-R)) ) ./ ((1-R).^2));
CC      = @(Bt,B,Ct,p)((Bt-B+(-1)^p*Ct)/(-1)^p);


Frames = moviein(nframes);
% UNCOMMENT THE NEXT TWO LINES IF YOU WANT TO OPEN THE
%  FIGURE WINDOW WITH THE FULL SCREEN.
%
if nargin == 4
 hand = figure(1); 
 set(hand,'position',[1           1        1024         648]); %get(0,'screensize')
  set(gca,'FontSize',18,'FontWeight','bold')
end



z=window(2);
%plot a bunch of frames
for q = 1:nframes
    
    if(r(q)<=0.5)
     h = get_trans_points(r(q));
    else
        h=0;
    end
     
     % COLOR SCHEME FAMILY SECTION
% ------------------------------
tol=.005;

    if(abs(r(q)-1/2)<tol || abs(r(q)-2/3)<tol || abs(r(q)-3/4)<tol)
        % R=1/2 color scheme-blue,green,black,red
        color_code = [0 0 1; 0 1 0];
        mod_number = 2; 
    elseif(abs(r(q)-0.4)<tol)
        % R=.4 
        color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0; 1 1 0; 1 .50196 0;];
        mod_number = 6; 
    elseif(abs(r(q)-1/3)<tol)
        % R=1/3 color scheme-blue,green,black,red
        color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0];
        mod_number = 4;               
    elseif(abs(r(q)-1/4)<tol)
        % R=0.25 color scheme-blue,green,black,red
        color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0; 0.584313725490196   0.388235294117647   0.3882352941; 1 .50196 0;];
        mod_number = 6;
    elseif(abs(r(q)-1/5)<tol)
        % R=0.2
        color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0; 1 1 0; 1 .50196 0;0 1 1;0.2314    0.4431    0.3373;];
        mod_number = 8; 
    elseif(abs(r(q)-1/6)<tol)
        % R=0.2
        color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0; 1 1 0; 1 .50196 0;...
            0 1 1;0.2314    0.4431    0.3373;1/2 1/2 1/2;  0.6000    0.2000    1.0000;];
        mod_number = 10; 
    else %old school color scheme
        mod_number=12;
        color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0; 1 1 0; 1 .50196 0;...
            0 1 1;0.2314    0.4431    0.3373;1/2 1/2 1/2; 0.6000    0.2000    1.0000;...
                       1 0 1; 0.8000    0.8000    1.0000];
    end
% %     old school colors
%  mod_number=12;
%     color_code = [0 0 1; 0 1 0; 0 0 0; 1 0 0;...
%                        1 1 0; 1 .50196 0; 0 1 1; ...
%                        0 0 1; 1/2 1/2 1/2; 0 1 0;...
%                        1 0 1; 1 0 0];
%     % end color scheme section. 
%------------------------------------------------

% compute coordinates and plot curves

    for i = 1:j  %number of bif curves determined by input var j
        w = linspace((i-1)*3.1416/(1-r(q)),i*3.14158/(1-r(q)),1100);
        if rotate_frame
                if isequal(mod(i,mod_number),0)
                    k = mod_number;
                else
                    k = mod(i,mod_number);
                end
                 plot(XX(w,a(q),r(q)),YY(w,a(q),r(q)),'color',color_code(k,:))
                 hold on
                 clear w;
        else
                if isequal(mod(i,mod_number),0)
                    k = mod_number;
                else
                    k = mod(i,mod_number);
                end
                 plot(B(w,a(q),r(q)),C(w,a(q),r(q)),'color',color_code(k,:))
                 hold on
                 clear w;
                 
                 
        end
    end
   
   
    
   % PLOT THE COORDINATE AXES
        ttt=max(window);
        b=linspace(-ttt,ttt,1100);
        plot(0,b,'black',b,0,'black')
        
        
        
        
 % PLOT DEGENERACY LINE - IF NEAR THE TRANSITION
 % PLOT PLANE
 % PLOT MRS
 
   if rotate_frame
       
        % PLOT DEGENERACY LINE
         if(~isequal(h,0))
           for k = 1:length(h(:,1))
               if(abs(a(q)-h(k,1))<10^(-1)) %then A is near Atrans
                   c=CC(Btrans(k,r(q)),b,Ctrans(k,r(q)),k);
                  plot((b+c)/sqrt(2),(b-c)/sqrt(2),'LineStyle','--','Color','black','LineWidth',1)
                  clear c 
               end
           end
          end   
       % PLOT rotated A+B+C=0 IN VIOLET 
        b;     c=-(a(q)+b);
        plot((b+c)/sqrt(2),(b-c)/sqrt(2),'color',[0.5 0 1]); %plot coordinate axes and the line C = -B
       % NOW PLOT THE MRS - which only exists when A > -(R-1)/R > 0 ...thm 2.2
         if ( a(q) > -(r(q)+1)/r(q) ) && ( a(q) > 0 ) 
            line([((0+a(q))/sqrt(2)), ((a(q)+0)/sqrt(2))],...
                 [((0-a(q))/sqrt(2)) ,((a(q)-0)/sqrt(2))],...
                 'Marker','.','LineStyle','--','Color','black')
            line([(a(q))/sqrt(2) (-a(q))/sqrt(2)],...
                [ (a(q))/sqrt(2) (a(q))/sqrt(2)],...
                'Marker','.','LineStyle','--','Color','black')
            line([-a(q)/sqrt(2) a(q)/sqrt(2)],...
                [-a(q)/sqrt(2) -a(q)/sqrt(2)],...
                'Marker','.','LineStyle','--','Color','black')   
         end
         % SET AXIS LABELS    
        xlabel('X axis','FontSize',18)
        ylabel('Y axis','FontSize',18)
    
   else
       
            %plot degeneracy line
          if(~isequal(h,0))
           for k = 1:length(h(:,1))
               if(abs(a(q)-h(k,1))<10^(-1)) %then A is near Atrans
              c=CC(Btrans(k,r(q)),b,Ctrans(k,r(q)),k);
              plot(b,c,'LineStyle','--','Color','black','LineWidth',1)
                clear  c
               end
           end
          end


        % PLOT A+B+C=0 IN VIOLET
            plot(b,-(a(q)+b),'color',[0.5 0 1]); %plot coordinate axes and the line C = -(A+B)
        % NOW PLOT THE MRS - which only exists when A > -(R-1)/R > 0 ...thm 2.2
         if ( a(q) > -(r(q)+1)/r(q) ) && ( a(q) > 0 ) 
            line([0 a(q)],[a(q) 0],'Marker','.','LineStyle','--','Color','black')
            line([a(q) 0],[ 0 -a(q)],'Marker','.','LineStyle','--','Color','black')
            line([-a(q) 0],[0 a(q)],'Marker','.','LineStyle','--','Color','black')
         end
         % SET AXIS LABELS    
        xlabel('B axis','FontSize',18)
        ylabel('C axis','FontSize',18)
   end
   
  

   
    

     
     
     
     
    % PRINT PARAMETER VALUES TO THE TITLE
        s = sprintf('R = %3.5f   and   A = %3.5f',r(q),a(q));  % %1.5f    ,r(q)
        %title(s,'BackgroundColor',[1 1 0.5020])%old setting for yellow
        title(s,'FontSize',18)
        
        
    
        
        
    % SET AXIS WINDOW - FOR A SCALAR OR AXIS VECTOR
        if isequal(length(window),1)
            axis(window*[-1 1 -1 1])
        else
            axis(window)
        end

     % SET AXIS WINDOW - with MRS Fixed   
        if fix_mrs

                     dA=(A(2)-A(1))/(nframes-1);
                     a_rate = dA/a(q);
                      axis(z*[-1 1 -1 1])
                      z=z+z*a_rate;
        end
                
                
    
    grid on;
    axis square

 
  %create a text box ui control to display mouse position
  hText = uicontrol(gcf,'Style','text',...
                        'BackgroundColor',[1 1 1],...
                        'ForegroundColor',[0 0 0]);
  % set figure props and define mouse motion callback function
  set( gcf,'pointer','crosshair','WindowButtonMotionFcn',{@hFigure_MotionFcn, hText});
  hold off
 Frames(:,q) = getframe(gcf);

end
  
% UNCOMMENT IF YOU WANT TO PLAY IT BACK - BUT YOU WILL LOSE THE AXIS
% NUMBERS, A MATLAB FLAW.
%
% movie(Frames,2,speed_control);  %movie(Frames,replay,speed) where 'speed' determines how fast
% the movie is played and 'replay' is the number of times the movie is replayed on the figure window
% see help 'movie'
return


function h = get_trans_points(R)
    % 
    % returns an array whose rows are (A,B,C) triplets (corresponding to equations 2.6, and 2.7
    % in Mahaffy, et.al pg 783) at which there exists a transition, for the
    % given value of R
  j = 1;
  [A,B,C] = get_min_stab_point(R);
     fprintf('\n')
    fprintf('The min point of stability is (%5.5f ,%5.5f, %5.5f)\n',A,B,C)
     fprintf('\n')
      fprintf('TRANSITION POINTS (A,B,C) for R = %1.4f \t\t DEGENERACY LINE\n',R)
    fprintf('---------------------------------------------------------------\n')
    tol=10^16;
    hh=zeros(200,3);
    Atrans =0;
    flag =0;
h=[];
    while (abs(Atrans) < tol) 
        Atrans = -j*pi/(1-R)*cot(j*R*pi/(1-R));

        Btrans = (-1)^j*( (1-R)*cos(j*pi*R/(1-R)) -j*pi*R*csc(j*R*pi/(1-R)) ) / ((1-R)^2);

        Ctrans = (-1)^(j+1)*( (1-R)*cos(j*pi/(1-R)) -j*pi*csc(j*pi/(1-R)) ) / ((1-R)^2);

        hh(j,1) = Atrans;
        hh(j,2) = Btrans;
        hh(j,3) = Ctrans;

            temp = sign(Atrans);
            if isequal(temp,1)
                flag = 1;
            end
            if isequal(temp,-1) && isequal(flag,1)
                fprintf('\n\n')
                return
            end
            if(abs(Atrans) < tol) 
                fprintf('j =%2d   (%5.5f ,%5.5f, %5.5f)    \t\tB %+2d\b C = %5.5f \n',...
                 j+1, Atrans,Btrans,Ctrans,(-1)^j,Btrans+(-1)^j*Ctrans)
            h(j,:)=hh(j,:);
            end


        j = j+1;
        
    end
    if isempty(h)
        h=0;
    end
    fprintf('\n\n')
    h;
 return
 
 
 
 function [A,B,C] = get_min_stab_point(R)

%returns minimum point of stability for the given delay ratio, R.
% according to thm 2.2, p782

A = -(R+1)/R;
B = R/(R-1);
C = 1/(R*(1-R));

return

    function hFigure_MotionFcn(hFigure,event_data,hText)
        %get mouse position
        set(gcf,'Units','Normalized')
        mouse_pos = get(gcf,'CurrentPoint');
%         pos(1) = mouse_pos(1)-341.25;
%          pos(2) = mouse_pos(2)-261.75; 
        set( hText,...
             'Position', [1 1 200 30],...
             'String',   ['B = ' num2str(mouse_pos(1)) ...
                           ', C = ' num2str(mouse_pos(2))]); 



