objref G1,vVec,indVec,nc, nil,apVec,G2,fn,voltMat,voltVec,axonVec,axonMat,dendVec,dendMat,stateMatNa,stateMatMut,stateMatNax

getcwd()

proc runSCN2A(){
voltMat = new Matrix(100.01/dt,600)
axonMat = new Matrix(100.01/dt,600)
dendMat = new Matrix(100.01/dt,600)
inactAtVinitRation = 1
simnum=0
counter = 0
soma nc = new NetCon(&v(.5), nil)
nc.threshold = -5 // watch out! only one threshold per presyn location
apVec = new Vector(1100)
indVec = new Vector()
nc.record(indVec)
G1 = new Graph()

na_somaOrig = na_soma1	
na_aisOrig     = na_ais1
	init()
	
	
for (aa=1;aa>=0;aa-=0.1){
	na_soma1 = na_somaOrig*aa
	na_ais1 = na_aisOrig*aa
	IClamp[0].amp = 0.5
	init_channels()
counter = 0
	printf("runing %f %% of Na density\n",aa*100)	

	//na_soma		=	0.5*500		

	tstop = 600
	IClamp[0].del = 500
	IClamp[0].dur = 100

		voltVec = new Vector(10001)
		axonVec = new Vector(10001)
		dendVec = new Vector(10001)
	p = 0
	
			flg0 = 1
		flg1 = 1
	while(IClamp[0].amp<=3){
		init()
		timeSteps = (tstop-t)/dt
	
		counter=0
		printf("runing %f",IClamp[0].amp)

		while (t<tstop) {
			G1.line(t,v(0.5))
			voltVec.x(counter) = v(0.5)
			axonVec.x(counter) = axon[0].v(1)
			dendVec.x(counter) = apic[37].v(0.5)
			if(v(0.5)>-10 && flg0){
				print "soma ap"
				flg0 = 0
			}
			if(axon[0].v(1)>-10&& flg1){
				print "ap axon"
				flg1 = 0
			}

			fadvance()
			counter+=1
		}
		IClamp[0].amp +=0.05
		printf("printing %f %f\n",aa,simnum)
		voltMat = voltMat.setcol(simnum,voltVec) 
		axonMat = axonMat.setcol(simnum,axonVec) 
		dendMat = dendMat.setcol(simnum,dendVec) 
		simnum +=1
		
	}
	G1.flush()
	
	//apVec.x(simnum) = IClamp[0].amp
	//printf("somaCurrent1.amp is %f\n",apVec.x(simnum))
		//apVec.x(simnum) = somaCurrent1.amp
		
}
fn = new File()
fn.wopen("Volts.csv")
voltMat.fprint(fn,"%f,","\n")
fn.close()
fn.wopen("VoltsAxon.csv")
axonMat.fprint(fn,"%f,","\n")
fn.close()
fn.wopen("VoltsDend.csv")
dendMat.fprint(fn,"%f,","\n")
fn.close()
}
/*