function V=expIF_neuron(time,varargin)
% V=IF_neuron(TIME,VPAR,MPAR,TYPE,CURRENT)
% VPAR, MPAR, TYPE, CURR all optional. TIME is in ms.
% VPAR=[Vrest,Vthresh,Vreset,Vpeak,Vinit]; in mV, respectively.
% MPAR=[Cm,Rm,Tref]; in nF, MOhms, and ms, respectively.
% TYPE=[1 or other]; 1 indicates that sineCurrent should be used, any other
% number indicates that stepCurrent should be used.
% CURR=[fsin,Imax,phi]; in kHz, nA, and rads, respectively, if sineCurrent
% and [T_STEP,I_STEP]; in ms and nA, respectively, if stepCurrent.
if nargin>1
    Vpar=varargin{1};
else
    Vpar=[-60 -50 -70 40 -60]; % default values
end
if nargin>2
    Mpar=varargin{2};
else
    Mpar=[1 1 10]; % default values
end
if nargin>3
    type=varargin{3};
else
    type=1; % default value is for sineCurrent
end
if nargin>4
    Curr=varargin{4};
else
    Curr=[.05 2 0]; % default values for sineCurrent
end

dt=(time(end)-time(1))/(numel(time)-1); % delta t
R=Mpar(2); % Resistance
C=Mpar(1); % Capacitance
Sref=round(Mpar(3)/dt); % Number of refractory steps

if type==1;
    for k=1:numel(time);
        I(k)=sineCurrent(time(k),Curr(1),Curr(2),Curr(3));
        V_inf(k)=Vpar(1)+(I(k)*R);
    end
else
    for k=1:numel(time);
        I(k)=stepCurrent(time(k),Curr(1),Curr(2));
        V_inf(k)=Vpar(1)+(I(k)*R);
    end
end

V=zeros(1,numel(time));
V(1)=Vpar(5);
counter=1;

for k=1:(numel(V)-1); % Finding the voltage for the steps past Sref.
    if V(k)==Vpar(3) & counter<Sref
        V(k+1)=Vpar(3);
        counter=counter+1;
    elseif V(k)>=Vpar(2) % Spiking threshold
        V(k)=Vpar(4);
        V(k+1)=Vpar(3);
        counter=1;
    else
        V(k+1)=V_inf(k)+(V(k)-V_inf(k))*exp(-dt/(R*C));
    end
end

Published with MATLAB® 7.10