%Hodgkin-Huxley:
%
a=-59%-105; %-88.18951%start membrane potential
for i=1:15
ySS=[0.3177 0.0530 0.5962 a];%[0.5 0.5 0.5 a];
[t,x]=ode45('hhode',[0 200],ySS); %%HH equation
subplot(2,1,1);
%
if max(x(:,4))>0 && max(x(:,4))<50 % AP
    plot(t, x(:,4),'r')
    drawnow;
    hold on;
    

elseif max(x(:,4))>50 % too high membrane potential
    plot(t, x(:,4),'b')
    drawnow;
    hold on;
    
else
    plot(t, x(:,4),'k') %no AP
    drawnow;
    hold on;
    
end
    title(['Threshold Behavior, voltage: ', num2str(a)]);
    xlabel('Time (ms)');
    ylabel('Transmembrane Voltage (mV)');
    xlim([-1 200]);

subplot(2,1,2);
%iongate activities
plot(t, x(:,3),'r')
hold on
plot(t, x(:,2),'k')
plot(t, x(:,1),'b')
legend('h: Na gate inactivity','m: K gate activity','n: Na gate activity');
title('Ion gate activities');
xlabel('Time (ms)');
ylabel('Activity');
ylim([0 1]);
xlim([-1 200]);
hold off;

waitforbuttonpress;
a=a+2;%0.000001; %

end
%}
%{
resting = [0.3177    0.0529    0.5962  -52];
ySS=resting;
[t,x]=ode45('hhode',[0 20],ySS);
figure
if max(x(:,4))>0 && max(x(:,4))<50
    plot(t, x(:,4),'r')
    drawnow;
    hold on;
    
elseif max(x(:,4))>50
    plot(t, x(:,4),'b')
    drawnow;
    hold on;
    
else
    plot(t, x(:,4),'k')
    drawnow;
    hold on;
end
title('Threshold Behavior');
    xlabel('Time (ms)');
    ylabel('Transmembrane Voltage (mV)');
    xlim([-1 20]);

figure;
plot(t, x(:,3),'r')
hold on
plot(t, x(:,2),'k')
plot(t, x(:,1),'b')
legend('h: Na gate inactivity','m: K gate activity','n: Na gate activity');
title('Ion gate activities');
xlabel('Time (ms)');
ylabel('Activity');
ylim([0 1]);
xlim([-1 20]);
%}