function xdot=mlds2cellfun(t,x,flag,tdel)
%% Carol Lucas
%% Cartoon shown in Figure 6.8

% a simple model for motion and directionally-selective cell
% Morris-Lecar cell (1st layer) that receives two separate inputs
%  -- one inhibitory and one excitatory; temporally offset.
% 2nd-layer ML cell receives same excitatory input as 1st-layer but
% also gets inhibitory input from 1st-layer neuron
% this ML model uses the standard ML params from JR/GBE chapt in 
% Koch and Segev and in JEK chapt%

% For inputs to (1st layer) target neuron N! fix the start time of
% excitation and delay the
% start of inhibition by tdel (if tdel<0; inhibition precedes excitation)

% delay the excitatory input from receptor R2 to neuron N2 by tdel2 
%   (a fixed propagation delay)
%
v1=x(1,1);
w1=x(2,1);
v2=x(3,1);
w2=x(4,1);
se=x(5,1);
si=x(6,1);
se2=x(7,1);
s21=x(8,1);

t0=100;teon=5;tion=5;vrest=-60;tdel2=10;
gsyne=2;gsyni=2.5;vsyne=100;vsyni=-84;
gsyne2=2;gsyn21=2.5;vsyn21=-84;alpha21=0.5;beta21=0.1;
vk=-84;vl=-60;vca=120;
i=0;gk=8;gl=2;c=20;
va=-1.2;vb=18;vc=2;vd=30;phi=.04;gca=4.;
vt=20;vshp=2;alphae=1;betae=0.3;alphai=0.5;betai=.1;

ve=100*(t>=t0)*(t<=teon+t0)+vrest;
vi=100*(t>=t0+tdel)*(t<=tion+t0+tdel)+vrest;
ve2=100*(t>=t0+tdel2)*(t<=teon+t0+tdel2)+vrest ;
%
isyne=gsyne*se*(vsyne-v1);
isyni=gsyni*si*(vsyni-v1);
isyne2=gsyne2*se2*(vsyne-v2);
isyn21=gsyn21*s21*(vsyn21-v2);

minf1=.5*(1+tanh((v1-va)/vb));
winf1=.5*(1+tanh((v1-vc)/vd));
tauw1=1/cosh((v1-vc)/(2*vd));

minf2=.5*(1+tanh((v2-va)/vb));
winf2=.5*(1+tanh((v2-vc)/vd));
tauw2=1/cosh((v2-vc)/(2*vd));

sinfvi=1/(1+exp(-(vi-vt)/vshp));
sinfve=1/(1+exp(-(ve-vt)/vshp));
sinfve2=1/(1+exp(-(ve2-vt)/vshp));
sinfv1=1/(1+exp(-(v1-vt)/vshp));

xdot(1,1) = (i-gca*minf1*(v1-vca)-gk*w1*(v1-vk)-gl*(v1-vl)+isyne+isyni)/c;
xdot(2,1) = phi*(winf1-w1)/tauw1;
xdot(3,1) = (i-gca*minf2*(v2-vca)-gk*w2*(v2-vk)-gl*(v2-vl)+isyne2+isyn21)/c;
xdot(4,1) = phi*(winf2-w2)/tauw2;
xdot(5,1)=alphae*sinfve*(1-se)-betae*se;
xdot(6,1)=alphai*sinfvi*(1-si)-betai*si;
xdot(7,1)=alphae*sinfve2*(1-se2)-betae*se2;
xdot(8,1)=alpha21*sinfv1*(1-s21)-beta21*s21;