load_file("stdlib.hoc")
load_file("ObliquePath.hoc")
tstop = 1500
steps_per_ms = 40
dt = 0.025
stamp=-0.1
stdur=1000
stdel=20
objref fchirp, ChirpStim
tstop = 25050
steps_per_ms = 40
dt = 0.025
ChirpAmp=0.05
rasoma = 50
raend = 35
raaxon = 50
rm = 200000
rmsoma = rm
rmdend = rm
rmend=12000
rmaxon = rm
c_m = 1
cmsoma = c_m
cmdend = c_m
cmaxon = c_m
v_init = -65
celsius = 34
Eh=-30
gh=0.0006
gh=3.4e-5
create axon[2]
objref sh, st, apc, vec, g
strdef str1, str2, str3
objref SD, AXON, SA, Basal, Trunk
objref s, ampasyn[1], nmdasyn[1]
objref bpropampasyn, bpropnmdasyn
objref LM, RAD, RADt, LUC, PSA, PSB, ORI
objref secref, f1
objref netlist
objref apamp[1], ampvec[1]
objref sapamp, sampvec
objref g1, nc
objref r, st1
objref Trial
objref local_trunk
objref pl[150],opl[150]
create soma[1], apical[1], basal[1]
objref v3
proc setpassive() {
forall {
insert pas
e_pas = v_init
Ra=rasoma
}
forsec SD {
cm=cmdend
g_pas=1/rmdend
}
forsec "soma" {
cm = cmsoma
g_pas=1/rmsoma
}
forsec Trunk {
for (x) {
rdist=raddist(x)
rm = rmsoma + (rmend - rmsoma)/(1.0 + exp((200-rdist)/50))
g_pas(x)=1/rm
Ra = rasoma + (raend - rasoma)/(1.0 + exp((210-rdist)/50))
}
}
for i=0,plcount {
seccount=0
forsec pl[i] {
if(!seccount){
trunk_pas=g_pas(1)
seccount=seccount+1
} else {
g_pas=trunk_pas
seccount=seccount+1
}
}
}
}
proc setactive () {local xdist, ndist
forsec Trunk {
insert hd
for (x) {
xdist=raddist(x)
ghdbar_hd(x) = gh*(1+100/(1+exp(-(xdist-380)/34.0)))
if (xdist > 100){
if (xdist>300) {
ndist=300
} else {
ndist=xdist
}
vhalfl_hd(x)=-82-8*(ndist-100)/200
} else {
vhalfl_hd(x)=-82
}
}
}
for i=0,plcount {
seccount=0
forsec pl[i] {
if(!seccount){
trunk_h=ghdbar_hd(1)
trunk_vhalf=vhalfl_hd(1)
seccount=seccount+1
} else {
insert hd
ghdbar_hd=trunk_h
vhalfl_hd=trunk_vhalf
seccount=seccount+1
}
}
}
forsec "soma" {
insert hd ghdbar_hd=gh vhalfl_hd=-82
}
forsec Basal {
insert hd ghdbar_hd=gh vhalfl_hd=-82
}
forall if (ismembrane("hd") ) ehd_hd=Eh
}
proc init_cell() {
setpassive()
addaxon()
setactive()
access soma[2]
distance()
print "Number of compartments: ", totcomp
finitialize(v_init)
fcurrent()
forall {
for (x) {
if (ismembrane("hd")) {
e_pas(x)=v(x)+i_hd(x)/g_pas(x)
} else {
e_pas(x)=v(x)
}
}
}
}
proc current_clamp(){
access apical[65]
st1=new IClamp(.5)
st1.dur = stdur
st1.del = stdel
st1.amp = stamp
access soma[0]
}
proc load_3dcell() {
forall delete_section()
xopen($s1)
SD = new SectionList()
SA = new SectionList()
Trunk = new SectionList()
Trial = new SectionList()
Basal = new SectionList()
forsec "soma" {
SD.append()
SA.append()
}
forsec "basal" {
SD.append()
Basal.append()
}
forsec "apical"{
SD.append()
SA.append()
}
soma[0] Trunk.append()
apical[0] Trunk.append()
apical[4] Trunk.append()
apical[6] Trunk.append()
apical[14] Trunk.append()
apical[15] Trunk.append()
apical[16] Trunk.append()
apical[22] Trunk.append()
apical[23] Trunk.append()
apical[25] Trunk.append()
apical[26] Trunk.append()
apical[27] Trunk.append()
apical[41] Trunk.append()
apical[42] Trunk.append()
apical[46] Trunk.append()
apical[48] Trunk.append()
apical[56] Trunk.append()
apical[58] Trunk.append()
apical[60] Trunk.append()
apical[62] Trunk.append()
apical[64] Trunk.append()
apical[65] Trunk.append()
apical[69] Trunk.append()
apical[71] Trunk.append()
apical[81] Trunk.append()
apical[83] Trunk.append()
apical[95] Trunk.append()
apical[103] Trunk.append()
apical[104] Trunk.append()
load_file("oblique-paths.hoc")
setpassive()
totcomp=0
forall {
soma[1] area(0.5)
nseg = int((L/(0.1*lambda_f(100))+.9)/2)*2 + 1
totcomp=totcomp+nseg
}
access soma[2]
distance()
init_cell()
}
proc addaxon() {
AXON = new SectionList()
for i = 30,34 basal[i] {
AXON.append()
Basal.remove()
}
for i = 18,22 basal[i] {
AXON.append()
Basal.remove()
}
forsec AXON {
e_pas=v_init
g_pas = 1/rmaxon
Ra=raaxon
cm=cmaxon
}
}
proc initsub() {
load_3dcell("n123.hoc")
finitialize(v_init)
fcurrent()
}
proc update_Gh() { local xdist
forsec "apical" {
if (ismembrane("hd") ) {
for (x) {
xdist=raddist(x)
ghdbar_hd(x)=gh*(1+2*xdist/100.0)
}
}
}
forsec "soma" {
ghdbar_hd=gh
}
forsec Basal {
ghdbar_hd=gh
}
update_init()
}
proc update_init(){
finitialize(v_init)
fcurrent()
forall {
for (x) {
if (ismembrane("hd")) {
e_pas(x)=v(x)+i_hd(x)/g_pas(x)
} else {
e_pas(x)=v(x)
}
}
}
}
somax=2.497
somay=-13.006
somaz=11.12
double distances[200]
func raddist() {
distn0=distance(0)
distances[0]=0
sum=0
for i=1,n3d()-1 {
xx=(x3d(i)-x3d(i-1))*(x3d(i)-x3d(i-1))
yy=(y3d(i)-y3d(i-1))*(y3d(i)-y3d(i-1))
zz=(z3d(i)-z3d(i-1))*(z3d(i)-z3d(i-1))
sum=sum+sqrt(xx+yy+zz)
distances[i]=sum
}
xval=$1
distn=distance(xval)
match=distn-distn0
matchptdist=100000
for i=0,n3d()-1 {
matptdist=(match-distances[i])*(match-distances[i])
if(matchptdist>matptdist){
matchptdist=matptdist
matchi=i
}
}
xx=(x3d(matchi)-somax)*(x3d(matchi)-somax)
yy=(y3d(matchi)-somay)*(y3d(matchi)-somay)
zz=(z3d(matchi)-somaz)*(z3d(matchi)-somaz)
return sqrt(xx+yy+zz)
}
proc GenerateChirp() {
ChirpStim=new Vector(25050/dt)
ChirpStim.fill(0)
for i=50/dt, ChirpStim.size()-1 {
xval=(i-50/dt)*dt/1000
ChirpStim.x[i]=ChirpAmp*sin(2*3.141592654*25*(i-50/dt)/(2*(25000/dt)) *xval)
}
fchirp=new File()
fchirp.wopen("ChirpI.txt")
ChirpStim.printf(fchirp,"%f\n")
}
proc Chirp_Trunk() {
print "\n\nComputing Trunk Compartment Responses to Local Chirp Stimuli.."
print "\n\nOutput will be saved in the Output directory..."
print "\n\nNote: These simulations take time, as each location requires 25050 ms of simulation data!"
initsub()
GenerateChirp()
print "Chirp Generated....\n"
count=0
update_init()
forsec Trunk {
for(x) {
if(x != 0) {
print count, " ", secname(), " ", x, raddist(x)
st=new IClamp(x)
st.dur=tstop
st.del=0
v3=new Vector()
v3.record(&v(x))
ChirpStim.play(&st.amp,dt)
f1=new File()
finitialize(v_init)
fcurrent()
while (t < tstop){
fadvance()
}
sprint(str1,"Output/ChirpOutput_%d.txt",count)
f1.wopen(str1)
v3.printf(f1,"%f\n")
f1.close()
count=count+1
}
}
}
}
objref st
objref f1
proc RN_Trunk () {
print "\n\nComputing Input Resistances Along Trunk..."
print "\n\nOutput will be saved in Trunk_RN.txt..."
initsub()
v3=new Vector()
v3.record(&soma.v(0.5))
f1=new File()
tstop=250
f1.wopen("Trunk_RN.txt")
update_init()
count=0
forsec Trunk {
for(x) {
if(x != 0) {
st=new IClamp(0.5)
st.dur=tstop
st.del=0
st.amp=-0.1
v3=new Vector()
v3.record(&v(x))
finitialize(v_init)
fcurrent()
while (t < tstop){
fadvance()
}
f1.printf("%f\t%f\n",raddist(x),-((v3.x[tstop/dt-1]-v_init)/100e-6)*1e-3)
print secname(), " ", x," ", raddist(x)," ", -((v3.x[tstop/dt-1]-v_init)/100e-6)*1e-3
}
}
count=count+1
}
f1.close()
}
xpanel("Chirp Stimulus Responses in a 3D model with h Channels")
xlabel("Select one of the simulations below")
xbutton("Compute Input Resistances Along Trunk", "RN_Trunk()")
xbutton("Save Local Chirp Responses for Locations Along Trunk", "Chirp_Trunk()")
xpanel()