//****************************************
// GC-Ball.hoc         Josef Bischofberger
// modified by JM Schulz 2020
//****************************************
// simple model of a young DG granule cell 
// Passive parameters from Heigele et al. 2016, Nature Neuroscience 19:263-70
// Cap = 25 pF
// Rin = 4 GOhm with Rm = 100000 Ohm cm^2
// 

debug_mode = 0

load_file("fixnseg.hoc")

somaSurface = 0
denSurface= 0
proxSurface = 0
midSurface = 0
distSurface = 0
axonSurface = 0
totalSurface = 0

axonGC = 1

denSomaL = 10 	//  um
somaLongL = 15
somaShortL = 10
axonL = 20		//  50 um

denD0 = 2.5	// um soma border
denD1 = 1.5 	// proximal diam um
denD2 = 1.5	// mid
denD3 = 1.2	// dist
denD4 = 0.75

fA1=0.4		// fraction of prox dendritic membrane
fA2=0.3		// mid
fA3=0.3		// distal

gen1 = ndend	// number of prox branches
gen2 = 2		// mid
gen3 = 2		// distal

//*********************************************
//			Topology
//*********************************************
//dendrites: prox +  middle   +   distal 
objref proxDen, middleDen, distalDen //, periSoma

//periSoma = new SectionList()
proxDen = new SectionList()
middleDen = new SectionList()
distalDen = new SectionList()
nSections =  gen1 + gen1*gen2 + gen1*gen2*gen3

create soma[4], dendrite[nSections],axon, AIS 

soma[0] connect soma[1](0), 1
soma[1] connect soma[2](0), 1
soma[2] connect soma[3](0), 1

//for (i=0; i<4; i=i+1){soma[i] periSoma.append()}

// proximal dendrites connect to soma
for i =0,gen1-1 {soma[3] connect dendrite[i](0), 1
			dendrite[i] proxDen.append() 
}

// middle dendrites
for i =0,gen1-1 {
	for j =0, gen2-1{ dendrite[i] connect dendrite[gen1+(gen2*i)+j](0), 1
				  dendrite[gen1+(gen2*i)+j] 	middleDen.append() 	
	}
}

// distal dendrites
for i =0, (gen1*gen2 -1) {
	for j =0, gen3-1{ dendrite[i+gen1] connect dendrite[gen1+gen1*gen2+(gen3*i)+j](0), 1
				 dendrite[gen1+gen1*gen2+(gen3*i)+j] distalDen.append() 	
	}
}

/*print "prox"
forsec proxDen print secname()
print "middle"
forsec middleDen print secname()
print "distal"
forsec distalDen print secname()*/

soma[0] connect AIS(0), 0
AIS connect axon(0), 1


//forall{Ra = Rint }
//*************************************************************
//         		Geometry
//**************************************************************
// mature cell: dendrite 250 µm, diam 1.5:0.7, AIS nseg= 11, axon nseg = 101
// Cap = 100 ms/ 4GOhm = 25 pF => 2500 um^2 => ;  L = 79.61 um
// Cm = 1 uF/cm^2 = 0.01 uF/mm^2 = 0.01 pF/um^2
// lenght Cylinder = Area / 2PI r 


proc adjustL() {
 denMemArea=0  
 somMemArea=0 

 //*************** soma ****************
 firstDiam =0
 lastDiam =0
 soma[1]  {
	pt3dclear()
 	for i = -99, 97{
 	   yi = sqrt(1-i*i/10000)  //circle x^2 + y^2 = 1
 	   pt3dadd(somaLongL/2*i/100,0,0,yi*somaShortL )
 	}
	firstDiam = sqrt(1-99*99/10000)*somaShortL
	lastDiam = yi*somaShortL 
 	nseg=11  	               
 }
 soma[0]  {nseg=1 L = 0.1    		diam(0:1) = firstDiam:firstDiam} //2.3  
 soma[2]  {nseg=3 L = 3   			diam(0:1) = lastDiam:denD0}
 soma[3]  {nseg=11 L = denSomaL		diam(0:1) = denD0:denD1}

 if (axonGC) {
 	AIS {nseg=11  L = 32  	diam(0:0.3) = firstDiam:1   //1.2
				 		diam(0.3:1) = 1:0.3	} 
	axon  {nseg=11 L = axonL	diam = 0.3} 
 }else{
 	AIS  {nseg =1 L=0.01 diam = 0.3}
	axon {nseg =1 L=0.01 diam = 0.3}
 }
 // to do : axon surface
 if (axonGC) {
 	AIS		{axonSurface = secArea() }
 	axon		{axonSurface = axonSurface + secArea()}
 }

 somaSurface  = 0
 for i = 0,3 {	soma[i] {	somaSurface = somaSurface +secArea()}
 }


 //*************** dendrites **************
 //CylinderL = (1-propDendrite)*100*Cap/(12.00 *PI)  //Length of a 12um Cylinder
 //calculate dendritic Area after substracting a constant soma-dend compartment

 denMemArea = 100*Cap - somaSurface - axonSurface 
 
  
 L1 = fA1*denMemArea/(PI*(denD1+denD2)/2) /gen1
 L2 = fA2*denMemArea/(PI*(denD2+denD3)/2) /(gen1*gen2)
 L3 = fA3*denMemArea/(PI*(denD3+denD4)/2) /(gen1*gen2*gen3)

 TDL = denSomaL + gen1*L1 + gen1*gen2*L2 + gen1*gen2*gen3*L3
// check better:
// TDL = 0
// forall {TDL=TDL+L}


 forsec proxDen   {nseg=11 L = L1 diam(0:1) = denD1 :denD2}
 forsec middleDen {nseg=11 L = L2 diam(0:1) = denD2 :denD3}
 forsec distalDen {nseg=11 L = L3 diam(0:1) = denD3 :denD4}
 //if (mature == 1){geom_nseg_shared()} // readjust nseg according to d_lambda(1kHz)

 denSurface= 0

 forsec proxDen   {proxSurface = proxSurface + secArea()
			   //print secname()
			   //print "Area =",  secArea()
 }
 print "prox denSurface = ", proxSurface

 forsec middleDen {midSurface = midSurface + secArea()}
 print "middle denSurface = ", midSurface

 forsec distalDen {distSurface= distSurface+ secArea()}
 print "distal denSurface = ", distSurface

 denSurface= proxSurface + midSurface + distSurface
 totalSurface = somaSurface+ axonSurface + denSurface 

 access soma[1]

 print ""
 print "relative contribution of dendrites to Cell membrane = ", denSurface/totalSurface 
 print ""
 printf("Dendrites L1=%.2f um, L2=%.2f um, L3=%.2f um, TDL=%.2f um \n", L1, L2, L3, TDL)

 print "\nmeasured  Surface of Dendrites = ", denSurface, "um^2" 
 print   "predicted Surface of Dendrites = ", denMemArea, "um^2" 

 if (axonGC) {
   print "measured  Surface of axon      = ", axonSurface, "um^2" 
 }
 print "\nmeasured Cell Surface = ", totalSurface , "um^2" 
 print   "required total Cap    = ", Cap, "pF" 

}

adjustL() // initialize Morphology which is updated later 



