L-type-Ca-Dend-Spike.hoc

strdef modelName, loadProgram, cellName, outputFile, cellPath, loadProgram
objref trunc, secR, fi,vC, mbSec
objref sh, axonal, dendritic, dendritic_only, initZone, stimSec
objref trunc,mbSec, middleSec
objref mydendrites

load_file("stdgui.hoc")
load_proc("nrnmainmenu")
objref sh, st, axonal, dendritic, dendritic_only
objref mydendrites,tuftS,temp,mytuft, mbSec
create soma
access soma

tstop = 500
steps_per_ms = 40
dt = 0.025


// --------------------------------------------------------------
// passive & active membrane 
// --------------------------------------------------------------

ra        = 150
global_ra = ra
rm        = 10000
c_m       = 1
cm_myelin = 0.04
g_pas_node = 0.02

v_init    = -75
celsius   = 37

Ek = -90
Ena = 60


gna_dend = 20
gna_node = 30000
gna_soma = 1500

gkv_axon = 2000
gkv_soma = 200

gca = .5
gkm = 2.2
gkca = 2.5

gca_soma = gca
gca_dend = gca
gkm_soma = gkm
gkca_soma = gkca
git_soma = 0.0000
gkv_dend = 20
git_dend = 0.0008
gil_dend = 0.00

// --------------------------------------------------------------
// Axon geometry
//
// Similar to Mainen et al (Neuron, 1995)
// --------------------------------------------------------------

n_axon_seg = 5

create soma,iseg,hill,myelin[2],node[2]

proc create_axon() {

  create iseg,hill,myelin[n_axon_seg],node[n_axon_seg]

  soma {
    equiv_diam = sqrt(area(.5)/(4*PI))

    // area = equiv_diam^2*4*PI
  }
  if (numarg()) equiv_diam = $1

  iseg {                // initial segment between hillock + myelin
     L = 15
     nseg = 5
     diam = equiv_diam/10        // see Sloper and Powell 1982, Fig.71
  }

  hill {                
    L = 10
    nseg = 5
    diam(0:1) = 4*iseg.diam:iseg.diam
  }

  // construct myelinated axon with nodes of ranvier

  for i=0,n_axon_seg-1 {
    myelin[i] {         // myelin element
      nseg = 5
      L = 100
      diam = iseg.diam         
    }
    node[i] {           // nodes of Ranvier
      nseg = 1
      L = 1.0           
      diam = iseg.diam*.75       // nodes are thinner than axon
    }
  }

  soma connect hill(0), 0.5
  hill connect iseg(0), 1
  iseg connect myelin[0](0), 1
  myelin[0] connect node[0](0), 1

  for i=0,n_axon_seg-2  { 
      node[i] connect myelin[i+1](0), 1 
      myelin[i+1] connect node[i+1](0), 1
  }
}

// --------------------------------------------------------------
// Spines
// --------------------------------------------------------------

      // Based on the "Folding factor" described in
      // Jack et al (1989), Major et al (1994)
      // note, this assumes active channels are present in spines 
      // at same density as dendrites

spine_dens = 1
      // just using a simple spine density model due to lack of data on some 
      // neuron types.

spine_area = 0.83 // um^2  -- K Harris

proc add_spines() { local a
  forsec $o1 {
    a =0
    for(x) a=a+area(x)

    F = (L*spine_area*spine_dens + a)/a

    L = L * F^(2/3)
    for(x) diam(x) = diam(x) * F^(1/3)
  }
}



proc init_cell() {

  // passive
  forall {
    insert pas
    Ra = ra 
    cm = c_m 
    g_pas = 1/rm
    e_pas = v_init
  }

  // exceptions along the axon
  forsec "myelin" cm = cm_myelin
  forsec "node" g_pas = g_pas_node

  // na+ channels
  forall insert na
  forsec dendritic gbar_na = gna_dend
  forsec "myelin" gbar_na = gna_dend
  hill.gbar_na = gna_node
  iseg.gbar_na = gna_node
  forsec "node" gbar_na = gna_node

  // kv delayed rectifier channels
  iseg { insert kv  gbar_kv = gkv_axon }
  hill { insert kv  gbar_kv = gkv_axon }
  soma { insert kv  gbar_kv = gkv_soma }

  // dendritic channels
  forsec dendritic {
    insert kv    gbar_kv = gkv_dend
    insert km    gbar_km  = gkm
    insert kca   gbar_kca = gkca
    insert ca    gbar_ca = gca_soma
    insert cad
    insert it2   gcabar_it2 = git_dend
    insert sca	 gbar_sca=gil_dend
  }

  soma {
	insert na          gbar_na = gna_soma
	insert kv          gbar_kv = gkv_soma
	insert km          gbar_km = gkm_soma
	insert kca         gbar_kca = gkca_soma
	insert ca          gbar_ca = gca_soma
	insert it          gbar_it = git_soma
  }

 
  forall if(ismembrane("k_ion")) ek = Ek
  forall if(ismembrane("na_ion")) {
    ena = Ena
    // seems to be necessary for 3d cells to shift Na kinetics -5 mV
    vshift_na = -5
  }
  forall if(ismembrane("ca_ion")) {
    eca = 140
    ion_style("ca_ion",0,1,0,0,0)
    vshift_ca = 0
  }
}

proc load_3dcell() {

// $s1 filename

  aspiny = 0
  forall delete_section()
  xopen($s1)
  access soma

  dendritic = new SectionList()

   // make sure no compartments exceed 50 uM length
  forall {
    diam_save = diam
    n = L/50
    nseg = n + 1
    if (n3d() == 0) diam = diam_save
    dendritic.append()
  }    

  dendritic_only = new SectionList()
  forsec dendritic dendritic_only.append()
  soma  dendritic_only.remove()

  create_axon()

  init_cell()

  if (!aspiny) add_spines(dendritic_only,spine_dens)
}
   load_3dcell("j8.hoc")

objectvar stim,stim2,stim3,stim4,stim5,stim6,stim7,stim8,stim9,stim10

I=10
dur=1
soma stim = new IClamp(0.5)

stim.del = 100
stim.dur = dur
stim.amp = I

soma stim2 = new IClamp(0.5)

stim2.del = 110
stim2.dur = dur
stim2.amp = I

soma stim3 = new IClamp(0.5)

stim3.del = 120
stim3.dur = dur
stim3.amp = I

soma stim4 = new IClamp(0.5)

stim4.del = 130
stim4.dur = dur
stim4.amp = I

soma stim5 = new IClamp(0.5)

stim5.del = 140
stim5.dur = dur
stim5.amp = I


soma stim6 = new IClamp(0.5)

stim6.del = 150
stim6.dur = dur
stim6.amp = I


/*
somaA stim7 = new IClamp(0.5)

stim7.del = 166
stim7.dur = dur
stim7.amp = I

somaA stim8 = new IClamp(0.5)

stim8.del = 170
stim8.dur = dur
stim8.amp = I

somaA stim9 = new IClamp(0.5)

stim9.del = 180
stim9.dur = dur
stim9.amp = I

somaA stim10 = new IClamp(0.5)

stim10.del = 190
stim10.dur = dur
stim10.amp = I

*/

/********************************************************************/

proc update_init(){
	finitialize(v_init)
	fcurrent()
	for (x) {
		e_pas(x)=v(x)
		e_pas(x)=e_pas(x)+(ina(x)+ik(x)+i_hd(x))/g_pas(x)
	}

}


objectvar g[20]         // max 20 graphs
ngraph = 0

proc addgraph() { local ii  // define subroutine to add a new graph
                // addgraph("variable", minvalue, maxvalue)
    ngraph = ngraph+1
    ii = ngraph-1
    g[ii] = new Graph()
    g[ii].size(0,tstop,$2,$3)
    g[ii].xaxis()
    g[ii].yaxis()
    g[ii].addvar($s1,1,0)
    g[ii].save_name("graphList[0].")
    graphList[0].append(g[ii])
}

if(ismenu==0) {
  nrnmainmenu()         // create main menu
  nrncontrolmenu()      // crate control menu
  ismenu=1
}



addgraph("soma.v(0.5)",-75,90)
addgraph("a3_12122222112.v(0.5)",-75,90)