{load_file("stdlib.hoc")}
celsius=36
{load_file("soma.hoc")}

proc default_var() {localobj s, pc
    s = new String()
    if (name_declared($s1) == 0) {
	if (argtype(2) == 2) {
	    sprint(s.s, "~strdef %s", $s1)
	    execute(s.s)
	    sprint(s.s, "%s = \"%s\"", $s1, $s2)
	}else{
	    hoc_ac_ = $2
	    sprint(s.s, "%s = hoc_ac_", $s1)
	}
	execute(s.s)
    }
    pc = new ParallelContext()
    if (pc.id == 0) {
	sprint(s.s, "print \"%s  \", %s", $s1, $s1)
	execute(s.s)
    }
}

///// Define default variables /////
default_var("gnav13",1000)
default_var("gnav17",1000)
default_var("gnav19",1000)
default_var("gkdr",1000)
default_var("gka",1000)
default_var("gcan",1000)
default_var("gkcaf",1000)
default_var("gcansc",1000)
default_var("gkcas",1000)
default_var("gih",1000)
default_var("vaih",0)
default_var("vs1ih",100)
default_var("vs2ih",100)
default_var("tc1ih",100)
default_var("tc2ih",100)
default_var("jtcmih",100)
default_var("stimamp",300)

///// conductance changes /////
access soma
gbar_kca_slow = 4e0
gbar_nav13 *= gnav13/1000.0
gnabar_nav17 *= gnav17/1000.0
gbar_nap *= gnav19/1000.0
gbar_can *= gcan/1000.0
gbar_kdr *= gkdr/1000.0
gbar_ka *= gka/1000.0
gbar_kca_fast *= gkcaf/1000.0
gbar_kca_slow *= gkcas/1000.0
gbar_cansc *= gcansc/1000.0
gbar_ih *= gih/1000.0
vhalf_ih += vaih
vslope1_ih *= vs1ih/100.0
vslope2_ih *= vs2ih/100.0
tmc1_ih *= tc1ih/100.0
tmc2_ih *= tc2ih/100.0
jtmc_ih *= jtcmih/100.0

//integration time step
dt = 0.01
tstop = 700

//inserting current clamp at soma
access soma
objectvar stim1
stim1 = new IClamp(0.5)
stim1.del = 100
stim1.dur = 500
stim1.amp = stimamp/1000.0

//for writing result files
objref vrec[1]
vrec[0] = new Vector()
{vrec[0].record(&soma.v(0.5), dt)}
objref tvec
tvec = new Vector()
{tvec.record(&t, dt)}
objref all
all = new File()
strdef datei

//graph
objref gmem
gmem = new Graph()
gmem.size(0, tstop, -100, 50)
{gmem.addvar("Membrane Potential", &soma.v(0.5), 2, 1)}
gmem.begin()

//defining some functions
proc run_sn() {
    sprint(datei, "./fig3currents.txt")
    all.wopen(datei)
    soma v = -55
    finitialize()
    while (t < tstop) {
	fadvance()
	all.printf("%g %g %g %g %g %g %g %g %g %g %g\n", t, cai, jina13_nav13, jina17_nav17, ica, jikdr_kdr, jika_ka, jikcaf_kca_fast, iother2, jikcas_kca_slow, jiih_ih)
	gmem.plot(t)
    }
    gmem.flush()
    all.close
    write_file()
}

proc write_file() {
    sprint(datei, "./fig3mem.txt")
    all.wopen(datei)
    for i=0, vrec[0].size()-1 {
        all.printf("%g %g\n", tvec.x(i), vrec[0].x(i))
    }
    all.close
}

run_sn()
//quit()