// $Id: shufmua.hoc,v 1.16 2012/01/11 19:39:35 samn Exp $ ncells = 300 nsec = 30 binsz = 5 refrac = 2.5 jitterdt = 0 declare("reg",0) declare("noisespks", 0) objref vs[ncells],vh[ncells],vmua,rdm,vy,vx vmua=new Vector() rdm=new Random() vy=new Vector() vx=new Vector() sz = (1e3/binsz) * nsec rate = 1 // rate in Hz sampr = 1e3 / binsz objref nqp[300] SPECTY=0 PRESM=0 //* checkrefrac(vector,refrac) proc checkrefrac () { local i localobj v1 v1=$o1 for i=1,v1.size-1 { if(v1.x(i)-v1.x(i-1)<$2) v1.x(i) = v1.x(i) + $2 } } //* myshuf(vec,nshuffles,rdm) proc myshuf () { local idx,i,j,k,n localobj v1,rdm,vm,vr,vr2 v1=$o1 n=$2 rdm=$o3 vm=new Vector(v1.size) vm.resize(0) for i=0,v1.size-1 if(v1.x(i)>0) vm.append(i) vr=new Vector(n) vr2=new Vector(n) rdm.discunif(0,vm.size-1) vr.setrand(rdm) rdm.discunif(0,v1.size-1) vr2.setrand(rdm) for i=0,n-1 { idx = vm.x(vr.x(i)) j = vr2.x(i) k = v1.x(idx) v1.x(idx) = v1.x(j) v1.x(j) = k } } //* applyjitter(vec,rdm,dt) proc applyjitter () { local i,jdt localobj v1,rdm,vj v1=$o1 rdm=$o2 jdt=$3 vj=new Vector(v1.size) rdm.uniform(-jdt,jdt) vj.setrand(rdm) v1.add(vj) for i=0,v1.size-1 if(v1.x(i)<0) v1.x(i)=0// make sure no neg #s v1.sort() } //* addnoise(numspikes,rdm) proc addnoise () { local i,ns localobj vec ns=$1 rdm.uniform(0,nsec*1e3) vec=new Vector(ns) for i=0,ncells-1 { vec.setrand(rdm) vs[i].append(vec) vs[i].sort() } } //* initcells proc initcells () { local tt,isi,i,nshuf localobj vr nspks = rate * nsec rdm.ACG(1234*nshuf) for i=0,ncells-1 if(vs[i]==nil) vs[i]=new Vector() else vs[i].resize(0) if(reg) { isi = 1e3 / rate tt = isi for i=0,nspks-1 { vs[0].append(tt) tt += isi } } else { vs[0].resize(nspks) rdm.uniform(0,nsec*1e3) vs[0].setrand(rdm) vs[0].sort() } checkrefrac(vs[0],refrac) for i=1,ncells-1 { vs[i].copy(vs[0]) applyjitter(vs[i],rdm,jitterdt) checkrefrac(vs[i],refrac) } if(noisespks) addnoise(noisespks,rdm) } //* histcells - make spike counts per time for each cell proc histcells () { local i for i=0,ncells-1 { if(vh[i]==nil) vh[i]=new Vector() vh[i].hist(vs[i],0,sz,binsz) } } //* shufhist(nshuf) proc shufhist () { local i,nshuf nshuf=$1 for i=1,ncells-1 { vh[i].copy(vh[0]) myshuf(vh[i],nshuf,rdm) if(refrac) checkrefrac(vh[i],refrac) } } //* mkmua proc mkmua () { local i vmua.resize(vh[0].size()) vmua.fill(0) for i=0,ncells-1 vmua.add(vh[i]) vmua.sub(vmua.mean()) } //* plotrast proc plotrast () { local i vrsz(0,vx,vy) for i=0,ncells-1 { for j=0,vs[i].size-1 { vx.append(vs[i].x(j)) vy.append(i) } } vy.mark(g,vx,"O",2,1) g.exec_menu("View = plot") } //* setjitter(jitterdt) proc setjitter () { jitterdt = $1 initcells() histcells() mkmua() } //* jittertest(maxjitter,jitterinc) proc jittertest () { local i for(jitterdt=0;jitterdt<=$1;jitterdt+=$2) { print "jitterdt is " , jitterdt setjitter(jitterdt) {nqsdel(nqp[i]) nqp[i]=getspecnq(vmua,sampr,SPECTY,PRESM)} // nqp.gr("pow","f",0,1,1) i += 1 } } proc nqpg () { nqp[$1].gr("pow","f",0,1,1) g.exec_menu("View = plot") } //* main setjitter(0) gg()