// By Jose Francisco Gomez Gonzalez, jfcgomez@ull.es, ULL, Spain
// jfcgomez@gmail.com
// last version: 28 June 2008
load_proc("nrnmainmenu") // load main NEURON library
//cvode.active(0)
//---------------------- xopen_filehoc() ------------------
//Outputs: Load Files.hoc https://correoweb.ccti.ull.es/horde/imp/mailbox.php?mailbox=INBOX&actionID=149
//Inputs: $s1 is a string root to "file.hoc" from the actual root
// $s2 is a string with the file name (not put .hoc)
proc xopen_filehoc() {
strdef tmp_str
sprint(tmp_str,"xopen(\"%s/%s.hoc\")",$s1,$s2)
execute1(tmp_str)
}
//------------------------------------------------------------------------------------------------
xopen_filehoc("../../template","ExperimentControl")
strdef accstr // not confuse experimental variable bindings with neurophysiological variable bindings
objref econ // Create an experiment object
show_errs=1
debug_lev=1
strdef textSegment
strdef cellname, cellfolder
cellfolder="youngcell" //"agedcell" //"youngcell"
cellname="n123" //"n123", "n125", "n128", "n129", "n130" YOUNG NEURONS
//"n170", "n172", "n175", "n178", "n180", "n182" AGED NEURONS
sprint(textSegment,"../../morphology/nXXX/%s",cellfolder)
econ=new ExperimentControl(show_errs,debug_lev)
econ.self_define(econ) // points the object at itself
econ.morphology_dir = textSegment // Setup morphology directory
econ.data_dir = "data" // Define directory to save produced data
sprint(econ.syscmd, "mkdir -p %s", econ.data_dir) // make the data directory
system(econ.syscmd)
//------------------------------------------- Setup cell------------------------------------------
econ.xopen_geometry_dependent(cellname)
cellfolder="nXXX"
sprint(textSegment,"../../morphology/%s",cellfolder)
econ.morphology_dir = textSegment // Setup morphology directory
xopen_filehoc(econ.morphology_dir,"apical-trunk-list")
xopen_filehoc(econ.morphology_dir,"apical-non-trunk-list")
xopen_filehoc(econ.morphology_dir,"soma-list")
xopen_filehoc(econ.morphology_dir,"basal-tree-list")
//------------------------------------------------------------------------------------------------
objref sr,distance_element
objref fileref //to save in file
//------------------------- Calculate the Center of Mass of the soma
proc CenterOfMass(){
xcg=0
ycg=0
zcg=0
ncp=0
Sum_diam=0
strdef temp
forsec "soma" {
print secname()," Area= ",area(0.5)
n=n3d()
nseg=n
ncp+=n
for i=0,n-1 {
print "# of segment",i,"\tDiam",diam3d(i),"\tx=",x3d(i),"\ty=",y3d(i),"\tz=",z3d(i)
xcg+=x3d(i)*diam3d(i)
ycg+=y3d(i)*diam3d(i)
zcg+=z3d(i)*diam3d(i)
Sum_diam+=diam3d(i)
}
}
xcg/=Sum_diam
ycg/=Sum_diam
zcg/=Sum_diam
print "Center of soma: ","Xcg= ", xcg ,"\tYcg= ", ycg,"\tZcg= ", zcg
}
//---------------------------- End Center of Mass
//---------------------------------- myunitvector() ----------------------------
//Inputs: $o1 is the SectionList called apical_trunk_list
// $o2 is the SectionLIst called apical_non_trunk_list
proc myunitvector(){local loop, leght
vector_ux=0
vector_uy=0
vector_uz=0
forsec $o1{
sr=new SectionRef()
lenght=sr.sec.L //I use the lenght as a weight to give more
//importance to longer dendrites
diameter=diam3d(0)
vector_ux+=lenght*diameter*(-xcg+x3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
vector_uy+=lenght*diameter*(-ycg+y3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
vector_uz+=lenght*diameter*(-zcg+z3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
i+=1
}
// forsec $o2{
// sr=new SectionRef()
// lenght=sr.sec.L
// diameter=diam3d(0)
// vector_ux+=lenght*diameter*(-xcg+x3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
// vector_uy+=lenght*diameter*(-ycg+y3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
// vector_uz+=lenght*diameter*(-zcg+z3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
// i+=1
// }
// forsec $o3{
// sr=new SectionRef()
// lenght=sr.sec.L
// diameter=diam3d(0)
// vector_ux+=lenght*diameter*(-xcg+x3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
// vector_uy+=lenght*diameter*(-ycg+y3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
// vector_uz+=lenght*diameter*(-zcg+z3d(0))/sqrt((xcg-x3d(0))^2+(ycg-y3d(0))^2+(zcg-z3d(0))^2)
// i+=1
// }
modulevector=sqrt((vector_ux)^2+(vector_uy)^2+(vector_uz)^2)
vector_ux=vector_ux/modulevector
vector_uy=vector_uy/modulevector
vector_uz=vector_uz/modulevector
print "FINAL vector normal = ",(vector_ux)^2+(vector_uy)^2+(vector_uz)^2
print "FINAL vector coord = ",vector_ux,vector_uy,vector_uz
}
//------------------------------END myunitvector()--------------------------------------------------
//---------------------------------- mydistance() ----------------------------
//Inputs: $o1 is the SectionList called apical_non_trunk_list
proc mydistance(){local loop, leght
i=0
forsec $o1 {i+=1} //to save in file
distance_element=new Vector(i)
fileref=new File() //to save in file
sprint(temp, "%s/distance%snewvector.dat", econ.data_dir,cellname) //to save in file
fileref.wopen(temp) //to save in file
fileref.printf("Apical Dendrite Dist.soma Lenght Diam.\n")
i=0
forsec $o1{
sr=new SectionRef()
lenght=sr.sec.L
distance_element.x[i]=(x3d(0)-xcg)*vector_ux+(y3d(0)-ycg)*vector_uy+(z3d(0)-zcg)*vector_uz
temp=secname()
fileref.printf( "%s %g %g %g\n", temp,distance_element.x[i],lenght,diam)
printf( " %s distance: %g Lenght: %g diam: %g\n", temp,distance_element.x[i],lenght,diam)
i+=1
}
fileref.close() //to save in file
}
//------------------------------END mydistance() ----------------------------
// -- Run the next functions ----
CenterOfMass()
//myunitvector(apical_trunk_list,apical_non_trunk_list,basal_tree_list)
myunitvector(apical_trunk_list,apical_non_trunk_list,soma_list)
mydistance(soma_list)
mydistance(apical_trunk_list)
mydistance(apical_non_trunk_list)
//------------------------------------------ Display NEURON--------------
objref display
display=new Shape()
display.show(0) // mode=0 diameters, =1 centroid, =2 schematic
create center_soma
center_soma{pt3dclear()
pt3dadd(xcg-1, ycg, zcg, 1) //-1 and + 1 are added to can see center_soma at display (shape)
pt3dadd(xcg, ycg, zcg, 1)
pt3dadd(xcg+1, ycg, zcg, 1)
}
access center_soma
print "Center of soma: ","Xcg= ", x3d(1) ,"\tYcg= ", y3d(1),"\tZcg= ", z3d(1)
create vector_cell
vector_cell{pt3dclear()
pt3dadd(xcg, ycg, zcg, 1) //-1 and + 1 are added to can see center_soma at display (shape)
pt3dadd(vector_ux*500+xcg, vector_uy*500+ycg, vector_uz*500+zcg, 1)
}
access center_soma
print "Center of soma: ","Xcg= ", x3d(1) ,"\tYcg= ", y3d(1),"\tZcg= ", z3d(1)
//center_list is to be able to display the center_soma
objref center_list
center_list=new SectionList()
center_soma center_list.append()
//vector_list is to be able to display the unit vector
objref vector_list
vector_list=new SectionList()
vector_cell vector_list.append()
//Origin_list is to be able to display the center_soma
create Origin
Origin{pt3dclear()
pt3dadd(0-1, 0, 0, 1) //-1 and + 1 are added to can see center_soma at display (shape)
pt3dadd(0, 0, 0, 1)
pt3dadd(0+1, 0, 0, 1)
}
objref Origin_list
Origin_list=new SectionList()
Origin center_list.append()
display.color_list(Origin_list,1)
display.color_list(center_list,1)
display.color_list(vector_list,6)
display.color_list(apical_trunk_list,6)
// display.color_list(axon_sec_list,3)
display.color_list(soma_list,2)
// display.color_list(dendrite_sec_list,4)
// display.color_list(apical_dendrite_sec_list,5)