// This file contains procedures used to set the membrane properties

proc init_params() {
	// Initialize user-defined membrane parameters

	//passive membrane parameters determined by fitting Cs data
	//nonuniform membrane resistance
	celsius = 35
	Rn = 50			// nodal resistivity
 	Vleak = -50		// leak reversal in Control
	Vrest = -60		// resting potential -63 in control
	target = 59.42      	// input resistance target in CsCl
	global_ra = 117.15	// axial resistance
	rmsoma = 2.8201e+05	// specific membrane resistance at soma
	rmend = 33938		// specific membrane resistance at distal dendrites
	rmhalfdis = 186.04	// inflection point of sigmoidal membrane resistance curve
	rmsteep = 31.681	// steepness of sigmoidal membrane resistance curve
	Vrest = -63		// resting potential
	Cm = 1.5657		// specific capacitance

	// nonuniform H-conductance distribution
	minq = 0.038794		// minimum H-conductance
	maxq = 20		// maximum H-conductance
	qhalfdis = 316.65	// inflection point of sigmoidal H-conductance curve
	qsteep = 20.001		// steepness of sigmoidal H-conductance curve
	erevh_h = -25		// These paramaters govern the reversal potential and temperature dependence of the H-conductance (see mod file) 
	zeta_h = 10.917	
	a0_h = 0.00028504	
	qten_h = 4.5
			
	// Active conductance values
	gkap = 0.1		//Proximal A-type conductance
	gkad = 0.1		//Distal A-type conductance
	dprox = 100             //Distance at which A-type switches from proximal to distal mod file
	dlimit = 450            //Cutoff for A-type conductance increase  
	dslope = 0.005		//Rate of increase of A-type conductance

	gkdr = 0.1		//KDR conductance
	gkdraxon = 0.05		//Axonal KDR conductance
				
	// User-defined properties for Na current
	gnabar = 0.05		// Total Na conductance
	gnanode = 1		// Na conducance at node	
	nalimit = 275         	// cut-off limit for increase of slow-type sodium conductance proportion
	distslowfrac = 1	// fraction of slow inactivating-channels in distal dendrites
	isegfactor = 10		
	isegfrac = 0.8		
	basalfrac=0.6		// fraction of somatic conductance in basal dendrites (as compared to apical)
}

// Insert and initialize passive membrane properties
proc insert_pass() {
	forsec all_dendrites {
		insert pas 
		Ra = global_ra 
		e_pas = Vleak
	}
	forall {insert spines scale_spines = 1.0}
	forsec all_dendrites 	{ insert spines scale_spines=1.0 }
	forsec no_spines 	{ scale_spines = 1.0 }
	forsec basals 		{ scale_spines = 2.7 }
	forsec med_spines_rad 	{ scale_spines = 1.2 }
	forsec max_spines_rad 	{ scale_spines = 1.6 }
	forsec thin_rad 	{ scale_spines = 2.6 }
	forsec med_spines_LM 	{ scale_spines = 1.1 }
	forsec thin_LM 		{ scale_spines = 1.2 }
	init_pass(Cm)
}

// Insert and initialize H-conductance
proc insert_h() {
	forsec all_dendrites { 
		insert h 
		gbar_h = 0 
	}
	init_h(a0_h*1000)	
}

// Insert and initialize A-type K conductance
proc insert_ka() {
	forall {
		insert kap
    		insert kad
    		gbar_kap = 0
    		gbar_kad = 0 
	}
	init_ka()
}

// Insert and initialize KDR-conductance
proc insert_kdr() {
	forall {
		insert kdr
		gbar_kdr = 0
	}
	init_kdr()
}

// Insert and initialize Na-conductance
proc insert_na() {
	forall {
		insert naslow
		insert nafast
    	}
	init_na()
}

//  Initialize parameters for nonuniform Rm
proc init_pass() {
	// sigmoidal decrease in Rm from rmsoma to rmend
	// rmsoma and rmend, and rmhalfdis and rmsteep are set in init_params
	// function applies to basal dendrites also
	Cm = $1
	soma.sec distance()
	forsec all_dendrites {
		for (x) {
			dis = distance(x)
			rmpoint = rmend + (rmsoma - rmend)/(1 + exp((dis - rmhalfdis)/rmsteep))
			g_pas(x) = 1/(rmpoint/scale_spines)
		}
		cm = Cm * scale_spines 
		Ra = global_ra
	}
	hill {
		cm = Cm
		for (x) {
			dis = distance(x)
			rmpoint = rmend + (rmsoma - rmend)/(1 + exp((dis - rmhalfdis)/rmsteep))
			g_pas(x) = 1/(rmpoint/scale_spines)
		}
	}
	iseg {
		cm = Cm
		for (x) {
			dis = distance(x)
			rmpoint = rmend + (rmsoma - rmend)/(1 + exp((dis - rmhalfdis)/rmsteep))
			g_pas(x) = 1/(rmpoint/scale_spines) 
		}
	}
	for i = 0,2 inode[i] { 
	   	cm = Cm
		for (x) {
			dis = distance(x)
			rmpoint = rmend + (rmsoma - rmend)/(1 + exp((dis - rmhalfdis)/rmsteep))
			g_pas(x) = 1/(rmpoint/scale_spines) 
		}
	}
	for i = 0,1 node[i] {
		cm = Cm
		for (x) {
			dis = distance(x)
			rmpoint = rmend + (rmsoma - rmend)/(1 + exp((dis - rmhalfdis)/rmsteep))
			g_pas(x) = 1/(rmpoint/scale_spines) 
		}
	}
	
}
			
//  Initialize nonuniform H conductance
proc init_h() {
	// scale up sag conductance as a sigmoidal function of distance from the soma
	// q = minq+(maxq-minq)/(1+exp((-dis-qhalfdis)/qsteep))
	// basals are set to the same value as the soma
	// set minq, maxq, and qhalfdis and qsteep in init_params
	a0_h = $1/1000
	soma.sec distance()
	forsec all_dendrites {
		for (x) {
			dis = distance(x)
			hpoint = minq + (maxq - minq)/(1 + exp(-(dis - qhalfdis)/qsteep))
			gbar_h(x) = hpoint * scale_spines*(1 - hblock)
		}
	}
	
	soma.sec distance()
	forsec basals {
		for (x) {
			gbar_h(x) = minq * scale_spines * (1 - hblock)
		}
	}
}

// Initialize parameters for Kdr channels
proc init_kdr() {
	forall {
		gbar_kdr = gkdr * (1 - kdrblock)
	}
	hill {
		gbar_kdr = gkdraxon * (1 - kablock)
	}
	iseg {
		gbar_kdr = gkdraxon * (1 - kablock)
	}
	for i = 0,2 inode[i] { 
	   	gbar_kdr = gkdraxon * (1 - kablock)  
	}
	for i = 0,1 node[i] {
		gbar_kdr = gkdraxon * 0.2 * (1 - kablock)
	} 
}


// Initialize parameters for Na channels
proc init_na() {
	nalimit = 300
	soma.sec distance()
	somaA {
		gbar_nafast = gnabar * (1 - nablock)
		gbar_naslow = 0
	}
	forsec new_basals {
		gbar_nafast = gnabar * (1 - nablock) * basalfrac
		gbar_naslow = 0    		
	}
	forsec apicals {
		for (x) { 
			xdist = distance(x)
			if (xdist>nalimit) {
				xdist=nalimit
			}
			gbar_nafast(x)=gnabar*(1-(distslowfrac*(xdist/nalimit)))
			gbar_naslow(x)=gnabar-gbar_nafast(x)
		}
    	}
	hill {
		gbar_nafast= gnabar * (1 - nablock)
		gbar_naslow = 0
	}
	iseg {
		for (x) { 
			gbar_nafast(x) = 0
			if (x < isegfrac) {
	        		gbar_nafast(x) = gnabar * (1 - nablock)
	       		} else {
	            		gbar_nafast(x) = isegfactor * gnabar * (1 - nablock)
	       		}
    		}
	}
	for i = 0,2 inode[i] { 
	   	gbar_nafast = 0  
	}
	for i = 0,1 node[i] {
		gbar_nafast = gnanode * (1 - nablock)
	} 
}

// Initialize parameters for A-type K channels
proc init_ka() {
	soma.sec distance()
	somaA {
		gbar_kap = gkap * (1 - kablock)
	}
	forsec new_basals {
		gbar_kap = gkap * (1 - kablock)
	}
	forsec apicals {
		for (x) { 
			xdist = distance(x)
			if (xdist > dlimit) {
			  	xdist = dlimit
			}
			if (x!=1) {
				gbar_kap(x) = 0
				gbar_kad(x) = 0
				if (xdist > dprox) {
				       	gbar_kad(x) = gkad * (1 + xdist * dslope) * (1 - kablock)
				} else {
				        gbar_kap(x) = gkap * (1 + xdist * dslope) * (1 - kablock)
				}
			}
    		}
	}
}

proc reinit_props() {
	init_pass(Cm)		
	init_h(a0_h*1000)	
	init_kdr()
	init_ka()
	init_na()
}