def fillHists(tree,hdict,hdict2=None,hdict3=None): event=0 for e in tree: event+=1 if args.nevts!=-1 and event>args.nevts: break if event%100==0: print event for tr in e.CutCosmics: theta=tr.GetTheta() phi=tr.GetPhi() mom=tr.GetMom() mcmom=tr.GetMCmom() mctheta=mcmom.Theta()*ROOT.TMath.RadToDeg() mcphi=mcmom.Phi()*ROOT.TMath.RadToDeg() mz=tr.GetMeanZ() hdict['theta'].Fill(theta-mctheta) hdict['mom'].Fill(mom.Mag()-mcmom.Mag()) hdict['phi'].Fill(phi-mcphi) hdict['phi2'].Fill(phi) hdict['phi z'].Fill(mz,phi-mcphi) if hdict2!=None: theta=tr.GetzAxisDir().Theta()*ROOT.TMath.RadToDeg() phi=tr.GetzAxisDir().Phi()*ROOT.TMath.RadToDeg() hdict2['theta'].Fill(theta-mctheta) hdict2['mom'].Fill(mom.Mag()-mcmom.Mag()) hdict2['phi'].Fill(phi-mcphi) hdict2['phi2'].Fill(phi) hdict2['phi z'].Fill(mz,phi-mcphi) if hdict3!=None: theta=tr.GetFirstMCDir().Theta()*ROOT.TMath.RadToDeg() phi=tr.GetFirstMCDir().Phi()*ROOT.TMath.RadToDeg() hdict3['theta'].Fill(theta-mctheta) hdict3['mom'].Fill(mom.Mag()-mcmom.Mag()) hdict3['phi'].Fill(phi-mcphi) hdict3['phi2'].Fill(phi) hdict3['phi z'].Fill(mz,phi-mcphi) return def print_info(d1,d2,var): print args.titles[0]+' rms',var,':',d1[var].GetRMS() print args.titles[1]+' rms',var,':',d2[var].GetRMS() print 'rms diff',var,':',100*(d1[var].GetRMS()-d2[var].GetRMS())/d2[var].GetRMS(),'%' print args.titles[0]+' mean',var,':',d1[var].GetMean() print args.titles[1]+' mean',var,':',d2[var].GetMean() print 'mean diff',var,':',100*(d1[var].GetMean()-d2[var].GetMean())/d2[var].GetMean(),'%' print '' def plotHists(hists,c1,c2,c3): maxi={} for f in hists: for hn in f: h=f[hn] if type(h)==ROOT.TH1D: h.Scale(1/h.GetEntries()) h.SetStats(0) maxi[hn]=max(maxi.get(hn,0),h.GetMaximum()) ''' for hn in hists: h=hists[hn] if type(h)==ROOT.TH1D: h.Scale(1/h.GetEntries()) h.SetStats(0) ''' opt='' counter=-1 canv=[c1,c2,c3] legends.append(ROOT.TLegend(0,0,1,1)) for f in hists: counter+=1 c1.cd(1) f['theta'].SetMaximum(maxi['theta']*1.1) f['theta'].Draw(opt) f['theta'].GetXaxis().SetTitle('Theta Residual (#circ)') f['theta'].SetLineColor(colors[counter]) legends[-1].AddEntry(f['theta'],args.titles[counter],'l') c1.cd(2) f['phi'].SetMaximum(maxi['phi']*1.1) f['phi'].Draw(opt) f['phi'].GetXaxis().SetTitle('Theta Residual (#circ)') f['phi'].SetLineColor(colors[counter]) c1.cd(3) f['mom'].SetMaximum(maxi['mom']*1.1) f['mom'].Draw(opt) f['mom'].GetXaxis().SetTitle('Theta Residual (#circ)') f['mom'].SetLineColor(colors[counter]) c1.Update() c3.cd(1) f['phi2'].SetLineColor(colors[counter]) f['phi2'].Draw(opt) opt='same' c1.Update() c3.Update() import ROOT, argparse,sys,math sys.path.append('$PANDAPATH/macro/tpc/FOPI/mberger') from functions import * from cosmic_histos import * parser=argparse.ArgumentParser(description='plot track residuals (theta, phi, mom)') parser.add_argument('filenames',help='the reco files',default=['None'],type=str,nargs='+') parser.add_argument('--titles',help='the titles',default=['1','2'],type=str,nargs='+') parser.add_argument('--nevts',help='number of events',type=int,default=-1) parser.add_argument('--pfile',help='the file to save all',type=str,default='None') parser.add_argument("--hlp",help="display help",action='store_true') args=parser.parse_args() if args.hlp: parser.print_help() exit(0) ROOT.gROOT.ProcessLine(".x rootlogon.C") ROOT.gROOT.ProcessLine('gROOT->SetStyle("Plain")') set_palette() mcangle=90 cv1=ROOT.TCanvas('c1','First Point',500,1000) cv1.Divide(1,3) cv2=ROOT.TCanvas('c2','First Point',500,1000) cv2.Divide(1,3) cv3=ROOT.TCanvas('c3','First Point',500,1000) cv3.Divide(1,3) cv4=ROOT.TCanvas('c4','zAxis',500,1000) cv4.Divide(1,3) cv5=ROOT.TCanvas('c5','zAxis',500,100) cv5.Divide(1,3) cv6=ROOT.TCanvas('c6','zAxis',500,1000) cv6.Divide(1,3) cv7=ROOT.TCanvas('c7','MC Point',500,1000) cv7.Divide(1,3) cv8=ROOT.TCanvas('c8','MC Point',500,1000) cv8.Divide(1,3) cv9=ROOT.TCanvas('c9','MC Point',500,1000) cv9.Divide(1,3) rfile=[] trees=[] hists=[] hists2=[] hists3=[] fcounter=-1 for f in args.filenames: fcounter+=1 print f rfile.append(ROOT.TFile(f,'read')) trees.append(rfile[-1].Get('cbmsim')) trees[-1].SetBranchStatus('*',0) trees[-1].SetBranchStatus('CutCosmics.*',1) hists.append({}) hists2.append({}) hists3.append({}) hists[-1]['theta']=ROOT.TH1D('htheta'+str(fcounter),'Theta-MCTheta '+args.titles[0],100,-1,1) hists[-1]['phi']=ROOT.TH1D('hphi'+str(fcounter),'Phi-MCPhi '+args.titles[0],200,-4,4) hists[-1]['phi2']=ROOT.TH1D('hphi2'+str(fcounter),'Phi '+args.titles[0],360*2,-180,180) hists[-1]['mom']=ROOT.TH1D('hmom'+str(fcounter),'Mom-MCmom '+args.titles[0],1000,-1,1) hists[-1]['phi z']=ROOT.TH2D('hphiZ'+str(fcounter),'Phi vs Z '+args.titles[0],80,0,80,100,-1,1) for h in hists[-1]: hists[-1][h].SetTitle(hists[-1][h].GetTitle()+' First Point') hists2[-1][h]=hists[-1][h].Clone(hists[-1][h].GetName()+'zAxis') hists2[-1][h].SetTitle(hists2[-1][h].GetTitle()+' Z-Axis') hists3[-1][h]=hists[-1][h].Clone(hists[-1][h].GetName()+'firstMC') hists3[-1][h].SetTitle(hists3[-1][h].GetTitle()+' first MC') fillHists(trees[-1],hists[-1],hists2[-1],hists3[-1]) if len(args.filenames)<10: colors=getColors2() else: colors=getColors() global legends legends=[] print 'First track point' plotHists(hists,cv1,cv2,cv3) print 'Z-Axis' plotHists(hists2,cv4,cv5,cv6) print 'first MC point' plotHists(hists3,cv7,cv8,cv9) cleg=ROOT.TCanvas('leg','leg',500,1000) legends[-1].Draw() if args.pfile!='None': cleg.SaveAs(args.pfile+'(') cv1.SaveAs(args.pfile) #cv2.SaveAs(args.pfile) #cv3.SaveAs(args.pfile) cv4.SaveAs(args.pfile) #cv5.SaveAs(args.pfile) #cv6.SaveAs(args.pfile) cv7.SaveAs(args.pfile) #cv8.SaveAs(args.pfile) #cv9.SaveAs(args.pfile+')') cv9.SaveAs(args.pfile+']') u=raw_input("done?")