MembraneUtilities.hoc

// proc MechList( List )	-	adds String with the name of all membrane mechanisms in cas to List
	// dependency: none

// proc getglist( List )	-	get the names of all membrane mechanisms in cas with a 'g' variable
	// dependency: MechList (local)

// func Rm( channel name(s), Rm/Gm flag )	-	returns membrane resistance or conductance
	// dependency: MakeStringList (local), getglist (local)

// func meanRm( sections, mechanism(s), Rm/Gm flag)
	// dependency: Rm (local), MakeStringList (local), getglist (local)

// proc HCNDist( strdef gradient name )	-	set the gradient of HCN distributions
	// dependency: issplit, StopPar, StartPar (parinit); meanRm (local)


load_file("chanAnalysis.hoc")


proc MechList() { local i, cas	localobj ml, mt, strobj
//MechList( string list, cas flag )
// creates a list of strings with the names of all membrane mechanisms in cas or all sections

	if (object_id($o1,1) == -1)	$o1 = new List()	// if object is null make a new List
	ml = $o1
	mt = new MechanismType(0)
	if (numarg()>1) cas = $2	else cas = 0

	for i=0,mt.count()-1 {
		strobj = new String()
		mt.select(i)
		mt.selected(strobj.s)
		if (cas) {
			if (MechChk(strobj.s) == 1)	ml.append(strobj)
		} else	ml.append(strobj)
		
	}
}

proc PPList() { local i	localobj ml, mt, strobj
//PPList( string list, type )
// creates a list of strings with the names of point processes in cas

	if (object_id($o1,1) == -1)	$o1 = new List()	// if object is null make a new List
	ml = $o1
	mt = new MechanismType(1)
// 	if (numarg()>1) cas = $2	else cas = 0

	strobj = new String()
	for i=0,mt.count()-1 {
		mt.select(i)
		mt.selected(strobj.s)
		if (MechChk(strobj.s) == 1)	ml.append(strobj)
	}
}

func MechChk() { local i,val localobj ms, strobj
//bool = MechChk( string )
// checks whether membrane mechanism named "string" is in current section. Returns 1 if yes, 0 if not

	ms = new MechanismStandard($s1,0)
	val = -1.314e20		// assuming no parameter will be equal to this value
	
	strobj = new String()
	if (ms.count() > 0) {
		ms.name(strobj.s,0)
		ms.set(strobj.s,val)
		ms.in()
		i = (val != ms.get(strobj.s))
	}
	return i
}


proc getglist() { local i,j,flag, cas localobj glist, ms, SF
// getglist( List, cas flag )
// generate a list of conductance mechanisms in cas or all sections (all membrane mechs with a 'g' value)

	if (numarg()>1) cas = $2	else cas = 0
	if (object_id($o1,1) == -1)	$o1 = new List()	// if object is null make a new List
	glist = $o1
	MechList(glist, cas)
	
	SF = new StringFunctions()

	for (i=glist.count()-1; i>=0; i=i-1) {
		ms = new MechanismStandard(glist.o(i).s, 0)
		flag=0
		for j=0,ms.count()-1 {
			ms.name(strtmp,j)
			if (SF.substr(strtmp,"g_") != -1) {
				flag=1
			}
		}
		if (flag == 0) glist.remove(i)
	}
}


proc getgvec() { local i,j localobj gvec, ml, ms, SF
// getgvec( Vector, List )
// getgvec( Vector )
// generate vector of conductance values of all mechanisms in List. If no List given, calls
//		MechList and uses all conductance mechanisms

	gvec = $o1
	if (numarg()>1) {
		ml = $o2
	} else {
		ml = new List()
		MechList(ml)
	}
	if (object_id(gvec,1) == -1) gvec=new Vector(ml.count(),0)
	
	SF = new StringFunctions()

	for i=0,ml.count()-1 {
		ms = new MechanismStandard(ml.o(i).s, 0)
		for j=0,ms.count()-1 {
			ms.name(strtmp,j)
			if (SF.substr(strtmp,"g_") != -1) {
				gvec.x[i]=1
			}
		}
	}
}


proc MakePassive() {local i,Gsum,argn,vi	localobj sl, gls

	argn = numarg()
	vi = v_init
	
	if (verbosity > 3) for i=1,argn print argtype(i)

	if (argn > 1) {
		vi = $2
	} else if (argn > 0) {
		if (argtype(1)==2)	{
			sl = new SectionList()
			forall ifsec $s1 sl.append()
			if (verbosity > 4) sl.printnames()
		} else if (argtype(1)==1)	sl = $o1
	} else {
		sl = new SectionList()
		sl.append()
	}
	
	init(v_init)

	forsec sl {
		gls = new List()
		getglist(gls,1)

		for (x,0) {
			RmGtot=0
			for i=0,gls.count()-1 {
				if (ismembrane(gls.o(i).s)) {
					sprint( strtmp, "RmGtot += g_%s", gls.o(i).s)
					execute(strtmp)
					sprint( strtmp, "uninsert %s", gls.o(i).s)
					execute(strtmp)
				}
			}
			insert pas
			g_pas = RmGtot
			e_pas = vi
		}
	}

}

func Rm() { local i,flag localobj gls
// val = Rm( strdef channel name )
// val = Rm( List channel name(s) )
// val = Rm( channel name(s), Rm/Gm flag )
// val = Rm( 0, Rm/Gm flag )
// measure membrane resistivity megaohm*cm^2 (flag == 0) or conductivity S/cm^2 (flag == 1)

	flag = 0
	RmGtot=0

	if (numarg() > 0) {
		if (argtype(1) == 1) {
			gls = $o1
		} else if (argtype(1)==2) {
			gls = new List()
			MakeStringList(gls,$s1)
		} else if (argtype(1)==0) {
			gls = new List()
			getglist(gls)
		}
	} else {
		gls = new List()
		getglist(gls)
	}
	
	if (numarg() > 1) flag = $2
	
	for i=0,gls.count()-1 {
		if (ismembrane(gls.o(i).s)) {
			sprint( strtmp, "RmGtot += g_%s", gls.o(i).s)
			//if (verbosity > 3) printf("%s", strtmp)
			execute(strtmp)
		}
	}
	if (flag ==1) {
		return RmGtot
	} else {
		if (RmGtot == 0) return 1e30 else return 1e-3/RmGtot
	}
}


func sum_state() { local err
// sum = sum_state( SectionList, strdef )

	Ssum=0
	err=-1
	sprint( strtmp, "Ssum += %s*area(0.5)", $s2)
	forsec $o1 {
		if (err==-1) {
			err = execute1(strtmp)
			if (err==0) return -1
		} else execute(strtmp)
	}
	return Ssum
}

func meanRm() { local argn,i,flag, Asum,Gsum localobj gls, sl
// val = meanRm( sections, mechanism(s), Rm/Gm flag)
// 'sections' can be a SectionList or string specifying sections to measure
// 'mechanisms' can be a List of membrane mechanisms, a string with name of mechanism, or 
// 		else a list of all membrane mechanisms is generated
// if Rm/Gm flag == 1 returms mean membrane conductivity (S/cm^2).
// measure membrane resistivity megaohm*cm^2 (flag == 0) or conductivity S/cm^2 (flag == 1)

	flag = 0
	Gsum=0
	Asum = 0
	argn = numarg()
	
	if (verbosity > 3) for i=1,argn print argtype(i)

	if (argn >= 1) {
		if (argtype(1)==2)	{
			sl = new SectionList()
			forall ifsec $s1 sl.append()
			if (verbosity > 4) sl.printnames()
		} else if (argtype(1)==1)	sl = $o1
	} else {
		if (verbosity > 1) printf("Must give at least 1 input arg: meanRm\n")
		return -1
	}

	if (argn >= 2) {
		if (argtype(2)==1) {
			gls = $o2
		} else if (argtype(2)==2) {
			gls = new List()
			MakeStringList(gls,$s2)
		} else if (argtype(2)==0) {
			gls = new List()
			getglist(gls)
		}
	} else {
		gls = new List()
		getglist(gls)
	}
	
	if ((argn >= 3) && (argtype(3)==0) ) flag = $3

	if ( SectionListCount(sl) == 0 )	return -1
	
	if (object_id(sl) != 0) {
		forsec sl for (x,0) {
			Asum += area(x)				// µm2
			Gsum += Rm(gls,1)*area(x)	// S/cm2*µm2 = S*10^-8
		}
	} else {
		forall for (x,0) {
			Asum += area(x)				// µm2
			Gsum += Rm(gls,1)*area(x)	// S/cm2*µm2 = S*10^-8
		}
	}

	if (flag ==1) {
		return Gsum/Asum		// S/cm2
	} else {
		if (Gsum == 0) return 1e30 else return Asum*1e-3/Gsum	// MΩ•cm2
	}
}


obfunc getchan() { local i,flag localobj glist, ms, SF
// mech = getchan( strdef )
// generate a list of conductance mechanisms in section (all membrane mechs with a 'g' value)

	glist = new List()	// if object is null make a new List
	MechList(glist)
	
	SF = new StringFunctions()

	for (i=glist.count()-1; i>=0; i=i-1) {
		if (SF.substr( glist.o(i).s, $s1) != -1) {
			ms = new MechanismStandard( glist.o(i).s, 0)
			flag=1
		}
	}
	return ms
}


func meanCaI() { local argn,i, A,ci localobj sl
	//meanCaI( sections )

	ci= 0
	A = 0
	argn = numarg()
	
	if (verbosity > 3) for i=1,argn print argtype(i)

	if (argn >= 1) {
		if (argtype(1)==2)	{
			sl = new SectionList()
			forall ifsec $s1 sl.append()
			if (verbosity > 2) sl.printnames()
		} else if (argtype(1)==1)	sl = $o1
	}
	
	if (object_id(sl) != 0) {
		forsec sl for (x,0) {
			if (ismembrane("ca_ion")) {
				A += area(x)
				ci += cai*area(x)
			}
		}
	} else {
		forall for (x,0) {
			if (ismembrane("ca_ion")) {
				A += area(x)
				ci += cai*area(x)
			}
		}
	}

	return ci/A
}

func ri2() {local val,i	localobj sref
// val = ri2( void )
// measures the sum axial resistance (MΩ) between cas and its adjacent section

	val=0
	
	sref = new SectionRef()
	if (sref.has_parent()) val=ri(0.5)
	for i= 0,sref.nchild-1 {
		sref.child[i] val+=ri(0.5)
	}
}

func lambda() { local val, rm

	if (numarg() > 0) rm = $1 else rm = Rm()
	val = 0.5*(diam*rm/Ra)^0.5

	return val
}

proc swapchans() {local g
	// swapchans( channel name1, channel name2 )
	// swap chans replaces one conductance with another maintaining gmax in each section
	
	forall {
		strtmp = $s1
		if (ismembrane(strtmp)) {
			sprint(tmpstr, "g = gmax_%s", strtmp)
			execute(tmpstr)
			sprint(tmpstr, "uninsert %s", strtmp)
			execute(tmpstr)
			strtmp = $s2
			sprint(tmpstr, "insert %s", strtmp)
			execute(tmpstr)
			sprint(tmpstr, "gmax_%s = g", strtmp)
			execute(tmpstr)
		}
	}

}

proc ZDirect_g() { local argn,zf,grel,reps,pwr,g,zo,i,max,glim	 localobj z1, z2, chan, sl
// ZDirect_g( section list, ref section, chan, grel, frequency, reps,pwr,glim,max change)

zf=0
grel=1e-5
reps=0
pwr=1
max=1e2
glim=1

finitialize(v_init)

argn=numarg()
chan = new String()
if (argtype(1)==2)	{
	sl = new SectionList()
	forall ifsec $s1 sl.append()
	if (verbosity > 2) sl.printnames()
} else if (argtype(1)==1)	sl = $o1

if (argn>1) {
	if (argtype(2)==2)	{
		sprint(tmpstr, "%s sref = new SectionRef()", $s2)
		if (verbosity > 2) printf("%s\n", tmpstr)
		execute(tmpstr)
	}
	if (argtype(2)==1)	sref = $o2
} else sref = new SectionRef()
if (argn>2) {
	if (argtype(3)==2)	chan.s = $s3
	if (argtype(3)==1)	chan = $o3
}
if (argn>3) grel =$4
if (argn>4) zf = $5
if (argn>5) reps =$6
if (argn>6) pwr = $7
if (argn>7) glim =$8
if (argn>8) max = $9

for i=0,reps {
	z1 = new Impedance()
	z2 = new Impedance()
	sref.sec { z1.loc(0.5) }
	z1.compute(zf)
	forsec sl {
		if (ismembrane(chan.s)) {
			z2.loc(0.5)
			z2.compute(zf)
			
			sref.sec zo = z2.ratio(0.5)
			if (grel == 0) {
				sprint(tmpstr, "val = gmax_%s", chan.s)
				execute(tmpstr)
				g = val+max*val*(zo-z1.ratio(0.5)^pwr)
			} else g = grel+max*grel*(zo-z1.ratio(0.5)^pwr)
			if (verbosity > 3) printf("%g", g)
			if (g<1e-9) g=0
			if (g>glim) g=glim
			sprint(tmpstr, "gmax_%s = %g", chan.s, g)
			execute(tmpstr)
		}
	}
}
}

proc Zratio_g() { local zf,grel,up,pwr,max,g localobj zz, chan, sl
// Zratio_g( section list, init section, chan, grel, frequency, increase/decrease, pwr, max )
// 
// sets conductance within sections based on the ratio of impedance between that section 
// and a reference location

zf=100
grel=1e-5
up=1
pwr=1
max=1

finitialize(v_init)

chan = new String()
if (argtype(1)==2)	{
	sl = new SectionList()
	forall ifsec $s1 sl.append()
	if (verbosity > 3) sl.printnames()
} else if (argtype(1)==1)	sl = $o1

if (numarg()>1) {
	if (argtype(2)==2)	{
		sprint(tmpstr, "%s sref = new SectionRef()", $s2)
		if (verbosity > 2) printf("%s\n", tmpstr)
		execute(tmpstr)
	}
	if (argtype(2)==1)	sref = $o2
} else sref = new SectionRef()
if (numarg()>2) {
	if (argtype(3)==2)	chan.s = $s3
	if (argtype(3)==1)	chan = $o3
	if (verbosity > 2) printf("%s\n", chan.s)
}
if (numarg()>3) grel = $4
if (numarg()>4) zf = $5
if (numarg()>5) up = $6
if (numarg()>6) pwr = $7
if (numarg()>7) max = $8

	//if (msplit) stopPar()
	
	if (verbosity > 2) printf("grel = %g\tzf = %g\tup = %g\tpwr = %g\tmax = %g\n", grel,zf,up,pwr,max)

	zz = new Impedance()
	sref.sec { zz.loc(0.5) }
	zz.compute(zf)
	forsec sl {
		if (ismembrane(chan.s)) {
			
			if (grel == 0) {
				if (strcmp(chan.s, "pas") ==0) { sprint(tmpstr, "val = g_%s", chan.s)
				} else sprint(tmpstr, "val = gmax_%s", chan.s)
				execute(tmpstr)
			} else val = grel

			if (up==0) {
				g = val*zz.ratio(0.5)^pwr
				if (g>max) g=max
				if (strcmp(chan.s, "pas") ==0) { sprint(tmpstr, "g_pas = %g", g)
				} else sprint(tmpstr, "gmax_%s = %g", chan.s, g)
			}
			if (up>0) {
				g = val/zz.ratio(0.5)^pwr
				if (g>max) g=max
				if (strcmp(chan.s, "pas") ==0) { sprint(tmpstr, "g_pas = %g", g)
				} else sprint(tmpstr, "gmax_%s = %g", chan.s, g)
			}
			execute(tmpstr)
		}
	}
}


proc handle_M_dist() { local D, Mhalf, MS, Mmax, Mmin, ms
        Mmax = SIZM
        Mmin = SIZM/10
        Mhalf = 150
        MS = -75
        
        ms = issplit()
        if (ms) stopPar()
        
        soma[0] {distance()}
        
        forsec "Handle" {
        	D = distance(0.5)
        	gmax_M = (Mmin + (Mmax-Mmin)/(1+exp((Mhalf-D)/MS)))
        }

        if (ms) startPar()
}

proc handle_Na_dist() { local D, d_half, Namax, S, Namin, ms
        Namin = SIZNa/50
        Namax = SIZNa
        d_half = 135
        S = -135
        
        ms = issplit()
        if (ms) stopPar()

		soma[0] {distance()}
		
		forsec "Handle" {
        	D = distance(0.5)
        	gmax_Na = (Namin + (Namax-Namin)/(1+exp((d_half-D)/S)))
        }
        if (ms) startPar()

}


proc gvecSetup() {local gi	localobj objtmp, gvec
// gvecSetup(gvec list, conductance/current flag)
// 1st input is list object that the recorded vectors will be added to
// 2nd input specifies whether to record conductance (0) or current (1)

	CleanStepEvent()
	
	if (numarg()>1) gi=$2 else gi=0
	idc_G_[1] = new IClamp()
	Tines[1] idc_G_[1].loc(0.5)
	idc_G_[1].amp=0

	if (numarg() > 0) {
		if (argtype(1)==1)	gvec = $o1 else gvec = gvec_G_
	} else gvec = gvec_G_

	if (SectionListCount(HList)>0) {
		if (strcmp(Hname,"hcn") == 0) {
			Hchan = new hcnAnalysis()
		} else if (strcmp(Hname,"h") == 0) {
			Hchan = new hAnalysis()
		}
		objtmp = new Vector()
		if (gi==0) {
			objtmp.record(idc_G_[1], &Hchan.nh, RecDt)
			objtmp.label("nh")
		} else if (gi==1) {
			objtmp.record(idc_G_[1], &Hchan.ih, RecDt)
			objtmp.label("Ih")
		}
		gvec.append(objtmp)
	}
	
	if (SectionListCount(MList)>0) {
		Mchan = new MAnalysis()
		objtmp = new Vector()
		if (gi==0) {
			objtmp.record(idc_G_[1], &Mchan.nM, RecDt)
			objtmp.label("nM")
		} else if (gi==1) {
			objtmp.record(idc_G_[1], &Mchan.iM, RecDt)
			objtmp.label("IM")
		}
		gvec.append(objtmp)
	}
	
	if (SectionListCount(KDList)>0) {
		if (strcmp(KDname,"KD") == 0) {
			KDchan = new KDAnalysis()
		} else if (strcmp(KDname,"KD_ca3") == 0) {
			KDchan = new KD3Analysis()
		}
		objtmp = new Vector()
		if (gi==0) {
			objtmp.record(idc_G_[1], &KDchan.nKD, RecDt)
			objtmp.label("nKD")
			gvec.append(objtmp)
			objtmp = new Vector()
			objtmp.record(idc_G_[1], &KDchan.lKD, RecDt)
			objtmp.label("lKD")
		} else if (gi==1) {
			objtmp.record(idc_G_[1], &KDchan.iKD, RecDt)
			objtmp.label("IKD")
		}
		gvec.append(objtmp)
	}
	
	if (SectionListCount(KAList)>0) {
		KAchan = new KAAnalysis()
		objtmp = new Vector()
		if (gi==0) {
			objtmp.record(idc_G_[1], &KAchan.nKA, RecDt)
			objtmp.label("nKA")
			gvec.append(objtmp)
			objtmp = new Vector()
			objtmp.record(idc_G_[1], &KAchan.lKA, RecDt)
			objtmp.label("lKA")
		} else if (gi==1) {
			objtmp.record(idc_G_[1], &KAchan.iKA, RecDt)
			objtmp.label("IKA")
		}
		gvec.append(objtmp)
	}

	if (SectionListCount(CaTList)>0) {
		CaTchan = new CaTAnalysis()
		objtmp = new Vector()
		if (gi==0) {
			objtmp.record(idc_G_[1], &CaTchan.sCaT, RecDt)
			objtmp.label("sCaT")
			gvec.append(objtmp)
	
			objtmp = new Vector()
			objtmp.record(idc_G_[1], &CaTchan.hCaT, RecDt)
			objtmp.label("hCaT")
		} else if (gi==1) {
			objtmp.record(idc_G_[1], &CaTchan.iCaT, RecDt)
			objtmp.label("ICaT")
		}
		gvec.append(objtmp)
	}
	
}

proc SaveSections() { local i, opt
// SaveSections( flag, sections )

	opt=0
	if ( (object_id(simsecs_G_,1) == -1) && (numarg()==0) )	{
		opt = 1
		seclist=1
	}
	
	if (numarg()>0) {
		if (argtype(1)==0)	{
		opt = 1
		seclist=$1
		}
	}
	
	if ( numarg() > 1 )	opt = 2
	
	if (object_id(simsecs_G_,1) == -1) simsecs_G_ = new SectionList()

	if ( opt == 0 ) {
		if (verbosity > 1) printf("SaveSections: simsecs_G_ unchanged\n")
		return
	} else if (opt == 1) {
	
		simsecs_G_ = new SectionList()
		if (seclist==1) {
			MakeSecList(simsecs_G_,"soma[0]","Handle[25]","MainTrunk[0]","FieldC[15]","FieldB[15]","Tines[287]",\
				"Tines[7]","Tines[621]","Tines[786]","Tines[1080]","Tines[895]")
		} else if (seclist==2) {
			MakeSecList(simsecs_G_,"Handle[40]","MainTrunk[0]","FieldC[25]","Tines[786]","Tines[0]")
		} else if (seclist==3) {
			MakeSecList( simsecs_G_, "MainTrunk[0]" )
		} else if (seclist==4) {
			MakeSecList(simsecs_G_,"MainTrunk[0]","Tines[287]","MainTrunk[40]",\
				"Tines[7]","Tines[786]","Tines[1080]","Handle[20]")
		} else if (seclist==5) {
			if (AEC) {
				MakeSecList(simsecs_G_,"soma[0]","Tine0[7]","Tine3[6]")//,"Tine5[5]")//,\
// 					"Tine8[8]","Tine17[4]")
			} else MakeSecList(simsecs_G_,"soma[0]","Handle[20]","MainTrunk[0]","Tines[530]","Tines[17]","Tines[621]",\
				"Tines[1101]","Tines[1200]","Tines[795]")
		}
	} else if (opt==2) {
		if (argtype(2)==1) {
			simsecs_G_ = $o2
		} else if (argtype(2)==2) {
			for i=2,numarg() {
				forall ifsec $si simsecs_G_.append()
			}
		} else {
			if (verbosity > 0) printf("SaveSections: Sections must be specified as string or SectionList\nsimsecs_G_ unchanged\n")
		}
	}

}