{load_file("nrngui.hoc") } //I assume these libraries will load from the first /nrn path
{load_file("stdrun.hoc")}
{load_file("stdlib.hoc")}
{load_file("netparmpi.hoc")}
{load_file("./externals.hoc")}
{load_file("./templates/syn.tem")}
{load_file("./templates/gap.tem")}
{load_file("./templates/iapp.tem")}
{load_file("./templates/pyrkop.tem")}
{load_file("./templates/bwb.tem")}
{load_file("./templates/ok.tem")}
{load_file("./parameters/synapses.tem")}
{load_file("./parameters/manycells.tem")}
{load_file("./templates/TGbignet2.tem")}
if(doextra){
//print "doing EXTRACELLULAR stuff..........."
{load_file("./interpxyz.hoc")} // only interpolates sections that have extracellular
{load_file("./setpointers.hoc")} // has the function grindaway() in interpxyz.hoc to set up pointers
{load_file("./field.hoc")} //functions to calculate the extracellular voltage
{load_file("./calcrxc_a.hoc")} //function to create the transfer resistance from each section
}
strdef cmd
numiters=1
objref pr, pc, nil, pyrspiketau_vec, baskspiketau_vec, baskconnvector, ipspg_vec, fih_progress, hines
hinest1=startsw()
hinest2=startsw()
gapstyle=96 //split baskets right, many tests
//19 lowered spikedur to 75
//20 from here onward the _b will actually be basketnoise, rather than conn %. I've switched conn to 100
//g97 is to signify switching to this new method
gapstyle=97
gapstyle=0 //this is to signify running it with the newer version, with the noisesyn mods added
//29 ungapped baskets, 2 rows of 10 inline gaps between pyr cells (.01 - high) 1-10 and 41-50
antennaDC=1.3 //constant current added to antenna cells
realrunFlag=1
shortspikeFlag=1
//celsius = 34 //for some reason this line gives an error for being out of range
iteration=0
{cvode.active(0)}
//the following variables determine the structure of the ramp in noise
/* ramp_on=100
ramp_off=1000
ramp_tau_start=0.4
ramp_tau_end=0.03
ramp_num_steps=100 */
num_trials=1 //set number of different trials for each combination of pyrspike_tau, baskspike_tau, baskconnvector, and ipspgmax. There will be a different global seed for each simulation
gind_start=0 //global index at which to start, in the parameter sweep
{pyrspiketau_vec = new Vector()}
{baskspiketau_vec = new Vector()}
{baskconnvector = new Vector()}
{ipspg_vec = new Vector()}
{pyrspiketau_vec.append(0.10)} //mean time between arrival of noise events to pyramidal cells (so smaller number implies more intense noise)
{baskspiketau_vec.append(6.0)} // 1e9 orig.: mean time between arrival of noise events to basket cells (so smaller number implies more intense noise)
{baskconnvector.append(0)} //controls what percentage of possible connections from basket cells to pyramidal cells are realized
{ipspg_vec.append(5.5)} //gmax for bask->pyr synapses, which goes into _f# in the filename
//29 was old method
//30 is adding in chris' new pieces one by one. First is having shortspikes and noisesyn. It worked
//31 is having baskets 80-99 and ant above. Working
//32 recording from ant cells better (was same recording until now). Working?
//33 gapped baskets to see ant cell response- gaps worked, but no real response.
//34 Will add bias current to ants, chose 1.5 to make it AP, not enough, made it 3 and now it works.
//35 adding noisesyn to baskets. worked
//36 adding noisesyn to ant, and activateAntSynapses.. worked. Makes normal HFO (high basket gaps, ant/bask/pyr all same noise (.33/ 3.333). Also two chains of pyr gaps
//37 trying to make abnormal HFO. (.2 to all no real change, .1 and basketstart 1e9 looks good)
//38 same as above without antenna's active (just antenna for baskets, not for spike). manymoreconn has one additional gap between 88 and 90 to remove the
// "split". 600 ant cells. baskets and pyr both active.
//39 3000 antenna, had to change the location style
pc = new ParallelContext()
pc.subworlds(1)
func getTstop() { return Tstop }
proc prinit() { //had to change name because non-pr functions couldn't address pr
{pr.setScatteredVoltages(-85, -60)}
{pr.connectNetwork($1,$2)} // took this out of init() in TGbignet2.tem
{pr.setSeed($3)} //set global index for Random123
// { pr.activeSynapsesZero()} //CF: This inactivates all connections. I have no idea why you would want to do this.
finitialize()
finitialize()
}
func onerun() {local id, num, ipspg, pyrthr, basketthr, pyr_spike_tau, bask_spike_tau, bask_perc, temp_time, temp_tau,ii localobj pc, fo, fo1, forast
id= hoc_ac_
pc = new ParallelContext()
{pr = new TGbignet2()}
if(doextra) {
print "Doing EXTRACELLULAR stuff..... ....."
setpointers()
setelec(-50, 0, 0)
}
//hines = new FInitializeHandler(2, "hinest1=startsw() hinest2=startsw() hines1()")
print "id_world number ", pc.id_world, " id_bbs ", pc.id_bbs, " id ", id, " pc.id ", pc.id
//fih_progress = new FInitializeHandler(2, "cvode.event(100, \"progress()\" )" ) //took out if (pc.id == 0)
//forall { for (x,0) print secname() }
{pr.recordVoltages()}
{pr.pnm.set_maxstep(0.01)}
{pr.pnm.want_all_spikes()}
runningTime = startsw()
stdinit()
bask_perc=$3 //percentage of inhibitory connections that are allowed to exist
ipspg=$4
g_ind=$5
prinit(bask_perc,ipspg,g_ind) //added this extra function to allow for non-pr functions. The input is passed to connectNetwork as the connthr
normmean=0 //can set this to $5 if useful in the future
pyr_spike_tau=$1
pyr_nospike_tau=1.0 // 1e9
bask_spike_tau=$2
bask_nospike_tau=6.0 // 1e9
{ pr.activatePyrSynapses(pyr_spike_tau,pyr_nospike_tau) }
{ pr.activateAntSynapses(pyr_spike_tau,pyr_nospike_tau) }
{ pr.activateBaskSynapses(bask_spike_tau,bask_nospike_tau) }
//we no longer need shortspikes, because noisesyn.mod provides stochastic stimulation; do need to pilfer 'addAntennaDC' from shortspikes, though
// if (shortspikeFlag) {
// {pr.shortspikes(Tstop,pyrthr,basketthr, spikedur, spikefreq, normmean, normvar)}
/* } else {
//{pr.shortnonrandomspikes(Tstop,pyrthr,basketthr, spikedur, spikefreq)} } else {
//{pr.singlecellnonrandom(Tstop,pyrthr,basketthr, spikedur, spikefreq)} } else {
{pr.activeSynapsesRandom(Tstop, pyrthr, basketthr)}
} */
{ pr.addAntennaDC(antennaDC) } //CF: add DC current to all antenna cells, so that their resting membrane potential is above the -80mV chloride reversal potential; spike threshold is slightly above 1.3
pc.post(id, pyrthr, ipspg, basketthr)
forast = new File()
sprint(cmd, "data/spikes_b%4.2f_p%4.2f_g%4.2f_f%d.dat", bask_spike_tau, pyr_spike_tau, ipspg,bask_perc)
{ forast.wopen(cmd) }
//now set up all the cvode events which will change the pyr_spike_tau values delivered to the pyramidal cells
/* for ii=1, ramp_num_steps {
temp_time=ramp_on + ii*(ramp_off-ramp_on)/ramp_num_steps //calculate the time at which the iith noise change should occur
temp_tau=ramp_tau_start + ii*(ramp_tau_end-ramp_tau_start)/ramp_num_steps //calculate the value of tau for the iith noise change
sprint(cmd,"pr.activatePyrSynapses(%5.3f,pyr_nospike_tau,normalstd)",temp_tau) //create the command to use 'activatePyrSynapses' to change the value of pyr_spike_tau
cvode.event(temp_time,cmd)
} */
//advance through simulation in increments of t_seg (defined in externals.hoc); after every t_seg, write voltage data to files, and delete vectors containing this data,
//so that program does not run out of memory
t_curr = 0
while (t_curr < Tstop-dt){ //include the '-dt' to account for rounding error; otherwise, may get error in writeVoltages
print "Time = ",t_curr
if(t_curr + t_seg < Tstop) {
{ pr.pnm.pcontinue(t_curr+t_seg)}
} else {
{pr.pnm.pcontinue(Tstop)}
}
for i=0, pr.pnm.spikevec.size-1 {
forast.printf("%-10.6lf, %d\n", pr.pnm.spikevec.x[i], pr.pnm.idvec.x[i])
}
pr.pnm.spikevec.resize(0)
pr.pnm.idvec.resize(0)
pr.writeVoltages(bask_spike_tau, pyr_spike_tau, ipspg, bask_perc,t_curr)
t_curr = t_curr + t_seg
}
{forast.close()}
//basketthr=$2 //reset in case I used it above differently
//if (realrunFlag) {pr.writeVoltages(basketthr, pyrthr, gapstyle, sigfreq)}
runningTime = startsw() - runningTime
iteration=iteration+1
print "Running Time: ", runningTime, "iteration: ",iteration
/* {pr.pnm.gatherspikes()}
fo = new File()
fo1= new File()
sprint(cmd, "data/spikes.dat")
{fo1.wopen(cmd)}
sprint(cmd, "data/spikes_b%d_p%d_g%d_f%d.dat", basketthr, pyrthr, gapstyle,sigfreq)
if (realrunFlag) {fo.wopen(cmd)}
for i=0, pr.pnm.spikevec.size-1 {
if (realrunFlag) fo.printf("%-10.6lf, %d\n", pr.pnm.spikevec.x[i], pr.pnm.idvec.x[i])
fo1.printf("%-10.6lf, %d\n", pr.pnm.spikevec.x[i], pr.pnm.idvec.x[i])
}
if (realrunFlag) {fo.close()}
{fo1.close()} */
//{pr.pnm.pc.done()} //remove this one?
if (realrunFlag){
{fo=new File()}
{sprint(cmd, "spikes_b%4.2f_p%5.3f_g%4.2f_f%d.dat", bask_spike_tau, pyr_spike_tau,ipspg,bask_perc)}
{fo.aopen("data/spikelog.dat")}
{fo.printf("%s\n",cmd)}
{fo.close()}
{fo=new File()}
{sprint(cmd, "sum_b%4.2f_p%5.3f_g%4.2f_f%d.dat", bask_spike_tau, pyr_spike_tau,ipspg,bask_perc)}
{fo.aopen("data/sumlog.dat")}
{fo.printf("%s\n",cmd)}
{fo.close()}
pr.writeParameters(bask_spike_tau, pyr_spike_tau, gapstyle, ipspg, bask_nospike_tau, pyr_nospike_tau, antennaDC, Tstop, t_seg) //write parameters to file
}
{pc.gid_clear()}
{pr.pnm.pc.gid_clear()}
pr=nil
return id
}
proc progress() {
print t
cvode.event(t+100, "progress()" )
}
proc hines1() {
printf("pc.id_world= %d pc.id_bbs= %d t= %g dt= %g dreal=%g treal=%g\n", pc.id_world, pc.id_bbs, t, dt, startsw()-hinest2, startsw()-hinest1)
hinest2 = startsw()
cvode.event(t + 10, "hines1()")
}
{pc.runworker()}
proc series() {local i, j, k, delay, tstop, id, spkcnt, tmax, gid, num
for i = 0, pyrspiketau_vec.size()-1 {
for j = 0, baskspiketau_vec.size()-1 {
for k = 0, baskconnvector.size()-1 {
for l = 0, ipspg_vec.size()-1 {
for m = 0, num_trials-1 {
//generate a different global index for each simulation
g_index = gind_start+i*baskspiketau_vec.size()*baskconnvector.size()*ipspg_vec.size()*num_trials + j*baskconnvector.size()*ipspg_vec.size()*num_trials + k*ipspg_vec.size()*num_trials + l*num_trials + m
{pc.submit("onerun", pyrspiketau_vec.x[i], baskspiketau_vec.x[j], baskconnvector.x[k], ipspg_vec.x[l],g_index)}
}
}
}
}
}
while ((id= pc.working())!=0) {
pc.take(id, &num)
//printf("num= %d", num)
}
}
series() //this line actually runs the simulation
{pc.done()}
quit()
//onerun(45,90,1, 400) //pyrthr, baskthr, f, normvar