function [t,y]=hh1_nk

% hodgkin-huxley model usando las convenciones originales, en particular V es la despolaricaciòn con respecto al valor de reposo
% por eso ploteamos V-60 pera llevarlo a las convenciones actuales
%Cómo habría que cambiar las ecuaciones si queremos llevar todo a la convención actual?


%------------- default values ------------- 
 
global Cm gNa gK gL VNa VK VL;

Cm = 1;             % membrane capacitance [?mF/cm2]
gNa = 120;          % Maximum possible Sodium Conductance (about 120 mOhms^-1/cm2) 
gK = 36.001;            % Maximum possible Potassium Conductance (about 36 mOhms^-1/cm2) 
gL = 0.3;           % Maximum possible Leakage Conductance (about .3 mOhms^-1/cm2) 


%VNa = 55;           % equilibrium potential [mV] for Na
%VK = -72;           %                        ... for K
%VL = -60;         %                        ... for leakage

VNa = 115;          % equilibrium potential [mV] for Na
VK = -12;           %                        ... for K
VL = 10.6;          %                        ... for leakage


%------------- time period------------- 

tspan = [0; 150];           

y0=[ -0;  0.3177  ;  0.0529   ; 0.5961]


[t,y] = ode45(@f,tspan,y0);


%------------------ hodgkin and huxley equations ------------------

function dydt = f(t,y)
     global Cm gNa gK gL VNa VK VL;
     V=y(1); % V: the membrane potential (voltage) 
     n=y(2); % n: the probability that 1 of 4 activation particles has influenced the state of the K gate.
     m=y(3); % m: the probability that 1 of the 3 required activation particles has contributed to the activation of the Na gate (m^3 : the probability that all 3 activation particles have produced an open channel) 
     h=y(4); % h: the probability that the 1 inactivation particle has not caused the Na gate to close 

     % CONVENTION: outward flux of positive charge is positive current.
 
     Inat = ...                          % the "natural" cross-membrane current =
              gNa * m^3*h * (V-VNa)...   % sodium current +
            + gK *  n^4 *   (V-VK) ...   % potassium current +
            + gL *          (V-VL)  ;    % leakage current +

     if t>50 || t < 45 %una manera de hacer un pulso rectangular 
       Itot=-1.*Inat;
     else
       Itot = 10 - Inat;     % applied current (experimental or synaptic input)
     end
     % change of the membrane potential V as a function of m, h, n
     dVdt = Itot/Cm;  % Cm*dVdt=Itot 
    
     % change of the activation particle probability n for potassium
     dndt = -1 * 0.125*exp( 0.0125*V) * n + (0.01*(10-V)/(exp(-0.1*V+1  )-1)) * (1-n);

     % change of the activation particle probability m for sodium
     dmdt = -1 * 4    *exp(-0.0556*V) * m +  (0.1 *(25-V)/(exp(-0.1*(V-25))-1)) * (1-m);

     % change of the inactivation particle probability (1-h) for sodium
     dhdt = -1* 1/(1+exp(-0.1*V+3))  * h +  0.07*exp(-0.05*V)                * (1-h);
     
     dydt = [ dVdt       ;       dndt         ;     dmdt         ;     dhdt];
     
     
function m_inf=mcalc(V)
      z= (0.1.*(25-V)./(exp(-0.1.*(V-25))-1))./(4.*exp(-0.0556*V));
      m_inf =z./(1+z);