class fitdata: def __init__(self): self.phi=-1 self.theta=-1 self.xpar=[] self.ypar=[] self.zpar=[] self.apar=[[],[],[]] def setPhi(self,phi): self.phi=phi def setTheta(self,theta): self.theta=theta def setXpar(self,par): self.xpar=self.addPar(par) def setYpar(self,par): self.ypar=self.addPar(par) def setZpar(self,par): self.zpar=self.addPar(par) def setApar(self,par,axis): self.apar[axis]=self.addPar(par) def addPar(self,inpar): par=[] for p in inpar: par.append(float(p)) return par import glob,math,sys,subprocess,copy,os,commands,time#,argcomplete sys.path.append('$PANDAPATH/macro/tpc/FOPI/mberger') pandapath=os.environ.get('PANDAPATH') sys.path.append(pandapath+'/macro/tpc/FOPI/python/argparse-1.2.1') import argparse from anaFile import anaFile from functions import * def setGraph(g,title,col=1): g.SetTitle(title) g.SetName(title) #col=col%len(colors) if title.find(args.special)!=-1 and args.special!="None": g.SetMarkerStyle(29) g.SetMarkerSize(1) g.SetMarkerColor(1) else: g.SetMarkerStyle(20) g.SetMarkerSize(0.5) g.SetMarkerColor(colors[ fcol[col]%len(colors) ]) return def fillGraph(hist,graph,fnum=0,mini={},maxi={},append=1): nbins=hist.GetNbinsX() projs=[] tempg={} tempg['RMS']=ROOT.TGraphErrors() setGraph(tempg['RMS'],'RMS '+hist.GetTitle()+titles[fnum],fnum) tempg['Mean']=ROOT.TGraphErrors() setGraph(tempg['Mean'],'Mean '+hist.GetTitle()+titles[fnum],fnum) if append==0: tempg=graph if graph.get('RMS',None)==None: graph['RMS']=ROOT.TGraphErrors() setGraph(graph['RMS'],'RMS',fnum) if graph.get('Mean',None)==None: graph['Mean']=ROOT.TGraphErrors() setGraph(graph['Mean'],'Mean',fnum) pcounter=0 i=args.binwidth while 1: projs.append(hist.ProjectionY(hist.GetName()+titles[fnum]+'_'+str(i),i-args.binwidth,i)) tempg['RMS'].SetPoint(pcounter,i-args.binwidth+1,projs[-1].GetRMS()) tempg['RMS'].SetPointError(pcounter,1./2.,projs[-1].GetRMSError()) tempg['Mean'].SetPoint(pcounter,i-args.binwidth+1,projs[-1].GetMean()) tempg['Mean'].SetPointError(pcounter,1./2.,projs[-1].GetMeanError()) mini['RMS']=min(mini.get('RMS',1000),projs[-1].GetRMS()) maxi['RMS']=max(maxi.get('RMS',0),projs[-1].GetRMS()) mini['Mean']=min(mini.get('Mean',1000),projs[-1].GetMean()) maxi['Mean']=max(maxi.get('Mean',0),projs[-1].GetMean()) pcounter+=1 i+=1 if (i>=nbins): break if append!=0: graph.append(tempg) return copy.deepcopy(projs) def fitSlicesFromFile(graph,name,fnum,rfile,mini={},maxi={}): graph.append({}) print name graph[-1]['Sigma']=rfile.Get(name+'_Sigma') graph[-1]['Mean']=rfile.Get(name+'_Mean') graph[-1]['Const']=rfile.Get(name+'_Const') setGraph(graph[-1]['Sigma'],graph[-1]['Sigma'].GetTitle()+'_'+titles[fnum],fnum) setGraph(graph[-1]['Mean'] ,graph[-1]['Mean'].GetTitle() +'_'+titles[fnum],fnum) setGraph(graph[-1]['Const'],graph[-1]['Const'].GetTitle()+'_'+titles[fnum],fnum) projs=[] nbins=graph[-1]['Sigma'].GetN() #x3, y3 = ROOT.Double(0), ROOT.Double(0) x1, y1 = ROOT.Double(0), ROOT.Double(0) x2, y2 = ROOT.Double(0), ROOT.Double(0) for ibin in range(nbins): projs.append(rfile.Get(name+str(ibin))) graph[-1]['Sigma'].GetPoint(ibin,x1,y1) mini['Sigma']=min(mini.get('Sigma',1000),y1) maxi['Sigma']=max(maxi.get('Sigma',0),y1) graph[-1]['Mean'].GetPoint(ibin,x2,y2) mini['Mean']=min(mini.get('Mean',1000),y2) maxi['Mean']=max(maxi.get('Mean',0),y2) return projs def fitslices(graph,hist,fnum=0,mini={},maxi={},funcf='gaus'): if funcf=='landau': func=ROOT.TF1("landau","landau",0,hist.GetRMS(2)*2) func.SetParLimits(1,0,10000) nbins=hist.GetNbinsX() projs=[] graph.append({}) graph[-1]['Sigma']=ROOT.TGraphErrors() setGraph(graph[-1]['Sigma'],'Sigma '+hist.GetTitle()+titles[fnum],fnum) graph[-1]['Mean']=ROOT.TGraphErrors() setGraph(graph[-1]['Mean'],'Mean '+hist.GetTitle()+titles[fnum],fnum) graph[-1]['Const']=ROOT.TGraphErrors() setGraph(graph[-1]['Const'],'Const '+hist.GetTitle()+titles[fnum],fnum) pcounter=0 binwidth=args.binwidth i=binwidth while 1: projs.append(hist.ProjectionY(hist.GetName()+str(i),i-binwidth,i)) if funcf=='landau': projs[-1].Fit(func,"QN") data={} data['nsig']=func.GetParameter(2) data['nsigerr']=func.GetParError(2) data['nmean']=func.GetParameter(1) data['nmeanerr']=func.GetParError(2) data['nconst']=func.GetParameter(0) data['nconsterr']=func.GetParError(0) else: data=fitres(projs[-1]) graph[-1]['Sigma'].SetPoint(pcounter,i-binwidth+1,data['nsig']) graph[-1]['Sigma'].SetPointError(pcounter,float(1./2.),data['nsigerr']) graph[-1]['Mean'].SetPoint(pcounter,i-binwidth+1,data['nmean']) graph[-1]['Mean'].SetPointError(pcounter,float(1/2.),data['nmeanerr']) graph[-1]['Const'].SetPoint(pcounter,i-binwidth+1,data['nconst']) graph[-1]['Const'].SetPointError(pcounter,float(1/2.),data['nconsterr']) mini['Sigma']=min(mini.get('Sigma',1000),data['nsig']) maxi['Sigma']=max(maxi.get('Sigma',0),data['nsig']) mini['Mean']=min(mini.get('Mean',1000),data['nmean']) maxi['Mean']=max(maxi.get('Mean',0),data['nmean']) pcounter+=1 i+=1 if (i>=nbins): break #graph[-1]['Sigma'].SetName(hist.GetTitle()+'_Sigma') #graph[-1]['Mean'].SetName(hist.GetTitle()+'_Mean') #graph[-1]['Const'].SetName(hist.GetTitle()+'_Const') # tmpcanv=ROOT.TCanvas() # graph[-1]['RMS'].Draw('AP') # u=raw_input('next dingfes') return projs def drawprojs(canv,histos,fitname,graph,fits,fitf='gaus'): for i in range(len(histos)): opt='' canv[i].Divide(int(math.ceil(math.sqrt(len(histos[i][0])))), int(math.ceil(math.sqrt(len(histos[i][0]))))) fits.append([]) for f in range(len(histos[i])): hcounter=-1 fits[-1].append([]) for h in histos[i][f]: hcounter+=1 fits[-1][-1].append(ROOT.TF1(fitname+str(i)+str(f)+str(hcounter), fitf, h.GetXaxis().GetXmin(), h.GetXaxis().GetXmax())) x1, y1 = ROOT.Double(0), ROOT.Double(0) x2, y2 = ROOT.Double(0), ROOT.Double(0) x3, y3 = ROOT.Double(0), ROOT.Double(0) graph[i][f]['Const'].GetPoint(hcounter,x1,y1) fits[-1][-1][-1].SetParameter(0,y1) graph[i][f]['Mean'].GetPoint(hcounter,x2,y2) fits[-1][-1][-1].SetParameter(1,y2) graph[i][f]['Sigma'].GetPoint(hcounter,x3,y3) fits[-1][-1][-1].SetParameter(2,y3) ''' fits[-1][-1][-1].SetParameter(0,graph[i][f]['Const'].GetBinContent(hcounter+int(args.binwidth/2))) fits[-1][-1][-1].SetParameter(1,graph[i][f]['Mean'].GetBinContent(hcounter+int(args.binwidth/2))) fits[-1][-1][-1].SetParameter(2,graph[i][f]['Sigma'].GetBinContent(hcounter+int(args.binwidth/2))) ''' if h.GetTitle().find(args.special)!=-1: fits[-1][-1][-1].SetLineColor(1) h.SetLineColor(1) else: fits[-1][-1][-1].SetLineColor(colors[ fcol[f]%len(colors) ]) h.SetLineColor(colors[ fcol[f]%len(colors) ]) canv[i].cd(hcounter+1) h.Draw(opt) fits[-1][-1][-1].Draw("same") opt="same" canv[i].Update() return def drawprojsfile(c,histos,fits,leg=None): for i in range(len(histos)): for ibin in range(len(histos[i][0])): opt='' c.cd(1).cd(ibin%9+1) maxi=0 for ifile in range(len(histos[i])): maxi=max(maxi,histos[i][ifile][ibin].GetMaximum()) histos[i][ifile][ibin].Draw(opt) opt='same' histos[i][0][ibin].SetMaximum(maxi*1.1) if ibin%9==8: c.Update() c.SaveAs(args.pfile) c.SaveAs(args.pfile) return def drawgraphs(c,g,gmin,gmax,leg,save=True): c.Divide(2,1) c.cd(1).SetPad(0,0,0.8,1) c.cd(1).Divide(2,3) c.cd(2).SetPad(0.8,0,1,1) c.cd(2) leg.Draw() for g1 in range(len(g)): opt='AP' for g2 in g[g1]: c.cd(1).cd(2*g1+1) minmult=0.9 maxmult=1.1 if gmin[g1]['RMS']<0: minmult=1.1 if gmax[g1]['RMS']<0: maxmult=0.9 g2['RMS'].GetYaxis().SetRangeUser(gmin[g1]['RMS']*minmult,gmax[g1]['RMS']*maxmult) g2['RMS'].Draw(opt) c.cd(1).cd(2*g1+2) minmult=0.9 maxmult=1.1 if gmin[g1]['Mean']<0: minmult=1.1 if gmax[g1]['Mean']<0: maxmult=0.9 g2['Mean'].GetYaxis().SetRangeUser(gmin[g1]['Mean']*minmult,gmax[g1]['Mean']*maxmult) g2['Mean'].Draw(opt) opt='P' # gfile.cd() # g2['RMS'].Write() # g2['Mean'].Write() c.Update() if args.pfile!='None' and save: c.SaveAs(args.pfile) outfile.cd() c.Write() return def drawgraphf(c,g,gmin,gmax,leg,save=True): c.Divide(2,1) c.cd(1).SetPad(0,0,0.8,1) c.cd(1).Divide(2,3) c.cd(2).SetPad(0.8,0,1,1) c.cd(2) leg.Draw() for g1 in range(len(g)): opt='AP' for g2 in g[g1]: c.cd(1).cd(2*g1+1) g2['Sigma'].SetMaximum(gmax[g1]['Sigma']*1.1) g2['Sigma'].SetMinimum(gmin[g1]['Sigma']*0.9) g2['Sigma'].Draw(opt) c.cd(1).cd(2*g1+2) g2['Mean'].SetMaximum(gmax[g1]['Mean']*1.1) g2['Mean'].SetMinimum(gmin[g1]['Mean']*0.9) g2['Mean'].Draw(opt) #gfile.cd() #g2['Sigma'].Write() #g2['Mean'].Write() opt='P' c.Update() if args.pfile!='None' and save: c.SaveAs(args.pfile) outfile.cd() c.Write() return def writeGraphToFile(graph,text,f): x,y=ROOT.Double(0),ROOT.Double(0) graph.GetPoint(0,x,y) if x!=0: print 'blah',x for missing in range(int(x)): f.write(text+' '+str(missing)+' '+str(y)+'\n') for ibin in range(graph.GetN()): #x,y=ROOT.Double(0),ROOT.Double(0) graph.GetPoint(ibin,x,y) f.write(text+' '+str(x)+' '+str(y)+'\n') return parser=argparse.ArgumentParser(description='calculate a parametrization for the cluster error for several files') parser.add_argument('files',help='all the files to analyze',type=str,nargs='+') parser.add_argument('--list',help='files given in list',action='store_true') parser.add_argument('--titles',help='the titles to use for the files',type=str,nargs='+',default=['None']) parser.add_argument('--hlp',help='print help',action='store_true') parser.add_argument('--new',help='new analysis. create the histofiles',action='store_true') parser.add_argument('--oldana',help='dont run ana of clerrcalib again',action='store_true') parser.add_argument('--basepath',help='the basepath to store all the files into',type=str,default='./') parser.add_argument('--fitfile',help='the file to write the fitinfo into',type=str,default='fitinfo.txt') parser.add_argument('--listfile',help='name of the new listfile',type=str,default='histlist.txt') parser.add_argument('--outfile',help='the root-file to write the canvas',type=str,default='clerr.root') parser.add_argument('--events',help='number of events',type=int, default=-1) parser.add_argument('--fit',help='plot only data from existing fitfile',action='store_true') parser.add_argument('--phi',help='enable phi dependent fitting',action='store_true') parser.add_argument('--nobatch',help='disable batch mode for debugging',action='store_true') parser.add_argument('--plot',help='plot all pictures into the respective folders',action='store_true') parser.add_argument('--pfile',help='the file to plot all the stuff into',type=str,default='None') parser.add_argument('--genfit',help='use genfit track residuals to cluster residuals',action='store_true') parser.add_argument('--track',help='use genfit track to mc residuals',action='store_true') parser.add_argument('--binwidth',help='width of bins',type=int,default=2) parser.add_argument('--check',help='only check all files in list',action='store_true') parser.add_argument('--batch',help='use batch farm',action='store_true') parser.add_argument('--multi',help='run multiple jobs',action='store_true') parser.add_argument('--createList',help='create only the histlist, dont run clerrcalib_v2',action='store_true') parser.add_argument('--entry',help='do only entry number',type=int,default=-1) parser.add_argument('--rmsfile',help='the rms file to use',type=str,default='./rms.txt') parser.add_argument('--cutCov',help='use cut not proj for scaling calculation',action='store_true') parser.add_argument('--depfile',help='the file with the batch job dependencies',type=str,default='NONE') parser.add_argument('--ax3Zax',help='cut on angle of axis 3 and z axis (inside window)',type=float,nargs=2,default=[0,180]) parser.add_argument('--clsize2D',help='cut on 2D clustersize (inside window)',type=int,nargs=2,default=[0,99999999]) parser.add_argument('--clsize3D',help='cut on 3D clustersize (inside window)',type=int,nargs=2,default=[0,99999999]) parser.add_argument('--special',help='special entry in list',type=str,default="None") parser.add_argument('--quiet',help='turn on root batch mode to plot everything more quietly',action='store_true') parser.add_argument('--runtime',help='the time to request at batch',type=str,default='02:00:00') parser.add_argument('--make_lookup',help='enables creation of lookup table',action='store_true') parser.add_argument('--stderr',help='disable cuts the std err cant provide',action='store_true') parser.add_argument('--sort',help='sort the axes according to their orientation',action='store_true') parser.add_argument('--subintervall',help='time intervall in which the jobs are submitted (min)',type=int,default=2) #argcomplete.autocomplete(parser) args=parser.parse_args() if args.quiet: sys.argv.append( '-b-' ) import ROOT if args.quiet: ROOT.gROOT.SetBatch(True) ROOT.gROOT.ProcessLine(".x rootlogon.C") ROOT.gROOT.ProcessLine('gROOT->SetStyle("Plain")') #ROOT.gStyle.SetPaperSize(30,30) set_palette() if args.hlp: parser.print_help() exit(0) dependency="" print 'The depfile used:',args.depfile if args.depfile!='NONE': depfile=open(args.depfile,'read') for line in depfile: dependency+=line.replace('\n','')+':' dependency=dependency[0:len(dependency)-1] print 'followind dependencies are used:',dependency fcol=[] if args.list: listfile=open(args.files[0],'r') files=[] titles=[] fcolcounter=-1 for line in listfile: fcolcounter+=1 if line[0]=='#' or line[0]=='\n': continue line=line.replace('\n','') words=line.split(';') files.append(str(words[0])) if len(words)<2: fname=words[0].split('/')[-1] fnamep=fname.split('_') fnamep_end=fnamep.pop().split('.') fnamep=fnamep+fnamep_end endindex=fnamep.index('reco') startindex=fnamep.index('1') tit='' for i in range(startindex+1,endindex): tit+='_{0}'.format(fnamep[i]) tit=tit[1:len(tit)] titles.append(tit) else: titles.append(str(words[1]).replace('\n','')) titles[-1]=titles[-1].replace(' ','') fcol.append(fcolcounter) listfile.close() if not args.new: args.basepath=os.path.dirname(args.files[0]) else: files=args.files fcol=range(len(files)) #exit() if len(files)<10: colors=getColors2() else: colors=getColors() if args.titles[0]=='None' and not args.list: titles=range(len(files)) for i in range(len(titles)): titles[i]=str(titles[i]) elif not args.list: titles=args.titles if args.check: for f in files: print '' rfile=ROOT.TFile(f,'read') tree=rfile.Get('cbmsim') if tree==None: print 'no tree in:',f continue print f,':',tree.GetEntries() exit() if not args.fit and args.new: counter=-1 fitfile=open(args.basepath+'/'+args.fitfile,'w') if args.new and args.entry==-1: if args.genfit: args.listfile=args.listfile.replace('.txt','_gf.txt') if args.track: args.listfile=args.listfile.replace('.txt','_tr.txt') histlist=open(args.basepath+'/'+args.listfile,'w') run_procs={} fitfile.close() for f in files: counter+=1 proc=[] proc.append('python') proc.append('FOPI/mberger/clerrcalib_v2.py') if args.oldana: if args.genfit: proc.append(args.basepath+'/'+titles[counter]+'_gf.root') elif args.track: proc.append(args.basepath+'/'+titles[counter]+'_tr.root') else: proc.append(args.basepath+'/'+titles[counter]+'_mc.root') else: proc.append(f) proc.append('--new') proc.append('--outfile') if args.genfit: outname=titles[counter]+'_gf.root' elif args.track: outname=titles[counter]+'_tr.root' else: outname=titles[counter]+'_mc.root' if args.sort: outname=outname.replace('.root','_sorted.root') proc.append(outname) proc.append('--fitfile') proc.append(args.fitfile) proc.append('--basepath') proc.append(args.basepath) if not args.nobatch: proc.append('--batch') if args.events!=-1: proc.append('--events') proc.append(str(args.events)) if args.phi: proc.append('--phi') proc.append(titles[counter].split('_')[1]) if args.plot: proc.append('--pfile') proc.append('{0}/{1}'.format(args.basepath,outname.replace('.root','.pdf') ) ) if args.genfit: proc.append('--genfit') if args.track: proc.append('--track') if f.find('cutCov')!=-1: proc.append('--cutCov') if args.sort: proc.append('--sort') proc.append('--title') proc.append(titles[counter]) proc.append('--rmsfile') proc.append(str(args.rmsfile)) proc.append('--ax3Zax') proc.append(str(args.ax3Zax[0])) proc.append(str(args.ax3Zax[1])) proc.append('--clsize2D') proc.append(str(args.clsize2D[0])) proc.append(str(args.clsize2D[1])) proc.append('--clsize3D') proc.append(str(args.clsize3D[0])) proc.append(str(args.clsize3D[1])) proc.append('--binwidth') proc.append(str(args.binwidth)) if args.stderr: proc.append('--stderr') # proc.append('-b') if args.new and args.entry==-1: histlist.write(args.basepath+'/'+outname+';'+titles[counter]+'\n') if args.createList: continue fitfile=open(args.basepath+'/'+args.fitfile,'a') fitfile.write('Name: '+titles[counter]+'\n') fitfile.close() # print 'phi=',titles[counter].split('_')[1] if args.entry!=-1 and args.entry!=counter: continue runfilename='clerrcalib_'+titles[counter] runfile=open(runfilename,'w') runfile.write('pwd\n') runfile.write('export PANDAPATH=/home/mberger/fopiroot/fopiroot_dev/ \n') runfile.write('cd $PANDAPATH\n') runfile.write('pwd\n') runfile.write('echo '+runfilename+'\n') runfile.write('. ./setenv\n') for p in proc: runfile.write(p+' ') runfile.write('\n') runfile.close() os.system('chmod a+x '+runfilename) if args.batch: outname=f.split('/')[-1] commStr='' commStr+='qsub -l nodes=x64 -l walltime='+args.runtime+' -l mem=1500m -q farmq -N clerr'+titles[counter]+' -o ~/sgeout/clerr_'+f.replace('/','_')+'_'+outname+'_'+titles[counter] commStr+=' -j oe ' if dependency!="": commStr+='-W depend=afterok:'+dependency+' ' commStr+=runfilename status,jobid=commands.getstatusoutput(commStr) print 'job:',jobid print outname elif args.multi: runCommand=['gnome-terminal','--title=clerr'+titles[counter],'--disable-factory','-e','./clerrcalib_'+titles[counter]] print runCommand run_procs[titles[counter]]=subprocess.Popen(runCommand) else: p=subprocess.call(proc) fitfile.close() if args.new and args.entry!=-1: histlist.close() running=True while running: running=False for p in run_procs: run_procs[p].poll() if run_procs[p].returncode!=0: running=True print 'job',p,'is still running' time.sleep(10) if args.batch or args.createList: exit() if args.fit : if not args.list: fitfile=open(args.basepath+'/'+args.files[0],'r') else: fitfile=open(args.basepath+'/'+args.fitfile,'r') fitdataArr=[] for line in fitfile: if line=='\n': continue if line.find('Name:')!=-1: words=line.split() phi=words[1].split('_')[1] fitdataArr.append(fitdata()) fitdataArr[-1].setPhi( float(phi) ) continue words=line.split(';') if words[0]=='X': print words fitdataArr[-1].setXpar(words[1:len(words)-1]) elif words[0]=='Y': fitdataArr[-1].setYpar(words[1:len(words)-1]) elif words[0]=='Z': fitdataArr[-1].setZpar(words[1:len(words)-1]) elif words[0][0]=='A': fitdataArr[-1].setApar(words[1:len(words)-1],int(words[0][1])-1) graphs={} graphsA={} counter=[0,0,0] coords=['X','Y','Z'] coordsA=['A1','A2','A3'] for data in fitdataArr: for i in range(len(coords)): if graphs.get(coords[i],None)==None: graphs[coords[i]]=[] for k in range(len(data.xpar)): graphs[coords[i]].append(ROOT.TGraph()) graphs[coords[i]][-1].SetTitle(coords[i]+': parameter '+str(k)) graphs[coords[i]][-1].SetMarkerStyle(20) graphs[coords[i]][-1].SetMarkerSize(.8) if graphsA.get(coordsA[i],None)==None: graphsA[coordsA[i]]=[] for k in range(len(data.apar[0])): graphsA[coordsA[i]].append(ROOT.TGraph()) graphsA[coordsA[i]][-1].SetTitle(coordsA[i]+': parameter '+str(k)) graphsA[coordsA[i]][-1].SetMarkerStyle(20) graphsA[coordsA[i]][-1].SetMarkerSize(.8) if i==0: par=data.xpar elif i==1: par=data.ypar elif i==2: par=data.zpar for j in range(len(par)): # print counter[i],par[j],data.phi,i,j,coords[i] graphs[coords[i]][j].SetPoint(counter[i],data.phi,par[j]) for j in range(len(data.apar[0])): # print counter[i],data.phi,data.apar[i][j] graphsA[coordsA[i]][j].SetPoint(counter[i],data.phi,data.apar[i][j]) counter[i]+=1 #if args.list and not args.new: if not args.new: hErrT2=[] hErrT=[] hErr=[] hErrR=[] hErrA=[] hErrAProj=[] hErrAr=[] hErrA3=[] hPullUV=[] hprojUV=[] hNum=[] hprojNum=[] hPullxyzC=[] hprojPullxyzC=[] gErrT2=[] gErrT=[] gfErrT=[] gErr=[] gErrR=[] gErrA=[] gErrAr=[] gErrA3=[] gfErrA3=[] gfErrAr=[] gAllErrA=[] gPullUV=[] gfPullUV=[] gNum=[] gfNum=[] gfPullxyzC=[] gPullxyzC=[] gT2min=[] gT2max=[] gTmin=[] gTmax=[] gmin=[] gmax=[] gRmin=[] gRmax=[] gAmin=[] gAmax=[] gArmin=[] gArmax=[] gfArmin=[] gfArmax=[] gA3min=[] gA3max=[] gfA3min=[] gfA3max=[] gAllAmin=[] gAllAmax=[] gPullxyzCmin=[] gPullxyzCmax=[] gfPullxyzCmin=[] gfPullxyzCmax=[] gPullUVmin=[] gPullUVmax=[] gfPullUVmin=[] gfPullUVmax=[] gNummin=[] gNummax=[] gfNummin=[] gfNummax=[] hNum.append([]) hprojNum.append([]) gNum.append([]) gfNum.append([]) gNummin.append({}) gNummax.append({}) gfNummin.append({}) gfNummax.append({}) for i in range(3): hErrT2.append([]) gErrT2.append([]) gT2min.append({}) gT2max.append({}) hErrT.append([]) gErrT.append([]) gfErrT.append([]) gTmin.append({}) gTmax.append({}) hErr.append([]) gErr.append([]) gmin.append({}) gmax.append({}) hErrR.append([]) gErrR.append([]) gRmin.append({}) gRmax.append({}) hErrA.append([]) hErrAProj.append([]) gErrA.append([]) hErrA3.append([]) gErrA3.append([]) gfErrA3.append([]) gfErrAr.append([]) gAllErrA.append({}) gAmin.append({}) gAmax.append({}) gfArmin.append({}) gfArmax.append({}) gAllAmin.append({}) gAllAmax.append({}) hErrAr.append([]) gErrAr.append([]) gArmin.append({}) gArmax.append({}) gA3min.append({}) gA3max.append({}) gfA3min.append({}) gfA3max.append({}) hPullxyzC.append([]) hprojPullxyzC.append([]) gPullxyzCmin.append({}) gPullxyzCmax.append({}) gfPullxyzCmin.append({}) gfPullxyzCmax.append({}) gfPullxyzC.append([]) gPullxyzC.append([]) for i in range(2): hPullUV.append([]) hprojUV.append([]) gPullUV.append([]) gPullUVmin.append({}) gPullUVmax.append({}) gfPullUV.append([]) gfPullUVmin.append({}) gfPullUVmax.append({}) outfile=ROOT.TFile(args.basepath+'/'+args.outfile,'RECREATE') cgErrT2=ROOT.TCanvas('cgErrT2','Error T2',1200,1000) cgErrT=ROOT.TCanvas('cgErrT','Error*MCRes T',1200,1000) cgErr=ROOT.TCanvas('cgErr','Error XYZ',1200,1000) cgErrR=ROOT.TCanvas('cgErrR','Error*MCRes XYZ',1200,1000) cgErrA=ROOT.TCanvas('cgErrA','Error A',1200,1000) cgfErrAr=ROOT.TCanvas('cgfErrAr','Error fA',1200,1000) cgErrAr=ROOT.TCanvas('cgErrAr','Error*MCRes A',1200,1000) cAllErrA=ROOT.TCanvas('cAllErrA','All ErrA',1200,1000) cgErrA3=ROOT.TCanvas('cgErrA3','ErrorA3*MCRes A',1200,1000) cgfErrA3=ROOT.TCanvas('cgfErrA3','ErrorA3*MCRes fA',1200,1000) cgPullUV=ROOT.TCanvas('cgPullUV','Pull UV',1200,1000) cgfPullUV=ROOT.TCanvas('cgfPullUV','Pull UV f',1200,1000) cgNum=ROOT.TCanvas('cgNum','Scaling Factor',1200,1000) cgfNum=ROOT.TCanvas('cgfNum','Scaling Factor Fit',1200,1000) cgPullxyzC=ROOT.TCanvas('cgPullxyzC','Cov XYZ pulls',1200,1000) cgfPullxyzC=ROOT.TCanvas('cgfPullxyzC','Cov XYZ pulls f',1200,1000) chprojA=[] chprojXYZC=[] for i in range(3): chprojA.append(ROOT.TCanvas('chprojA'+str(i),'Error*MCRes A'+str(i)+' projections',1200,1000)) chprojXYZC.append(ROOT.TCanvas('chprojXYZC'+str(i),'Cov xyz-axis cut '+str(i)+' projections',1200,1000)) chprojUV=[] chprojUV.append(ROOT.TCanvas('chprojU','Pull U',1200,1000)) chprojUV.append(ROOT.TCanvas('chprojV','Pull V',1200,1000)) chprojNum=[] chprojNum.append(ROOT.TCanvas('chprojNum','Num',1200,1000)) if args.pfile!='None': cproj=ROOT.TCanvas('cproj','Projection',1200,1000) cproj.Divide(2,1) cproj.cd(1).SetPad(0,0,0.8,1) cproj.cd(1).Divide(3,3) cproj.cd(2).SetPad(0.8,0,1,1) cproj.SaveAs(args.pfile+'[') ctmp=ROOT.TCanvas('ctmp','tmp',1200,1000) clerrCname=['cerr','cerr2','cerrT','cerrT2','cerrA','cerrA2', 'cerrA3','cerrNum','cScale','cErrUV','cPullUV','cerrC','cPullC'] fcounter=-1 lookup_file=open(args.basepath+'/error_lookup.txt','w') lookup_file_num=open(args.basepath+'/error_lookup_num.txt','w') for f in files: fcounter+=1 print f rfile=ROOT.TFile(f,'read') if args.pfile!='None': ctmp.cd() text=ROOT.TText(0.5,0.5,titles[fcounter]) text.SetTextAlign(10) text.SetTextSize(0.05) text.Draw() ctmp.Update() ctmp.SaveAs(args.pfile) for cname in clerrCname: print cname clerrTmp=(rfile.Get(cname)) clerrTmp.cd() clerrTmp.Draw() clerrTmp.SaveAs(args.pfile) print 'all written' #hErrT2[0].append(rfile.Get('hErrLAmpVsZ2')) #hErrT2[1].append(rfile.Get('hErrP1AmpVsZ2')) #hErrT2[2].append(rfile.Get('hErrP2AmpVsZ2')) #hErrT[0].append(rfile.Get('hErrLrAmpVsZ')) #hErrT[1].append(rfile.Get('hErrP1rAmpVsZ')) #hErrT[2].append(rfile.Get('hErrP2rAmpVsZ')) hErr[0].append(rfile.Get('hErrXAmpVsZ2')) hErr[1].append(rfile.Get('hErrYAmpVsZ2')) hErr[2].append(rfile.Get('hErrZAmpVsZ2')) hErrR[0].append(rfile.Get('hErrXrAmpVsZ')) hErrR[1].append(rfile.Get('hErrYrAmpVsZ')) hErrR[2].append(rfile.Get('hErrZrAmpVsZ')) hErrA[0].append(rfile.Get('hErrA1AmpVsZ2')) hErrA[1].append(rfile.Get('hErrA2AmpVsZ2')) hErrA[2].append(rfile.Get('hErrA3AmpVsZ2')) hErrAr[0].append(rfile.Get('hErrA1rAmpVsZ')) hErrAr[1].append(rfile.Get('hErrA2rAmpVsZ')) hErrAr[2].append(rfile.Get('hErrA3rAmpVsZ')) hErrA3[0].append(rfile.Get('hErrA1AmpVsZ3')) hErrA3[1].append(rfile.Get('hErrA2AmpVsZ3')) hErrA3[2].append(rfile.Get('hErrA3AmpVsZ3')) hPullUV[0].append(rfile.Get('hPullUvsZ')) hPullUV[1].append(rfile.Get('hPullVvsZ')) hNum[0].append(rfile.Get('hErrNumVsZ')) hPullxyzC[0].append(rfile.Get('hPullXCov')) hPullxyzC[1].append(rfile.Get('hPullYCov')) hPullxyzC[2].append(rfile.Get('hPullZCov')) #fillGraph(hErrT2[0][-1],gErrT2[0],fcounter,gT2min[0],gT2max[0]) #fillGraph(hErrT2[1][-1],gErrT2[1],fcounter,gT2min[1],gT2max[1]) #fillGraph(hErrT2[2][-1],gErrT2[2],fcounter,gT2min[2],gT2max[2]) for i in range(3): #fillGraph(hErrT[i][-1],gErrT[i],fcounter,gTmin[i],gTmax[i]) #fitSlicesFromFile(gfErrT[i],hErrT[i][-1].GetName(),fcounter,rfile) hErrAProj[i].append(fillGraph(hErrAr[i][-1],gErrAr[i],fcounter,gArmin[i],gArmax[i])) fitSlicesFromFile(gfErrAr[i],hErrAr[i][-1].GetName(), fcounter,rfile,gfArmin[i],gfArmax[i]) fillGraph(hErrAr[i][-1],gAllErrA[i],fcounter,gAllAmin[i],gAllAmax[i],0) fillGraph(hErrA3[i][-1],gErrA3[i],fcounter,gA3min[i],gA3max[i]) fitSlicesFromFile(gfErrA3[i],hErrA3[i][-1].GetName(), fcounter,rfile,gfA3min[i],gfA3max[i]) hprojPullxyzC[i].append(fillGraph(hPullxyzC[i][-1],gPullxyzC[i],fcounter,gPullxyzCmin[i],gPullxyzCmax[i])) fitSlicesFromFile(gfPullxyzC[i],hPullxyzC[i][-1].GetName(),fcounter, rfile,gfPullxyzCmin[i],gfPullxyzCmax[i]) for i in range(2): hprojUV[i].append(fillGraph(hPullUV[i][-1],gPullUV[i],fcounter,gPullUVmin[i],gPullUVmax[i])) #fitslices(gfPullUV[i],hPullUV[i][-1],fcounter,gfPullUVmin[i],gfPullUVmax[i]) fitSlicesFromFile(gfPullUV[i],hPullUV[i][-1].GetName(),fcounter, rfile,gfPullUVmin[i],gfPullUVmax[i]) fillGraph(hErr[0][-1],gErr[0],fcounter,gmin[0],gmax[0]) fillGraph(hErr[1][-1],gErr[1],fcounter,gmin[1],gmax[1]) fillGraph(hErr[2][-1],gErr[2],fcounter,gmin[2],gmax[2]) fillGraph(hErrR[0][-1],gErrR[0],fcounter,gRmin[0],gRmax[0]) fillGraph(hErrR[1][-1],gErrR[1],fcounter,gRmin[1],gRmax[1]) fillGraph(hErrR[2][-1],gErrR[2],fcounter,gRmin[2],gRmax[2]) fillGraph(hErrA[0][-1],gErrA[0],fcounter,gAmin[0],gAmax[0]) fillGraph(hErrA[1][-1],gErrA[1],fcounter,gAmin[1],gAmax[1]) fillGraph(hErrA[2][-1],gErrA[2],fcounter,gAmin[2],gAmax[2]) hprojNum[0].append(fillGraph(hNum[0][-1],gNum[0],fcounter,gNummin[0],gNummax[0])) fitSlicesFromFile(gfNum[0],hNum[0][-1].GetName(),fcounter,rfile,gfNummin[0],gfNummax[0]) if args.make_lookup: if fcounter==0: dtheta=0 if len(titles)>1: dtheta=int(titles[1].split('_')[-1])-int(titles[0].split('_')[-1]) x,y=ROOT.Double(0),ROOT.Double(0) gfErrA3[1][-1]['Sigma'].GetPoint(0,x,y) offset=int(x) print 'offset:',x lookup_file.write(str(len(files))+' '+str(dtheta)+' '+str(gfErrA3[1][-1]['Sigma'].GetN()+offset)+' 1\n') lookup_file_num.write(str(len(files))+' '+str(dtheta)+' '+str(gfErrA3[1][-1]['Sigma'].GetN()+offset)+' 1\n') print 'writing lookup ',titles[fcounter] writeGraphToFile(gfErrA3[1][-1]['Sigma'],'0 '+titles[fcounter].replace('theta_',''),lookup_file) writeGraphToFile(gfErrA3[2][-1]['Sigma'],'1 '+titles[fcounter].replace('theta_',''),lookup_file) writeGraphToFile(gfNum[0][-1]['Mean'],'0 '+titles[fcounter].replace('theta_',''),lookup_file_num) writeGraphToFile(gfNum[0][-1]['Mean'],'1 '+titles[fcounter].replace('theta_',''),lookup_file_num) rfile.Close() lookup_file.close() lookup_file_num.close() gfile=ROOT.TFile(args.basepath+'/clerr_graph.root','RECREATE') leg=ROOT.TLegend(0,0,1,1) leg.SetFillColor(0) leg.SetTextSize(0.1) fcounter=-1 legp=[] for g in gErrT2[0]: fcounter+=1 legp.append(ROOT.TPolyMarker()) legp[-1].SetMarkerStyle(g['RMS'].GetMarkerStyle()) legp[-1].SetMarkerSize(3) legp[-1].SetMarkerColor(g['RMS'].GetMarkerColor()) leg.AddEntry(legp[-1],titles[fcounter],"P") cgErrT2.Update() cLegend=ROOT.TCanvas('cleg','Legend',1200,1000) leg.Draw() if args.pfile!='None': cLegend.SaveAs(args.pfile) drawgraphs(cgErrT2,gErrT2,gT2min,gT2max,leg) drawgraphs(cgErrT,gErrT,gTmin,gTmax,leg) drawgraphs(cgErr,gErr,gmin,gmax,leg) drawgraphs(cgErrR,gErrR,gRmin,gRmax,leg) drawgraphs(cgErrA,gErrA,gAmin,gAmax,leg) drawgraphs(cgErrAr,gErrAr,gArmin,gArmax,leg) drawgraphf(cgfErrAr,gfErrAr,gfArmin,gfArmax,leg) drawgraphs(cgErrA3,gErrA3,gA3min,gA3max,leg) drawgraphf(cgfErrA3,gfErrA3,gfA3min,gfA3max,leg) drawgraphs(cgPullUV,gPullUV,gPullUVmin,gPullUVmax,leg) drawgraphf(cgfPullUV,gfPullUV,gfPullUVmin,gfPullUVmax,leg) drawgraphs(cgNum,gNum,gNummin,gNummax,leg) drawgraphf(cgfNum,gfNum,gfNummin,gfNummax,leg) drawgraphs(cgPullxyzC,gPullxyzC,gPullxyzCmin,gPullxyzCmax,leg) drawgraphf(cgfPullxyzC,gfPullxyzC,gfPullxyzCmin,gfPullxyzCmax,leg) ''' gcounter=-1 for g in gAllErrA: gcounter+=1 cAllErrA.cd(2*gcounter+1) g['RMS'].Draw('AP') cAllErrA.cd(2*gcounter+2) g['Mean'].Draw('AP') cAllErrA.Update() gfile.Close() ''' UVFits=[] drawprojs(chprojUV,hprojUV,"UV",gfPullUV,UVFits) AxisFits=[] drawprojs(chprojA,hErrAProj,"Ax",gfErrAr,AxisFits) NumFits=[] drawprojs(chprojNum,hprojNum,'Num',gfNum,NumFits,'landau') XyzCFits=[] drawprojs(chprojXYZC,hprojPullxyzC,'XYZ C',gfPullxyzC,XyzCFits) if args.pfile!='None': cproj.cd(2) leg.Draw() drawprojsfile(cproj,hErrAProj,AxisFits) drawprojsfile(cproj,hprojUV,UVFits) drawprojsfile(cproj,hprojNum,NumFits) drawprojsfile(cproj,hprojPullxyzC,XyzCFits) if args.pfile!='None': cproj.SaveAs(args.pfile+']') c1=ROOT.TCanvas('cfitpars','Fit Pars') if args.phi: c1.Divide(5,3) else: c1.Divide(4,3) if args.fit: counter=-1 for i in coords: for g in graphs[i]: counter+=1 c1.cd(counter+1) g.Draw('AP') c1.Update() capars=ROOT.TCanvas('capars',' A Fit Pars',1000,700) capars.Divide(4,3) counter=-1 for i in coordsA: for g in graphsA[i]: counter+=1 capars.cd(counter+1) g.Draw('AP') capars.Update() if not args.quiet: u=raw_input("done?") if not args.new: outfile.Close()