// mit4.hoc

// 4 + 3 compartment mitral cell model. One compartment represents the 
// glomerulus, one the primary dendrite, one the secondary dendrites 
// and one the soma. The other three are small linking compartments
// which represent almost all the axial resistance in the model.

// Andrew Davison, The Babraham Institute
// 10-Jun-1999
// modified 28-Feb-2002: somastim and glomstim added to allow somatic or
// glomerular stimulation to be chosen without having to change alphas/alphag

// Note that although in theory the output is independent of total area,
//  due to rounding errors this only holds for Atotal > ~ 1e4

xopen("$(NEURONHOME)/lib/hoc/noload.hoc") // standard run tools
xopen("tabchannels.hoc")

// Set parameters
Atotal		= 100000	// um2
Len		= 100		// um
RM		= 100000	// ohm.cm2
Erest		= -65		// mV
p 		= 0.051
q 		= 0.084
r 		= 0.328
gpg 		= 5.86e-5	// S.cm-2
gsp 		= 5.47e-5	// S.cm-2
gsd 		= 1.94e-4	// S.cm-2
alphas 		= 1.37
alphag		= 1.85
Ifull		= 0.4		// nA
somastim	= 1
glomstim	= 0

// Create cell
create soma, glom, prim, dend, s2d, s2p, p2g
access soma

soma connect s2p(0),0
s2p connect prim(0),1
prim connect p2g(0),1
p2g connect glom(0),1
soma connect s2d(0),1
s2d connect dend(0),1

objref cvode
cvode  = new CVode(1)

// Insert stimuli
objref sstim, gstim
injcurrdens = 0.4/100072		// injected current density (nA.um-2)
soma {
   sstim = new IClamp(0.5)
   sstim.del = 50
   sstim.dur = 15000
}
glom {
   gstim = new IClamp(0.5)
   gstim.del = 50
   gstim.dur = 15000
}

// Define procedures

proc set_ra() {  			// the argument is the conductance (S.cm-2)
   Ra = (PI*1e4)/(4*Atotal) * ( 1/$1 )	// ohm.cm
}

proc set_size() {			// the argument is the membrane area (um2)
   diam = $1/(PI*Len)			// um
}

proc change_params() {

      injcurrdens = Ifull/100072	//note that 100072 um2 is the total membrane area of the full model
      Asoma = p*Atotal
      Aglom = q*Atotal
      Aprim = r*Atotal
      Adend = Atotal - Asoma - Aglom - Aprim

      soma {
        set_size(Asoma)
        sstim.amp = somastim*alphas*injcurrdens*Atotal
      }

      glom {
        set_size(Aglom)
	gstim.amp = glomstim*alphag*injcurrdens*Atotal
      }

      prim {
        set_size(Aprim)
      }  

      dend {
        set_size(Adend)
      }

      s2d { set_ra(gsd) }
      s2p { set_ra(gsp) }
      p2g { set_ra(gpg) }
}


// Set cell properties

forsec "*2*" {
   L = 1
   diam = 1
}

soma {

   insert pas
   insert nafast
   insert kfasttab
   insert kslowtab
   insert kA
   insert kca3
   insert lcafixed
   insert cad
   depth_cad = 8

   L			= Len
   Ra			= 1e-7		// ohm.cm This value can be set to any small number
					// although if it too small rounding errors are a problem
					// The range 1e-7 to 100 is basically OK
   e_pas 		= Erest		// reversal potential mV
   g_pas 		= 1/RM		// membrane conductance S.cm-2

   gnabar_nafast 	= 0.1532	// S.cm-2
   gkbar_kfasttab 	= 0.1956
   gkbar_kslowtab 	= 0.0028
   gkbar_kA 		= 0.00587
   gkbar_kca3 		= 0.0142
   gcabar_lcafixed 	= 0.0040
}

glom {

   insert pas
   insert kslowtab
   insert lcafixed
   insert cad
   depth_cad = 8

   L			= Len
   Ra			= 1e-7
   e_pas 		= Erest
   g_pas		= 1/RM

   gkbar_kslowtab	= 0.020
   gcabar_lcafixed	= 0.0095
}

prim {
   
   insert pas
   insert nafast
   insert kfasttab
   insert kslowtab
   insert lcafixed
   insert cad
   depth_cad = 8

   L			= Len
   Ra			= 1e-7
   e_pas 		= Erest
   g_pas		= 1/RM

   gkbar_kfasttab	= 0.00123
   gnabar_nafast	= 0.00134
   gkbar_kslowtab	= 0.00174
   gcabar_lcafixed	= 0.0022
}

dend {

   insert pas
   insert kfasttab
   insert nafast
   L			= Len
   Ra			= 1e-7
   e_pas 		= Erest
   g_pas		= 1/RM

   gkbar_kfasttab	= 0.0330
   gnabar_nafast	= 0.0226
}

change_params()				// set areas and link resistances


// set reversal potentials, etc.

forall if (ismembrane("ca_ion")) {
        eca = 70	// mV
	cai = 0.00001	// mM 
	cao = 2		// mM
	ion_style("ca_ion",3,2,0,0,1)
}

forall if (ismembrane("na_ion")) {
	ena = 45	// mV
}

forall if (ismembrane("k_ion")) {
	ek  = -70	// mV
}