% Morris_Lecar.m
% Carol Lucas
% Simulations for Figure 2.9
% Parameters come from Table 2.4 with starting V of -20 mv to match figure
% User should be certain parameters are initialized with GMorris_Lecar_Init.m
% Figures are not cleared so overlapping plots can be made

global Vk Vl Vca gca gk gl V1 V2 V3 V4 c phi tend Vinit Winit I

clear t;

C=input('Hit any number to erase figures: ');
if length(C)~=0; 
	figure(1);clf;figure(2);clf;
	end

tstart=0;
k=1;t=[tstart:.25:tend];
xint=[Vinit;Winit];

[t,x]=ode45('Gmlfun',t,xint);

figure(1);subplot(211);plot(t,x(:,1));axis([0 tend -75 75]);
xlabel('time (ms)'),ylabel('voltage (mv)');hold on;
subplot(212);plot(t,x(:,2));axis([0 tend 0 1]);
xlabel('time (ms)'),ylabel('W=fraction of open channels');hold on;
pause(5);
mess=' Off to get phase plane and nullclines'

%% Plot phase plane and nullclines
figure(3);clf;    %% Subroutine plots independently
[Vs,Ws,lams,VRAY,nullV,nullw]=GMorris_Lecar_Jac;

figure(2);plot(x(:,1),x(:,2),x(1,1),x(1,2),'md',VRAY,nullV,':',VRAY,nullw,'--');
axis([-75 75 0 1]);
ylabel('W=fraction of open channels'),xlabel('voltage (mv)');
title('Phase Plane and Nullclines');hold on;

C=input('Hit any number to get vector field:  ');
if length(C)~=0
 	Vray=[-60:10:60];Nray=[0:.1:1];figure(2);
	vector_field('Gmlfun',Vray,Nray,2);
	end;