import ROOT, glob, math from ROOT import std dir = "/nfs/hicran/data/tpc/fopi/2010/decoded" #dir = "/nfs/nas/user/sneubert/PANDA/" numfiles =100 preliminary = 1 outfile = ROOT.TFile("anaOut.root", "recreate") ROOT.gROOT.ProcessLine(".x rootlogon.C") ROOT.gROOT.ProcessLine('gSystem->Load("libPhysics")') ROOT.gROOT.ProcessLine('gStyle->SetPalette(1)') failed = ROOT.TH1D("Failed", "Failed Hits", 30,0,30) recoMom = ROOT.TH1D("recoMom", "Rec. Momenta", 500,900,1100) sampTimes = ROOT.TH1D("SamT", "Drift Time Distribution (samples)",511,0,511) occXY = ROOT.TH2D("OccXY", "Cluster XY occupancy created from cosmic tracks", 200,-15,15,200,-15,15) occZ = ROOT.TH1D("OccZ", "Cluster Z occupancy created from cosmic tracks", 200,0,75) files = glob.glob(dir + "/*.reco.root") files.sort() # -------------------- driftVel correction -------------------------------- #driftVelSim = 0.00055192 #(ArCo2 70/30 @ 250 V/cm) #driftVelReal = 0.00284 #(ArCo2 90/10 @ 350 V/cm) #velCorr = driftVelReal / driftVelSim # define cuts in Z cuts = (0,10,20,30,40,50,60) ampcuts = (0,400,800,1200,1600,2000,4000) ncuts = (0,1,2,3,5,7,10,15,20,30,40) # define offset in Z 1 mu s * 2.8 OFFSET = -2.8 resUs = dict([(i, ROOT.TH1D("StatsResX"+str(cuts[i]), "Cosmic Residuals U (Z')", 500,-1,1)) for i in range(len(cuts))]) pullsUs = dict([(i, ROOT.TH1D("StatsPullX"+str(cuts[i]), "Cosmic Pulls U (Z')", 100,-3,3)) for i in range(len(cuts))]) resVIDs = dict([(i, ROOT.TH1D("StatsResVID"+str(cuts[i]), "Cosmic Residuals V (X') single track events", 500,-1,1)) for i in range(len(cuts))]) resVs = dict([(i, ROOT.TH1D("StatsResV"+str(cuts[i]), "Cosmic Residuals V (X')", 500,-1,1)) for i in range(len(cuts))]) resVza = dict([(j,dict([(i,ROOT.TH1D("StatsResV"+str(cuts[i])+str(ampcuts[j]), "Cosmic Residuals V (X') z="+str(cuts[i])+" amp="+str(ampcuts[j]), 500,-1,1)) for i in range(len(cuts))])) for j in range(len(ampcuts))]) resVzn = dict([(j,dict([(i,ROOT.TH1D("StatsResVzn"+str(cuts[i])+str(ncuts[j]), "Cosmic Residuals V (X') z="+str(cuts[i])+" n="+str(ncuts[j]), 500,-1,1)) for i in range(len(cuts))])) for j in range(len(ncuts))]) resUzn = dict([(j,dict([(i,ROOT.TH1D("StatsResUzn"+str(cuts[i])+str(ncuts[j]), "Cosmic Residuals U (Z') z="+str(cuts[i])+" n="+str(ncuts[j]), 500,-1,1)) for i in range(len(cuts))])) for j in range(len(ncuts))]) gresVa = dict([(i, ROOT.TGraph(len(ampcuts)-1)) for i in range(len(cuts))]) gresVz = dict([(j, ROOT.TGraph(len(cuts)-1)) for j in range(len(ampcuts))]) pullsVs = dict([(i, ROOT.TH1D("StatsPullV"+str(cuts[i]), "Cosmic Pulls V (X')", 100,-3,3)) for i in range(len(cuts))]) pullsVa = dict([(i, ROOT.TH1D("StatsPullV"+str(ampcuts[i]), "Cosmic Pulls V (X') amp="+str(ampcuts[i]), 100,-3,3)) for i in range(len(ampcuts))]) cl2Ds = dict([(i, ROOT.TH1D("cl2Ds"+str(cuts[i]), "Cluster Size Distributions as Function of Z", 50,0,50)) for i in range(len(cuts))]) #resUVs = dict([(i, ROOT.TH1D("StatsResXY"+str(cuts[i]), # "Cosmic Residuals X'Y'", 500,-1,1)) for i in range(len(cuts))]) resVsUs = dict([(i, ROOT.TH2D("StatsResXsYs"+str(cuts[i]), "Cosmic Residuals U(Z') vs V(X')", 500,-0.5,0.5, 500,-0.5,0.5)) for i in range(len(cuts))]) clSize2D = ROOT.TH1D("clSize2D", "2D Cluster Size Distribution", 100,0,100) clSize2D.SetFillColor(ROOT.kBlack) clSize = ROOT.TH1D("clSize", "Cluster Size Distribution", 100,0,100) clSize.SetLineColor(ROOT.kRed+3) clSizevsDrift = ROOT.TH2D("clnSizevsDrift", "Cluster Size vs. Drift", 100,0,75,60,0,60) chi2prob = ROOT.TH1D("chi2prob", "ChiSqu Probability Distribution", 10000,0,1) chi2probID = ROOT.TH1D("chi2probID", "ChiSqu Probability Distribution", 10000,0,1) chi2probID.SetLineColor(ROOT.kRed+2) chi2raw = ROOT.TH1D("chi2raw", "ChiSqu2 / NDF", 300,0,10) phiDist = ROOT.TH1D("phiDist", "Phi Distribution of fitted Tracks", 300,-3.5,3.5) thetaDist = ROOT.TH1D("thetaDist", "Theta Distribution of fitted Tracks", 300,-3.5,3.5) sigxDist = ROOT.TH1D("sigxDist", "Cluster SigmaX Distribution", 500, 0, 0.3) sigyDist = ROOT.TH1D("sigyDist", "Cluster SigmaY Distribution", 500, 0, 0.3) sigXvsClSize = ROOT.TH2D("sigXClSize", "Cluster SigmaX vs Cluster Size", 80,0,40,500,0,0.3) sigXvsDrift = ROOT.TH2D("sigXvsZ", "Cluster SigmaX vs Drift Length Z", 400,0,65,400,0,3000) sigXvsDrift.GetYaxis().SetTitle("#sigma_X (#mu m)") nTracksDist = ROOT.TH1D("nTracksDist", "Distribution of track multiplicity",10,0,10) pullAllV = ROOT.TH1D("pullAllV", "V (X') Pull Distribution", 200,-3,3) DIFFX = 0.0227*10000 #(mu m/ mu s) diffT = ROOT.TF1("diffT","[0]*TMath::Sqrt(x)",0,75) diffT.SetNpx(1000) diffT.SetLineWidth(2) diffT.SetLineStyle(4) diffT.SetParameter(0,DIFFX) resDX = ROOT.TF1("resDX", "TMath::Sqrt([0]/12 + [1]*x)", 0, 75) resDX.SetNpx(1000) resDX.SetParameter(0,3000**2) resDX.SetParameter(1,DIFFX**2) resDX.SetLineWidth(2) chi2func = ROOT.TF1("meh", "[0]*x*TMath::Exp(-2*x)",0,10) chi2func.SetParameter(0,2500) #for file in files : for f in range(numfiles) : if f==62 : continue file = files[f] #file = "/nfs/hicran/data/tpc/fopi/2010/decoded/runC_1735.reco.root" print(file) Rfile = ROOT.TFile(file, "read") if (not Rfile.IsOpen()) or Rfile.IsZombie() : continue tree = Rfile.Get("cbmsim") tree.SetBranchStatus("*", 0) #tree.SetBranchStatus("PndTpcSLResiduals.*", 1) tree.SetBranchStatus("TrackFitStat.*", 1) tree.SetBranchStatus("PndTpcCluster.*", 1) #tree.SetBranchStatus("TrackPostFit.*", 1) for e in tree : #for trk in e.TrackPostFit : # trk.Print() # cand = trk.getCand() # #cand.Print() # res = std.vector('double')() # trk.getResiduals(2,0,0,res) # # LOOP OVER CANDIDATE AND PLOT PULLS # hitIDs = cand.getHitIDs # #for id in hitIDs : # # GET CLUSTER, SIGMA AND RESIDUAL for cl in e.PndTpcCluster : sig = cl.sig() sigxDist.Fill(sig.X()) sigyDist.Fill(sig.Y()) size = cl.size() pos = cl.pos() sigXvsClSize.Fill(size, sig.X()) sigXvsDrift.Fill(pos.Z(), sig.X()*10000) clSizevsDrift.Fill(pos.Z(),size) nTracks = e.TrackFitStat.GetEntriesFast() if nTracks > 0 : nTracksDist.Fill(nTracks) # if nTracks > 1 : # continue for tfs in e.TrackFitStat : chi2 = tfs.getChi2() #redChi2 = tfs.getRedChi2() NDF = tfs.getNDF() NDIM=4 numHits = tfs.GetHitPositionsZ().size() if numHits < 10 : print("Short track. Skipping") continue chi2Prob = ROOT.TMath.Prob(chi2, NDF) chi2prob.Fill(chi2Prob) chi2raw.Fill(chi2/(NDF)) if nTracks == 1: chi2probID.Fill(chi2Prob) mom = tfs.GetMom() mom.SetMag(1.) phiDist.Fill(mom.Phi()) thetaDist.Fill(mom.Theta()) #build orthogonal vector to fix plane - res x mom vecX = ROOT.TVector3(1,0,0) vecZ = ROOT.TVector3(0,0,1) #define new coordinate system for this track u = mom.Cross(vecX) v = mom.Cross(u) for p in range(numHits) : z = tfs.GetHitPositionsZ().at(p) x = tfs.GetHitPositionsX().at(p) y = tfs.GetHitPositionsY().at(p) amp = tfs.GetAmps().at(p) xyRad = math.sqrt(x**2 + y**2) resX = tfs.GetResX().at(p) resY = tfs.GetResY().at(p) resZ = tfs.GetResZ().at(p) res = ROOT.TVector3(resX,resY,resZ) sigX = tfs.GetSigX().at(p) sigY = tfs.GetSigY().at(p) sigZ = tfs.GetSigZ().at(p) sig = ROOT.TVector3(sigX,sigY,sigZ) resU = res.Dot(u) resV = res.Dot(v) sigU = sig.Dot(u) sigV = sig.Dot(v) cl2Dsize = tfs.Get2DClSizes().at(p) clTotSize = tfs.GetClSizes().at(p) clSize2D.Fill(cl2Dsize) clSize.Fill(clTotSize) for i in range(len(cuts)) : if z < cuts[i] : cl2Ds[i-1].Fill(cl2Dsize) resUs[i-1].Fill(resU) resVs[i-1].Fill(resV) pullsUs[i-1].Fill(resU/sigU) pullsVs[i-1].Fill(resV/sigV) #resUVs[i-1].Fill(math.sqrt(resX**2 + resY**2)) resVsUs[i-1].Fill(resV, resU) #compare with these filters: if chi2Prob > 1E-4 : resVIDs[i-1].Fill(resV) #if xyRad > 8 and xyRad < 12: #if numHits > 20 : for j in range(len(ncuts)) : if (cl2Dsize < ncuts[j]) : resVzn[j-1][i-1].Fill(resV) resUzn[j-1][i-1].Fill(resU) break for j in range(len(ampcuts)) : if ( amp < ampcuts[j]) and (cl2Dsize > 1) : resVza[j-1][i-1].Fill(resV) break break for i in range(len(ampcuts)) : if (amp < ampcuts[i]) and (cl2Dsize > 1) and ( cuts[5]< z < cuts[6]): pullsVa[i-1].Fill(resV/sigV) break x = tfs.GetHitPositionsX().at(p) y = tfs.GetHitPositionsY().at(p) occXY.Fill(x,y) occZ.Fill(z) recoMom.Fill(tfs.GetP()) outfile.cd() c1 = ROOT.TCanvas() c1.Divide(3,2) for i in range(6) : c1.cd(i+1) resVsUs[i].GetXaxis().SetTitle("Residual V (cm)") resVsUs[i].GetYaxis().SetTitle("Residual U (cm)") resVsUs[i].Draw("COLZ") resVsUs[i].Write() c2 = ROOT.TCanvas() occXY.Draw("COLZ") occXY.Write() c3 = ROOT.TCanvas() occZ.Draw() occZ.Write() diffV = ROOT.TGraph(6) diffV.SetName("diffX") diffV.SetTitle("V(X') Resolution as function of Drift Length") diffVID = ROOT.TGraph(6) diffVID.SetMarkerColor(ROOT.kRed+2) bckgrShare = ROOT.TGraph(6) bckgrShare.SetName("bckgShare") bckgrShare.SetTitle("Share of background (from fits)") c4 = ROOT.TCanvas() c4.Divide(3,2) for i in range(6) : c4.cd(i+1) resVs[i].SetFillColor(ROOT.kAzure-8) resVs[i].GetXaxis().SetTitle("Residual V (cm)") resVs[i].Draw() testfit = ROOT.TF1("testfitX"+str(1),"gaus",-1,1) resVs[i].Fit(testfit, "N+", "", -1,1) fit = ROOT.TF1("fitfuncX"+str(1),"gaus + gaus(3)",-1,1) fit.SetNpx(1000) fit.SetParameter(0,testfit.GetParameter(0)) fit.SetParLimits(0,testfit.GetParameter(0)*0.2, testfit.GetParameter(0)*5) fit.SetParameter(1,testfit.GetParameter(1)) fit.SetParLimits(1,testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(2,testfit.GetParameter(2)) fit.SetParLimits(2,testfit.GetParameter(2)*0.1, testfit.GetParameter(2)*20) fit.SetParameter(3, 100) fit.SetParLimits(3, 0, testfit.GetParameter(0)/4.) fit.SetParameter(4, 0) fit.SetParLimits(4, testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(5, testfit.GetParameter(2)) fit.SetParLimits(5, testfit.GetParameter(2), testfit.GetParameter(2)*30) resVs[i].Fit(fit, "+", "", -1,1) diffV.SetPoint(i,cuts[i]+5,fit.GetParameter(2)*10000) # preliminary if preliminary : histo = resVs[i]; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray resVs[i].Write() #calculate ratio of central and background integrals ratio = fit.GetParameter(3)*fit.GetParameter(5)/(fit.GetParameter(0)*fit.GetParameter(2)) bckgrShare.SetPoint(i,cuts[i]+5,ratio) c5 = ROOT.TCanvas() c5.Divide(3,2) for i in range(6) : c5.cd(i+1) ROOT.gPad.SetLogy(1) cl2Ds[i].Draw() c6 = ROOT.TCanvas() c6.Divide(3,2) for i in range(6) : c6.cd(i+1) resUs[i].SetFillColor(ROOT.kAzure-8) resUs[i].GetXaxis().SetTitle("Residual U (cm)") resUs[i].Draw() testfit = ROOT.TF1("testfitY"+str(1),"gaus",-1,1) resUs[i].Fit(testfit, "N+", "", -1,1) fit = ROOT.TF1("fitfuncY"+str(1),"gaus + gaus(3)",-1,1) fit.SetNpx(1000) fit.SetParameter(0,testfit.GetParameter(0)) fit.SetParLimits(0,testfit.GetParameter(0)*0.5, testfit.GetParameter(0)*2) fit.SetParameter(1,testfit.GetParameter(1)) fit.SetParLimits(1,testfit.GetParameter(1)-0.03, testfit.GetParameter(1)+0.03) fit.SetParameter(2,testfit.GetParameter(2)) fit.SetParLimits(2,testfit.GetParameter(2)*0.2, testfit.GetParameter(2)) fit.SetParameter(3, 10) fit.SetParLimits(3, 0, testfit.GetParameter(0)/1.5) fit.SetParameter(4, 0) fit.SetParLimits(4, testfit.GetParameter(1)-0.03, testfit.GetParameter(1)+0.03) fit.SetParameter(5, testfit.GetParameter(2)*3) fit.SetParLimits(5, testfit.GetParameter(2)*1.4, testfit.GetParameter(2)*20) resUs[i].Fit(fit, "+", "", -1,1) # preliminary if preliminary : histo = resUs[i]; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray resUs[i].Write() c7 = ROOT.TCanvas() c7.Divide(3,2) for i in range(6) : c7.cd(i+1) resVIDs[i].SetFillColor(ROOT.kAzure-8) resVIDs[i].GetXaxis().SetTitle("Residual V (cm)") resVIDs[i].Draw() testfit = ROOT.TF1("testfitX"+str(1),"gaus",-1,1) resVIDs[i].Fit(testfit, "N+", "", -1,1) fit = ROOT.TF1("fitfuncX"+str(1),"gaus + gaus(3)",-1,1) fit.SetNpx(1000) fit.SetParameter(0,testfit.GetParameter(0)) fit.SetParLimits(0,testfit.GetParameter(0)*0.2, testfit.GetParameter(0)*5) fit.SetParameter(1,testfit.GetParameter(1)) fit.SetParLimits(1,testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(2,testfit.GetParameter(2)) fit.SetParLimits(2,testfit.GetParameter(2)*0.1, testfit.GetParameter(2)*20) fit.SetParameter(3, 10) fit.SetParLimits(3, 0, testfit.GetParameter(0)/3.) fit.SetParameter(4, 0) fit.SetParLimits(4, testfit.GetParameter(1)-0.8, testfit.GetParameter(1)+0.8) fit.SetParameter(5, testfit.GetParameter(2)*2) fit.SetParLimits(5, testfit.GetParameter(2), testfit.GetParameter(2)*50) resVIDs[i].Fit(fit, "+", "", -1,1) # preliminary if preliminary : histo = resVIDs[i]; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray diffVID.SetPoint(i,cuts[i]+5,fit.GetParameter(2)*10000) resVIDs[i].Write() c8 = ROOT.TCanvas() bckgrShare.SetMarkerStyle(21) bckgrShare.GetYaxis().SetRangeUser(0,1) bckgrShare.GetXaxis().SetTitle("Drift Length Z (cm)") bckgrShare.Draw("ALP") bckgrShare.Write() c9 = ROOT.TCanvas() diffV.SetMarkerStyle(21) diffV.Draw("ALP") diffV.GetXaxis().SetTitle("Drift Length Z (cm)") diffV.GetYaxis().SetTitle("X Resolution (#mu m)") diffV.GetYaxis().SetRangeUser(0,2000) diffVID.SetMarkerStyle(21) diffVID.Draw("LP") diffV.Write() diffVID.Write() diffT.Draw("same") # preliminary if preliminary : histo = diffV; xaxis = histo.GetXaxis() yaxis = histo.GetYaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() ycenter=yaxis.GetBinCenter(int(yaxis.GetNbins()*0.5)) yrange=yaxis.GetXmax()-yaxis.GetXmin() x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c10 = ROOT.TCanvas() c10.Divide(2,1) c10.cd(1) clSize2D.Draw() clSize.Draw("same") # preliminary if preliminary : histo = clSize2D xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c10.cd(2).SetLogz(1) clSizevsDrift.Draw("colz") # preliminary if preliminary : histo = clSizevsDrift xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c11 = ROOT.TCanvas() c11.cd() ROOT.gPad.SetLogy(1) chi2prob.Draw() chi2probID.Draw("same") # preliminary if preliminary : histo = chi2prob; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=math.log(ycenter-0.3*yrange) prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c12 = ROOT.TCanvas() chi2raw.Draw() chi2func.SetNpx(1000) # preliminary if preliminary : histo = chi2raw; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray #chi2func.Draw("same") c13 = ROOT.TCanvas() c13.Divide(2,1) c13.cd(1) phiDist.Draw() c13.cd(2) thetaDist.Draw() c14 = ROOT.TCanvas() c14.Divide(2,1) c14.cd(1) sigxDist.Draw() # preliminary if preliminary : histo = sigxDist; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c14.cd(2) sigyDist.Draw() # preliminary if preliminary : histo = sigyDist; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c15 = ROOT.TCanvas() c15.Divide(2,1) c15.cd(1).SetLogz() sigXvsClSize.Draw("colz") # preliminary if preliminary : histo = sigXvsClSize; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c15.cd(2).SetLogz() sigXvsDrift.Draw("colz") diffT.Draw("same") #resDX.Draw("same") # preliminary if preliminary : histo = sigXvsDrift; xaxis = histo.GetXaxis() yaxis = histo.GetYaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() ycenter=yaxis.GetBinCenter(int(yaxis.GetNbins()*0.5)) yrange=yaxis.GetXmax()-yaxis.GetXmin() x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray c16 = ROOT.TCanvas() nTracksDist.Draw() c17 = ROOT.TCanvas() ROOT.gStyle.SetOptFit(1111) c17.Divide(3,2) for i in range(6) : c17.cd(i+1) pullsVs[i].SetFillColor(ROOT.kAzure-8) pullsVs[i].GetXaxis().SetTitle("Pulls V") pullsVs[i].Draw() testfit = ROOT.TF1("testfitX"+str(1),"gaus",-3,3) pullsVs[i].Fit(testfit, "N+", "", -3,3) fit = ROOT.TF1("fitfuncX"+str(1),"gaus + gaus(3)",-3,3) fit.SetNpx(1000) fit.SetParameter(0,testfit.GetParameter(0)) fit.SetParLimits(0,testfit.GetParameter(0)*0.2, testfit.GetParameter(0)*5) fit.SetParameter(1,testfit.GetParameter(1)) fit.SetParLimits(1,testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(2,testfit.GetParameter(2)) fit.SetParLimits(2,testfit.GetParameter(2)*0.1, testfit.GetParameter(2)*20) fit.SetParameter(3, 100) fit.SetParLimits(3, 0, testfit.GetParameter(0)/4.) fit.SetParameter(4, 0) fit.SetParLimits(4, testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(5, testfit.GetParameter(2)) fit.SetParLimits(5, testfit.GetParameter(2), testfit.GetParameter(2)*30) pullsVs[i].Fit(fit, "+", "", -3,3) # preliminary if preliminary : histo = pullsVs[i]; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray pullsVs[i].Write() c18 = ROOT.TCanvas() ROOT.gStyle.SetOptFit(1111) c18.Divide(3,2) for i in range(6) : c18.cd(i+1) pullsVa[i].SetFillColor(ROOT.kAzure-8) pullsVa[i].GetXaxis().SetTitle("Pulls V") pullsVa[i].Draw() testfit = ROOT.TF1("testfitX"+str(1),"gaus",-3,3) pullsVa[i].Fit(testfit, "N+", "", -3,3) fit = ROOT.TF1("fitfuncX"+str(1),"gaus + gaus(3)",-3,3) fit.SetNpx(1000) fit.SetParameter(0,testfit.GetParameter(0)) # fit.SetParLimits(0,testfit.GetParameter(0)*0.2, testfit.GetParameter(0)*5) fit.SetParameter(1,testfit.GetParameter(1)) #fit.SetParLimits(1,testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(2,testfit.GetParameter(2)) #fit.SetParLimits(2,testfit.GetParameter(2)*0.1, testfit.GetParameter(2)*20) fit.SetParameter(3, 100) #fit.SetParLimits(3, 0, testfit.GetParameter(0)/4.) fit.SetParameter(4, 0) #fit.SetParLimits(4, testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(5, testfit.GetParameter(2)) #fit.SetParLimits(5, testfit.GetParameter(2), testfit.GetParameter(2)*30) pullsVa[i].Fit(fit, "+", "", -3,3) # preliminary if preliminary : histo = pullsVa[i]; xaxis = histo.GetXaxis() xcenter=xaxis.GetBinCenter(int(xaxis.GetNbins()*0.5)) xrange=xaxis.GetXmax()-xaxis.GetXmin() yrange=-histo.GetMinimum()+histo.GetMaximum() ycenter=(histo.GetMinimum()+histo.GetMaximum())*0.5 x=xcenter-0.3*xrange y=ycenter-0.3*yrange prelim= ROOT.TLatex(x,y,"preliminary") prelim.SetTextColor(ROOT.kGray) prelim.SetTextSize(0.1) prelim.SetTextAngle(20) prelim.DrawLatex(x,y,"preliminary") # end preliminray pullsVa[i].Write() c19 = ROOT.TCanvas() ROOT.gStyle.SetOptFit(1111) c19.Divide(6,6) # loop over z-bins : i # loop over amp-bins : j for i in range(6) : for j in range(6) : c19.cd((i)*6+j+1) resVza[j][i].SetFillColor(ROOT.kAzure-8) resVza[j][i].GetXaxis().SetTitle("Residiual V") resVza[j][i].Draw() testfit = ROOT.TF1("testfitX"+str(1),"gaus",-0.5,0.5) resVza[j][i].Fit(testfit, "N+", "", -0.5,0.5) fit = ROOT.TF1("fitfuncX"+str(1),"gaus + gaus(3)",-0.5,0.5) fit.SetNpx(1000) fit.SetParameter(0,testfit.GetParameter(0)) fit.SetParLimits(0,testfit.GetParameter(0)*0.2, testfit.GetParameter(0)*5) fit.SetParameter(1,testfit.GetParameter(1)) fit.SetParLimits(1,testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(2,testfit.GetParameter(2)) fit.SetParLimits(2,testfit.GetParameter(2)*0.1, testfit.GetParameter(2)*20) fit.SetParameter(3, 100) fit.SetParLimits(3, 0, testfit.GetParameter(0)/4.) fit.SetParameter(4, 0) fit.SetParLimits(4, testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) fit.SetParameter(5, testfit.GetParameter(2)) fit.SetParLimits(5, testfit.GetParameter(2), testfit.GetParameter(2)*30) resVza[j][i].Fit(fit, "+", "", -0.5,0.5) gresVz[j].SetPoint(i,cuts[i]+5,fit.GetParameter(2)); gresVa[i].SetPoint(j,ampcuts[j],fit.GetParameter(2)); #c17.Update() c20 = ROOT.TCanvas() c20.Divide(2,1) mgresVz = ROOT.TMultiGraph() mgresVa = ROOT.TMultiGraph() for i in range(6) : mgresVz.Add(gresVz[i]) mgresVa.Add(gresVa[i]) c20.cd(1) mgresVz.Draw("ALP") c20.cd(2) mgresVa.Draw("ALP") c21 = ROOT.TCanvas() ROOT.gStyle.SetOptFit(1111) c21.Divide(len(ncuts),6) # loop over z-bins : i # loop over amp-bins : j for i in range(6) : for j in range(len(ncuts)) : c21.cd((i)*len(ncuts)+j+1) print i print j print resVzn[j][i].GetTitle() resVzn[j][i].SetFillColor(ROOT.kAzure-8) resVzn[j][i].GetXaxis().SetTitle("Residiual V") resVzn[j][i].Draw() # testfit = ROOT.TF1("testfitX"+str(1),"gaus",-0.5,0.5) #resVzn[j][i].Fit(testfit, "N+", "", -0.5,0.5) # # fit = ROOT.TF1("fitfuncX"+str(1),"gaus + gaus(3)",-0.5,0.5) # fit.SetNpx(1000) # fit.SetParameter(0,testfit.GetParameter(0)) # fit.SetParLimits(0,testfit.GetParameter(0)*0.2, testfit.GetParameter(0)*5) # fit.SetParameter(1,testfit.GetParameter(1)) # fit.SetParLimits(1,testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) # fit.SetParameter(2,testfit.GetParameter(2)) # fit.SetParLimits(2,testfit.GetParameter(2)*0.1, testfit.GetParameter(2)*20) # fit.SetParameter(3, 100) # fit.SetParLimits(3, 0, testfit.GetParameter(0)/4.) # fit.SetParameter(4, 0) # fit.SetParLimits(4, testfit.GetParameter(1)-0.08, testfit.GetParameter(1)+0.08) # fit.SetParameter(5, testfit.GetParameter(2)) # fit.SetParLimits(5, testfit.GetParameter(2), testfit.GetParameter(2)*30) # resVzn[j][i].Fit(fit, "+", "", -0.5,0.5) # gresVz[j].SetPoint(i,cuts[i]+5,fit.GetParameter(2)); # gresVa[i].SetPoint(j,ampcuts[j],fit.GetParameter(2)); c22 = ROOT.TCanvas() ROOT.gStyle.SetOptFit(1111) c22.Divide(len(ncuts),6) # loop over z-bins : i # loop over amp-bins : j for i in range(6) : for j in range(len(ncuts)) : c22.cd((i)*len(ncuts)+j+1) resUzn[j][i].SetFillColor(ROOT.kAzure-8) resUzn[j][i].GetXaxis().SetTitle("Residiual U") resUzn[j][i].Draw() input() outfile.Close()