# $Id: network.py,v 1.125 2011/06/10 15:10:05 samn Exp $ from pyinit import * from geom import * import random gGID = 0 # global ID for cells class Population: "Population of cells" # cell_type -- pyr, bas, olm # n -- number of cells in the population # x, y, z -- initial position for the first Cell # dx -- an increment of the x-position for the cell location # amp, dur, delay -- parameters for the IClamp in the soma # Spikes are stored in ltimevec (times) and lidvec (cell # within the population) def __init__(self, cell_type, n, x, y, z, dx, amp, dur, delay): global gGID self.cell = [] # List of cells in the population self.nc = [] # NetCon list for recording spikes self.n = n # number of cells self.x = x self.y = y self.z = z self.ltimevec = h.List() # list of Vectors for recording spikes, one per cell self.lidvec = h.List() self.nssidx = {} self.nseidx = {} self.ncsidx = {} self.nceidx = {} for i in range(n): self.cell.append(cell_type(x+i*dx,y,z,gGID)) self.cell[-1].somaInj.amp = amp self.cell[-1].somaInj.dur = dur self.cell[-1].somaInj.delay = delay self.nc.append(h.NetCon(self.cell[-1].soma(0.5)._ref_v, None, sec=self.cell[-1].soma)) self.ltimevec.append(h.Vector()) #NB: each NetCon gets own Vectors for recording. needed to avoid multithreading crash self.lidvec.append(h.Vector()) self.nc[-1].record(self.ltimevec[-1],self.lidvec[-1],gGID) # record cell spikes with gGID gGID = gGID + 1 # inc global cell ID def set_r(self, syn, r): for c in self.cell: c.__dict__[syn].syn.r = r class MSpec: # this class uses matlab to make a spectrogram def __init__(self,vlfp,maxfreq,nsamp,dodraw): #make a spectrogram using matlab h("jjj=name_declared(\"nqspec\")") h("if(jjj){nqsdel(nqspec) print \"deleted nqspec\"}") h("objref nqspec") vslfp = h.Vector() vslfp.copy(vlfp) vslfp.sub(vlfp.mean()) h.nqspec = h.matspecgram(vslfp,1e3/h.dt,maxfreq,nsamp,dodraw) self.nqspec = h.nqspec def powinrange(self,minf,maxf): # get scalar power in range of frequencies nn = self.nqspec.select(-1,"f","[]",minf,maxf) if nn == 0: return 0 h("jnk = 0") h("vec.resize(0)") for i in self.nqspec.ind: mystr = "vec.copy(nqspec.get(\"pow\"," mystr += str(int(i)) mystr += ").o)" h(mystr) h("jnk += vec.sum()") jnk = h.jnk return jnk / nn def vecinrange(self,minf,maxf): # get vector of power in range of frequencies vs time nn = self.nqspec.select(-1,"f","[]",minf,maxf) if nn == 0: return None h("objref vjnk") h("vjnk=new Vector()") h.vec.resize(0) for i in self.nqspec.ind: mystr = "vec.copy(nqspec.get(\"pow\"," mystr += str(int(i)) mystr += ").o)" h(mystr) if h.vjnk.size()==0: h.vjnk.copy(h.vec) else: h.vjnk.add(h.vec) h.vjnk.div(self.nqspec.ind.size()) return h.vjnk class Network: def __init__(self,noise=True,connections=True,DoMakeNoise=True,iseed=1234,UseNetStim=True,wseed=4321,scale=1.0,MSGain=1.0,SaveConn=False): import math print("Setting Cells") self.pyr = Population(cell_type=PyrAdr,n=int(math.ceil(800*scale)), x= 0, y=0, z=0, dx=50, amp= 50e-3, dur=1e9, delay=2*h.dt) self.bas = Population(cell_type=Bwb, n=int(math.ceil(200*scale)), x=10, y=0, z=0, dx=50, amp= 0, dur= 0, delay=2*h.dt) self.olm = Population(cell_type=Ow, n=int(math.ceil(200*scale)), x=20, y=0, z=0, dx=50, amp=-25e-3, dur=1e9, delay=2*h.dt) # psr = sensor cell to estimate the E->E connections self.psr = Population(cell_type=PyrAdr,n=1, x= 0, y=0, z=0, dx=50, amp= 50e-3, dur=1e9, delay=2*h.dt) self.cells = [self.pyr, self.bas, self.olm, self.psr] self.iseed = iseed # seed for noise inputs self.noise = noise self.DoMakeNoise = DoMakeNoise self.UseNetStim = UseNetStim self.wseed = wseed # seed for 'wiring' self.MSGain = MSGain # gain for MS weights self.RecPyr = False self.SaveConn = SaveConn if connections: print("Setting Connections") self.set_all_conns() def set_noise_inputs(self,simdur): #simdur only used for make_all_noise if self.DoMakeNoise: if self.UseNetStim: self.make_all_NetStims(simdur,self.iseed) else: self.make_all_noise(simdur,self.iseed) else: self.load_all_noise() print("Done!") def load_all_noise(self): #load noise from data files print("Loading Noise") print("to PYR") self.b_pyr_somaAMPAf=self.load_spikes("spike_noise_pyr_800_soma_AMPA_ISI_1_N_10000_noise_1.npy",self.pyr,"somaAMPAf",0.05e-3) self.b_pyr_Adend3AMPAf=self.load_spikes("spike_noise_pyr_800_Adend3_AMPA_ISI_1_N_10000_noise_1.npy",self.pyr,"Adend3AMPAf",0.05e-3) self.b_pyr_somaGABAf=self.load_spikes("spike_noise_pyr_800_soma_GABA_ISI_1_N_10000_noise_1.npy",self.pyr,"somaGABAf",0.012e-3) self.b_pyr_Adend3GABAf=self.load_spikes("spike_noise_pyr_800_Adend3_GABA_ISI_1_N_10000_noise_1.npy",self.pyr,"Adend3GABAf",0.012e-3) self.b_pyr_Adend3NMDA=self.load_spikes("spike_noise_pyr_800_Adend3_NMDA_ISI_100_N_100_noise_1.npy", self.pyr,"Adend3NMDA",6.5e-3) print("to BAS") self.b_bas_somaAMPAf=self.load_spikes("spike_noise_bas_200_soma_AMPA_ISI_1_N_10000_noise_1.npy",self.bas,"somaAMPAf",w=0.02e-3) self.b_bas_somaGABA=self.load_spikes("spike_noise_bas_200_soma_GABA_ISI_1_N_10000_noise_1.npy",self.bas,"somaGABAf",w=0.2e-3) self.b_bas_somaGABAf=self.load_spikes("spike_noise_bas_200_soma_GABAf_ISI_150_N_65_noise_0.npy",self.bas,"somaGABAss",w=1.6e-3) print("to OLM") self.b_olm_somaAMPAf=self.load_spikes("spike_noise_olm_200_soma_AMPA_ISI_1_N_10000_noise_1.npy",self.olm,"somaAMPAf",w=0.02e-3) self.b_olm_somaGABAf=self.load_spikes("spike_noise_olm_200_soma_GABA_ISI_1_N_10000_noise_1.npy",self.olm,"somaGABAf",w=0.2e-3) self.b_olm_somaGABAss=self.load_spikes("spike_noise_olm_200_soma_GABAss_ISI_150_N_65_noise_0.npy",self.olm,"somaGABAss",w=1.6e-3) #this should be called @ beginning of each sim - done in an FInitializeHandler in run.py def init_NetStims(self): # h.mcell_ran4_init(self.iseed) for i in range(len(self.nrl)): rds = self.nrl[i] sead = self.nrlsead[i] rds.MCellRan4(sead,sead) rds.negexp(1) # print i,rds,sead #creates NetStims (and associated NetCon,Random) - provide 'noise' inputs #returns next useable value of sead def make_NetStims(self,po,syn,w,ISI,time_limit,sead): po.nssidx[syn] = len(self.nsl) #index into net.nsl po.ncsidx[syn] = len(self.ncl) #index into net.ncl for i in range(po.n): cel = po.cell[i] ns = h.NetStim() ns.interval = ISI ns.noise = 1 ns.number = (1e3 / ISI) * time_limit ns.start = 0 nc = h.NetCon(ns,cel.__dict__[syn].syn) nc.delay = h.dt * 2 # 0 nc.weight[0] = w rds = h.Random() rds.negexp(1) # set random # generator using negexp(1) - avg interval in NetStim rds.MCellRan4(sead,sead) # seeds are in order, shouldn't matter ns.noiseFromRandom(rds) # use random # generator for this NetStim #ns.start = rds.discunif(0,1e3) # start inputs random time btwn 0-1e3 ms to avoid artificial sync #rds.MCellRan4(sead,sead) # reinit rand # generator self.nsl.append(ns) self.ncl.append(nc) self.nrl.append(rds) self.nrlsead.append(sead) sead = sead + 1 po.nseidx[syn] = len(self.nsl)-1 po.nceidx[syn] = len(self.ncl)-1 return sead # setup recording of pyramidal cell inputs, assumes using NetCon,NetStims def RecPYRInputs(self): self.RecPyr = True self.NCV = {} self.sys = ['somaAMPAf', 'Adend3AMPAf', 'somaGABAf', 'Adend3GABAf'] sys=self.sys for s in sys: self.NCV[s] = [] sidx = self.pyr.ncsidx[s] eidx = self.pyr.nceidx[s] for i in range(sidx,eidx+1): self.NCV[s].append(h.Vector()) self.ncl[i].record(self.NCV[s][-1]) # make an NQS with pyramidal cell input times def setnqin(self): try: h.nqsdel(self.nqin) except: pass self.nqin = h.NQS("id","sy","vt") nqin=self.nqin nqin.odec("vt") jdx = 0 for s in self.sys: sidx = self.pyr.ncsidx[s] eidx = self.pyr.nceidx[s] idx = 0 for i in range(0,len(self.NCV[s])): nqin.append(idx,jdx,self.NCV[s][i]) idx = idx + 1 jdx = jdx + 1 # make a histogram of pyramidal cell spike outputs def mkspkh(self,binsz): snq=self.snq snq.verbose = 0 self.spkh = h.List() for i in range(0,800): if snq.select("id",i) > 0: vt = snq.getcol("t") self.spkh.append(vt.histogram(0,h.tstop,binsz)) else: self.spkh.append(h.Vector()) snq.verbose=1 def make_all_NetStims(self,simdur,rdmseed): print("Making NetStims") # h.mcell_ran4_init(self.iseed) self.nsl = [] #NetStim List self.ncl = [] #NetCon List self.nrl = [] #Random List for NetStims self.nrlsead = [] #List of seeds for NetStim randoms # numpy.random.seed(rdmseed) # initialize random # generator print("Making Noise") print("to PYR") rdtmp = rdmseed # starting sead value - incremented in make_NetStims rdtmp=self.make_NetStims(po=self.pyr, syn="somaAMPAf", w=0.05e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.pyr, syn="Adend3AMPAf", w=0.05e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.pyr, syn="somaGABAf", w=0.012e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.pyr, syn="Adend3GABAf", w=0.012e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.pyr, syn="Adend3NMDA", w=6.5e-3, ISI=100,time_limit=simdur, sead=rdtmp) print("to BAS") rdtmp=self.make_NetStims(po=self.bas, syn="somaAMPAf", w=0.02e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.bas, syn="somaGABAf", w=0.2e-3, ISI=1, time_limit=simdur, sead=rdtmp) print("to OLM") #rdtmp=self.make_NetStims(po=self.olm, syn="somaAMPAf", w=0.02e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.olm, syn="somaAMPAf", w=0.0625e-3, ISI=1, time_limit=simdur, sead=rdtmp) rdtmp=self.make_NetStims(po=self.olm, syn="somaGABAf", w=0.2e-3, ISI=1, time_limit=simdur, sead=rdtmp) #setup medial septal inputs to OLM and BASKET cells, note that MSGain can be 0 == no effect ns = h.NetStim() ns.interval = 150 ns.noise = 0 # NO randomness for the MS inputs ns.number = (1e3 / 150.0) * simdur self.nsl.append(ns) for i in range(self.bas.n): # MS inputs to BASKET cells nc = h.NetCon(ns,self.bas.cell[i].__dict__["somaGABAss"].syn) nc.delay = 2*h.dt nc.weight[0] = 1.6e-3 * self.MSGain self.ncl.append(nc) for i in range(self.olm.n): # MS inputs to OLM cells nc = h.NetCon(ns,self.olm.cell[i].__dict__["somaGABAss"].syn) nc.delay = 2*h.dt nc.weight[0] = 1.6e-3 * self.MSGain self.ncl.append(nc) def make_all_noise(self,simdur,rdmseed): # create noise for simdur milliseconds numpy.random.seed(rdmseed) # initialize random # generator import math print("Making Noise") fctr = (simdur+simdur/2) / 10000.0 print("to PYR") self.b_pyr_somaAMPAf=self.make_spikes(self.pyr,"somaAMPAf",0.05e-3,self.pyr.n,"soma",1,math.ceil(10000*fctr),1,simdur) self.b_pyr_Adend3AMPAf=self.make_spikes(self.pyr,"Adend3AMPAf",0.05e-3,self.pyr.n,"Adend3",1,math.ceil(10000*fctr),1,simdur) self.b_pyr_somaGABAf=self.make_spikes(self.pyr,"somaGABAf",0.012e-3,self.pyr.n,"soma",1,math.ceil(10000*fctr),1,simdur) self.b_pyr_Adend3GABAf=self.make_spikes(self.pyr,"Adend3GABAf",0.012e-3,self.pyr.n,"Adend3",1,math.ceil(10000*fctr),1,simdur) self.b_pyr_Adend3NMDA=self.make_spikes(self.pyr,"Adend3NMDA",6.5e-3,self.pyr.n,"Adend3",100,math.ceil(100*fctr),1,simdur) print("to BAS") self.b_bas_somaAMPAf=self.make_spikes(self.bas,"somaAMPAf",0.02e-3,self.bas.n,"soma",1,math.ceil(10000*fctr),1,simdur) self.b_bas_somaGABA=self.make_spikes(self.bas,"somaGABAf",0.2e-3,self.bas.n,"soma",1,math.ceil(10000*fctr),1,simdur) self.b_bas_somaGABAf=self.make_spikes(self.bas,"somaGABAss",1.6e-3,self.bas.n,"soma",150,math.ceil(65*fctr),0,simdur) print("to OLM") self.b_olm_somaAMPAf=self.make_spikes(self.olm,"somaAMPAf",0.02e-3,self.olm.n,"soma",1,math.ceil(10000*fctr),1,simdur) self.b_olm_somaGABAf=self.make_spikes(self.olm,"somaGABAf",0.2e-3,self.olm.n,"soma",1,math.ceil(10000*fctr),1,simdur) self.b_olm_somaGABAss=self.make_spikes(self.olm,"somaGABAss",1.6e-3,self.olm.n,"soma",150,math.ceil(65*fctr),0,simdur) def make_conn(self, preN, postN, conv): conn = numpy.zeros((postN,conv),dtype=numpy.int16) for i in range(postN): conn[i,:]=random.sample(list(range(preN)),conv) return conn def set_all_conns(self): random.seed(self.wseed) # initialize random # generator for wiring print("PYR -> X , NMDA") # src, trg, syn, delay, weight, conv self.pyr_bas_NM=self.set_connections(self.pyr,self.bas, "somaNMDA", 2, 1.15*1.2e-3, 100) self.pyr_olm_NM=self.set_connections(self.pyr,self.olm, "somaNMDA", 2, 1.0*0.7e-3, 10) self.pyr_pyr_NM=self.set_connections(self.pyr,self.pyr, "BdendNMDA",2, 1*0.004e-3, 25) print("PYR -> X , AMPA") self.pyr_bas_AM=self.set_connections(self.pyr,self.bas, "somaAMPAf",2, 0.3*1.2e-3, 100) self.pyr_olm_AM=self.set_connections(self.pyr,self.olm, "somaAMPAf",2, 0.3*1.2e-3, 10) self.pyr_pyr_AM=self.set_connections(self.pyr,self.pyr, "BdendAMPA",2, 0.5*0.04e-3, 25) print("BAS -> X , GABA") #self.bas_bas_GA=self.set_connections(self.bas,self.bas, "somaGABAf",2, 1.0e-3, 60)#orig 1 #self.bas_bas_GA=self.set_connections(self.bas,self.bas, "somaGABAf",2, 2 * 1.5*1.0e-3, 60)#new 1 self.bas_bas_GA=self.set_connections(self.bas,self.bas, "somaGABAf",2, 3 * 1.5*1.0e-3, 60)#new 2 self.bas_pyr_GA=self.set_connections(self.bas,self.pyr, "somaGABAf",2, 2 * 2*0.18e-3, 50)#new 1 print("OLM -> PYR , GABA") #self.olm_pyr_GA=self.set_connections(self.olm,self.pyr, "Adend2GABAs",2, 3*6.0e-3, 20)#original weight value self.olm_pyr_GA=self.set_connections(self.olm,self.pyr, "Adend2GABAs",2, 4.0 * 3*6.0e-3, 20)#new weight value #pyramidal to PSR cell -- for testing only print("PYR -> PSR, AMPA/NMDA") self.pyr_psr_NM=self.set_connections(self.pyr,self.psr, "BdendNMDA",2, 1*0.004e-3, 25) self.pyr_psr_AM=self.set_connections(self.pyr,self.psr, "BdendAMPA",2, 0.5*0.04e-3, 25) def set_conn_weight(self, conn, weight): for nc in conn: nc.weight[0] = weight def set_connections(self,src,trg,syn,delay,w,conv): conn = self.make_conn(src.n,trg.n,conv) nc = [] for post_id, all_pre in enumerate(conn): for j, pre_id in enumerate(all_pre): nc.append(h.NetCon(src.cell[pre_id].soma(0.5)._ref_v, trg.cell[post_id].__dict__[syn].syn, 0, delay, w, sec=src.cell[pre_id].soma)) if self.SaveConn: try: print(self.nqcon.size()) except: self.nqcon = h.NQS("id1","id2","w","syn") self.nqcon.strdec("syn") for post_id, all_pre in enumerate(conn): for j, pre_id in enumerate(all_pre): self.nqcon.append(src.cell[pre_id].id,trg.cell[post_id].id,w,syn) return nc def load_spikes(self,fn,po,syn,w,time_limit=10000): fn = os.path.join("data",fn) events = numpy.load(fn) print("Begin setting events...", po) print(events.shape) for i,ii in enumerate(events): ii=ii[ii<=time_limit] po.cell[i].__dict__[syn].append(ii) po.cell[i].__dict__[syn].syn.Vwt = w print("End setting events") return events def make_spikes(self,po,syn,w,cellN,comp,ISI,eventN,noise,time_limit): events = numpy.random.exponential(ISI, (cellN,eventN))*noise+numpy.repeat(ISI,cellN*eventN).reshape((cellN,eventN))*(1-noise) events = numpy.cumsum(events,axis=1) print("Begin setting events...", po) print(events.shape) for i,ii in enumerate(events): ii=ii[ii<=time_limit] po.cell[i].__dict__[syn].append(ii) po.cell[i].__dict__[syn].syn.Vwt = w print("End setting events") return events def rasterplot(self,sz=2): pon = 0 if h.g[0] == None: h.gg() col = [2, 4, 3, 1] for po in self.cells: id = h.Vector() tv = h.Vector() for i in range(po.n): id.append(po.lidvec[i]) tv.append(po.ltimevec[i]) id.mark(h.g[0],tv,"O",sz,col[pon],1) pon += 1 h.g[0].exec_menu("View = plot") def setrastervecs(self): self.myidvec = h.Vector() #IDs and firing times for ALL cells self.mytimevec = h.Vector() for po in self.cells: for i in range(po.n): self.myidvec.append(po.lidvec[i]) self.mytimevec.append(po.ltimevec[i]) # setsnq - make an NQS with ids, spike times, types def setsnq(self): try: h.nqsdel(self.snq) except: pass self.snq = h.NQS("id","t","ty") ty = 0 vec = h.Vector() for po in self.cells: for i in range(po.n): self.snq.v[0].append(po.lidvec[i]) self.snq.v[1].append(po.ltimevec[i]) vec.resize(po.lidvec[i].size()) vec.fill(ty) self.snq.v[2].append(vec) ty += 1 # setfnq - make an NQS with ids, firing rates, types def setfnq(self,skipms=200): try: self.snq.tog("DB") except: self.setsnq() try: h.nqsdel(self.fnq) except: pass self.fnq = h.NQS("id","freq","ty") tf = h.tstop - skipms ty = 0 for po in self.cells: for i in range(po.n): id = po.cell[i].id n = float( self.snq.select("t",">",skipms,"id",id) ) self.fnq.append(id, n*1e3/tf, ty) ty += 1 # pravgrates - print average firing rates using self.fnq def pravgrates(self,skipms=200): try: self.fnq.tog("DB") except: self.setfnq(skipms) ty = 0 tf = float( h.tstop - skipms ) for po in self.cells: self.fnq.select("ty",ty) vf = self.fnq.getcol("freq") if vf.size() > 1: print("ty: ", ty, " avg rate = ", vf.mean(), "+/-", vf.stderr(), " Hz") else: print("ty: ", ty, " avg rate = ", vf.mean(), "+/-", 0.0 , " Hz") ty += 1 def calc_lfp(self): # lfp is modeled as a difference between voltages in distal apical and basal compartemnts self.vlfp = h.Vector(self.pyr.cell[0].Adend3_volt.size()) #lfp in neuron Vector for cell in self.pyr.cell: self.vlfp.add(cell.Adend3_volt) self.vlfp.sub(cell.Bdend_volt) self.vlfp.div(len(self.pyr.cell)) # normalize lfp by amount of pyr cells self.lfp=numpy.array(self.vlfp.to_python()) # convert to python array (so can do PSD) def calc_specgram(self,maxfreq,nsamp,dodraw,skipms=0): self.calc_lfp() if skipms > 0: vtmp = h.Vector() vtmp.copy(self.vlfp,skipms/h.dt,self.vlfp.size()-1) self.MSpec = MSpec(vtmp,maxfreq,nsamp,dodraw) else: self.MSpec = MSpec(self.vlfp,maxfreq,nsamp,dodraw) def calc_psd(self,fig=3): self.calc_lfp() t0 = 200 # reject first ms of the signal fmax = 200 # upper limit for a periodogram frequency div = int(1000/h.dt/(2*fmax)) # downsample the signal tr = [3, 12] # Theta frequency range gr = [30, 80] # Gamma frequency range t0i = int(t0/h.dt) if t0i > len(self.lfp): print("LFP is too short! (<200 ms)") return 0,0,0,0,0,0 pylab.figure(fig) pylab.clf() pylab.subplot(2,1,1) # plot LFP pylab.plot(numpy.array(list(range(len(self.lfp))))*h.dt, self.lfp) pylab.subplot(2,1,2) # plot periodogram data = self.lfp[t0i::div] # downsample data Pxx, freqs = pylab.psd(data-data.mean(), Fs=1000/h.dt/div) # calculate FFT tind = numpy.where((freqs>=tr[0]) & (freqs<=tr[1]))[0] # index where for theta frequences gind = numpy.where((freqs>=gr[0]) & (freqs<=gr[1]))[0] # index where for gamma frequences self.tp = Pxx[tind].mean() * numpy.diff(tr) # integral over theta power self.gp = Pxx[gind].mean() * numpy.diff(gr) # integral over gamma power self.ti = self.get_lim_max(Pxx, tind) # index of the frequency with a maximal power in theta range self.gi = self.get_lim_max(Pxx, gind) # index of the frequency with a maximal power in gamma range self.tf = freqs[self.ti] self.gf = freqs[self.gi] pylab.scatter(self.tf, 10*numpy.log10(Pxx[self.ti]), 100, 'b','o') pylab.scatter(self.gf, 10*numpy.log10(Pxx[self.gi]), 100, 'r','o') pylab.xlim(0,fmax) def get_lim_max(self, data, ind): # return the position of the maximal element in data located in the postion indexed by ind return ind[data[ind].argmax()] #make the Network - use params in rseed.txt if the file exists -- makes it easier to run a batch #if rseed.txt doesn't exist, the Network is created with default params try: fp = open("./rseed.txt","r") ls = fp.readlines() ISEED = int(ls[0]) WSEED = int(ls[1]) MSG = 1.0 if len(ls) > 2: MSG = float(ls[2]) fp.close() #create the network net = Network(noise=True,connections=True,DoMakeNoise=True,iseed=ISEED,UseNetStim=True,wseed=WSEED,scale=1.0,MSGain=MSG) print("set network from rseed.txt : iseed=",ISEED,", WSEED=",WSEED,", MSG = ",MSG) except: net = Network() print("set network from default constructor") #setup some variables in hoc def sethocix(): h("PYRt=0") h("BASKETt=1") h("OLMt=2") h("PSRt=3") h("CTYP.o(PYRt).s=\"PYRt\"") h("CTYP.o(BASKETt).s=\"BASKETt\"") h("CTYP.o(OLMt).s=\"OLMt\"") h("CTYP.o(PSRt).s=\"PSRt\"") h("ix[PYRt]=0") h("ixe[PYRt]=799") h("ix[BASKETt]=800") h("ixe[BASKETt]=999") h("ix[OLMt]=1000") h("ixe[OLMt]=1199") h("ix[PSRt]=1200") h("ixe[PSRt]=1200") h("numc[PYRt]=800") h("numc[BASKETt]=200") h("numc[OLMt]=200") h("numc[PSRt]=1") sethocix()