// This hoc file reproduces the simulation of a border cell in a circular environment from Fig. 9C of Welday et al.
// The border cell is formed by a target neuron that receives inhibitory input from twelve theta cells.  
//
// Before running this simulation, the theta cell spike trains from which the border cell is formed must first be created
// and saved in files. This is done by running the MATLAB script 'boundary_thetaspikes.m' from within the simulation directory.
// Note that since the spike trains are generated stochastically, the MATLAB script will not generate identical 
// spike trains every time it is run (unless the random number generator is identically seeded).
// 
// Since we are modeling a 1 hour recording session (3600 seconds), the simulation runs for a long time. To maximize speed 
// and performance, you may wish to close the graph windows that plot the voltage trace (Graph[2]) and the spike raster
// (SpikePlot[0] for NetData[0]) before hitting the 'Init & Run' button to start the simulation.
//
// To save the place cell's spike times when the simulation finishes, type:
//
//   >>load_file("savespikes.hoc")
//
// at the interpreter prompt after the simulation has finished running. Spikes will be saved in a file called 'SPIKOUT.dat',  
// and you can then create a path plot of the simulated grid cell's spike output by running the MATLAB script 'plotnrn.m'.
//

load_file("nrngui.hoc")

// -----------------------------------------------------------------
// read theta cell spike trains from disk files into NEURON vectors
// -----------------------------------------------------------------

load_file("boundary_vecstims.ses") 	//create vecstim objects for delivering theta spike trains to the model neuron

numinputs = 12				//number of inhibitory theta inputs to the model neuron

objref evec[numinputs] 			//array of event vectors (VecStim.mod) for reading theta spike trains
objref spikefile       			//file object through which to read in theta cell spike times from disk
objref pplist          			//list of point processes containing theta spike trains
objref vmvec
objref vmfile

spikefile = new File()
pplist = new List()

vmfile = new File()


load_file("read_boundary_theta_spikes.hoc")  //read in theta cell spike timestamps from files
simlength=1000*3600


// --------------------------------------------------------------
// build a single-compartment postsynaptic cell
// --------------------------------------------------------------

ra        = 150 	// axial resistance through cytoplasm (ohms)
global_ra = ra		
rm        = 15000 	// passive membrane resistance (ohms) 
c_m       = 1   	// membrane capacitance (microFarads per centimeter squared, uF/cm^2)

create soma		//single somatic compartment
access soma

PI=3.14159
{L=10/PI  diam=10}
nseg=1
Ra = ra
cm = c_m 

vmvec = new Vector()
vmvec.record(&v(.5))

//Insert voltage-activated persistent sodium current (Nap.mod)
    insert nap
    gbar_nap=0.00005	//Nap conductance
    sh_nap=-16		//Nap voltage activation threshold shift parameter (mV)

//Insert Hodkin-Huxley kinetics (hh.mod, standard NEURON mechanism)
    insert hh
    gkbar_hh=0.005	//delayed rectifier K+ conductance
    gnabar_hh=0.05	//voltage-gated Na+ conductance
    el_hh=-65		//leak reversal potential 
    gl_hh=1/rm		//leak conductance

// --------------------------------------------------------------
// connect theta spike train inputs to the model neuron
// --------------------------------------------------------------

objref nclist          			//list of netcon objects for synaptic connections
nclist = new List()

load_file("place_GABA_inputs.ses") //GABA-A synapses are simulated by modifying an AMPA current (found in ampa.mod)
Erev_AMPA_S=-80			//set reversal potential to -80 mV to convert the AMPA synapse to GABA
synw=.001			//conductance of GABA synapses

for i=0,numinputs-1 {

nclist.append(new NetCon(pplist.o(i), AMPA_S[i], -20, 1, synw)) //make the input connections

}



tstop=simlength

// ------------------------------------------------------------------------
// store model neuron's spike times in a vector to be saved later if needed
// ------------------------------------------------------------------------

objref spiketimes	//vector in which to store timestamps of spikes generated by the model neuron
objref spikenc		//netcon object through which spikes are passed into the 'spiketimes' vector
objref null		//null object

spiketimes = new Vector()
spikenc = new NetCon(&v(0.5), null)
spikenc.threshold = -25    // store a spike timestamp if the postsynaptic membrane voltage exceeds -20 mV
spikenc.record(spiketimes) // record the spike times to the 'spiketimes' vector

//To save the postsynaptic neuron's spike times, type:
//
//   >>load_file("savespikes.hoc")
//
//at the interpreter prompt after the simulation has finished running.
//Spikes will be saved in a file called 'SPIKOUT.dat'.  MATLAB code 
//(plotnrn.m) is provided for generating a "path plot" of the model 
//neuron's spatial tuning function from the 'SPIKEOUT.dat' file.