FSROOT Advanced Tutorial#
Note
This tutorial is a PyROOT port of GlueX 2022 Tutorial Session 2d by Malte Albrecht. You must be on the JLab farm to run this tutorial to access the ROOT data files. Translating form c++ FSROOT to PyROOT FSROOT is as simple as changing how attributes are accessed (i.e. -> to .). FSROOT syntax is pythonic already but now you have access to python string manipulation.
from pyamptools import atiSetup
import ROOT
import os
atiSetup.setup(globals(), use_fsroot=True)
Welcome to JupyROOT 6.28/06
atiSetup| jupyter-book called python3.9
------------------------------------------------
atiSetup| MPI is disabled
atiSetup| GPU is disabled
------------------------------------------------
atiSetup| Loading library libIUAmpTools.so ............ ON
atiSetup| Loading library libAmpTools.so .............. ON
atiSetup| Loading library libAmpPlotter.so ............
ON
atiSetup| Loading library libAmpsDataIO.so ............ ON
atiSetup| Loading library libFSRoot.so ................ ON
atiSetup| Loading library libAmpsGen.so ............... OFF
------------------------------------------------
------------------------------------------------
atiSetup| Saved aliases found in /d/grid17/ln16/PyAmpTools/src/pyamptools/.aliases.txt, attempting to load...
atiSetup| minor warning: Unable to alias omegapiAngles - doesn't exist under ROOT namespace
atiSetup| minor warning: Unable to alias URConfig - doesn't exist under ROOT namespace
atiSetup| minor warning: Unable to alias URtypes - doesn't exist under ROOT namespace
atiSetup| minor warning: Unable to alias FSFitFunctions - doesn't exist under ROOT namespace
atiSetup| minor warning: Unable to alias FSFitPrivate - doesn't exist under ROOT namespace
(False, False, 0)
General Specifications
FND = "/work/halld/gluex_workshop_data/tutorial_2022/session2d/skim/tree_pi0eta__B4_M17_M7_DATA_sp17_*_GENERAL_SKIM.root"
NT = "ntFSGlueX_MODECODE"
bggen = False
Define event selections in FSROOT
FSCut.defineCut("unusedE","EnUnusedSh<0.1")
FSCut.defineCut("unusedTracks","NumUnusedTracks<1")
FSCut.defineCut("zProton","ProdVz>=52&&ProdVz<=78")
FSCut.defineCut("protMom","MOMENTUM([p+])>=0.3")
FSCut.defineCut("cet0103","OR(abs(-1*MASS2(GLUEXTARGET,-[p+])-0.2)<0.1)")
FSCut.defineCut("e8288","(EnPB>8.2&&EnPB<8.8)")
FSCut.defineCut("chi2","Chi2DOF<3.3","Chi2DOF>10.0&&Chi2DOF<15.0")
FSCut.defineCut("photFiducialA","(acos(COSINE([eta]a))*180/3.141>2.5 && acos(COSINE([eta]a))*180/3.141<10.3) || (acos(COSINE([eta]a))*180/3.141>11.9)")
FSCut.defineCut("photFiducialB","(acos(COSINE([eta]b))*180/3.141>2.5 && acos(COSINE([eta]b))*180/3.141<10.3) || (acos(COSINE([eta]b))*180/3.141>11.9)")
FSCut.defineCut("photFiducialC","(acos(COSINE([pi0]a))*180/3.141>2.5 && acos(COSINE([pi0]a))*180/3.141<10.3) || (acos(COSINE([pi0]a))*180/3.141>11.9)")
FSCut.defineCut("photFiducialD","(acos(COSINE([pi0]b))*180/3.141>2.5 && acos(COSINE([pi0]b))*180/3.141<10.3) || (acos(COSINE([pi0]b))*180/3.141>11.9)")
FSCut.defineCut("rejectOmega","!((MASS([pi0]a,[eta]a)<0.15 && MASS([pi0]b,[eta]b)<0.15) || (MASS([pi0]a,[eta]b)<0.15 && MASS([pi0]b,[eta]a)<0.15) || (MASS([pi0]a,[eta]a)<0.12 && MASS([pi0]b,[eta]a)<0.12) || (MASS([pi0]a,[eta]b)<0.12 && MASS([pi0]b,[eta]b)<0.12))")
FSCut.defineCut("delta","MASS([p+],[pi0])>1.4")
FSCut.defineCut("rf","OR(abs(RFDeltaT)<2.0)", "abs(RFDeltaT)>2.0",0.125)
FSCut.defineCut("eta","abs(MASS([eta])-0.548)<0.05","(abs(MASS([eta])-0.548+0.107)<0.025 || abs(MASS([eta])-0.548-0.103)<0.025)", 1.0)
FSCut.defineCut("pi0","abs(MASS([pi0])-0.135)<0.025","(abs(MASS([pi0])-0.135+0.0425)<0.0125 || abs(MASS([pi0])-0.135-0.0425)<0.0125)", 1.0)
FSHistogram.readHistogramCache();
FSModeCollection.addModeInfo("101_1").addCategory("m101_1");
READING CACHE:
root file = /d/grid17/ln16/PyAmpTools/demos/.cache.root
data file = /d/grid17/ln16/PyAmpTools/demos/.cache.dat
# prepare and cleanup previous plots
os.system("rm -rf plots")
os.system("mkdir plots")
0
Draw first set of plots#
c1 = ROOT.TCanvas("c1","c1",1200,800)
c1.Divide(3,2)
c1.cd(1)
hM1 = FSModeHistogram.getTH1F(FND,NT,"m101_1","EnUnusedSh","(100,0.0,1.0)","CUT(unusedTracks,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
hM1.SetTitle("MC: E_{unused} for -t in (0.1,0.3)")
hM1.SetXTitle("E_{unused} [GeV/c^{2}]")
hM1.SetYTitle("Events")
hM1.Draw()
cutUnusedE = ROOT.TLine(0.1,0,0.1,hM1.GetMaximum())
cutUnusedE.SetLineColor(ROOT.kRed)
cutUnusedE.Draw("same")
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","EnUnusedSh","(100,0.0,1.0)","CUT(unusedTracks,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
c1.cd(2)
hM6 = FSModeHistogram.getTH1F(FND,NT,"m101_1","ProdVz","(100,0.,100.0)","CUT(unusedTracks,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
hM6.SetTitle("MC: ProdVz for -t in (0.1,0.3)")
hM6.SetXTitle("ProdVz [GeV/c^{2}]")
hM6.SetYTitle("Events")
hM6.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","ProdVz","(100,0.,100.0)","CUT(unusedTracks,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
cutVz_low = ROOT.TLine(52,0,52,hM6.GetMaximum())
cutVz_low.SetLineColor(ROOT.kRed)
cutVz_low.Draw("same")
cutVz_hi = ROOT.TLine(78,0,78,hM6.GetMaximum())
cutVz_hi.SetLineColor(ROOT.kRed)
cutVz_hi.Draw("same")
c1.cd(3)
hM4a = FSModeHistogram.getTH1F(FND,NT,"m101_1","abs(-1*MASS2(GLUEXTARGET,-[p+]))","(100,0,1)","CUT(unusedTracks,zProton,chi2,unusedE,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
hM4a.SetTitle("MC: |-t|")
hM4a.SetXTitle("|-t| [GeV^{2}]")
hM4a.SetYTitle("Entries")
cutt_low = ROOT.TLine(0.1,0,0.1,hM4a.GetMaximum())
cutt_low.SetLineColor(ROOT.kRed)
hM4a.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","abs(-1*MASS2(GLUEXTARGET,-[p+]))","(100,0,1)","CUT(unusedTracks,zProton,chi2,unusedE,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
cutt_low.Draw("same")
cutt_hi = ROOT.TLine(0.3,0,0.3,hM4a.GetMaximum())
cutt_hi.SetLineColor(ROOT.kRed)
cutt_hi.Draw("same")
c1.cd(4)
hM4b = FSModeHistogram.getTH1F(FND,NT,"m101_1","EnPB","(125,5,12)","CUT(unusedTracks,zProton,chi2,unusedE,cet0103,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
hM4b.SetTitle("MC: E_{beam} for -t in (0.1,0.3)")
hM4b.SetXTitle("E_{beam} [GeV]")
hM4b.SetYTitle("Entries")
cutEb_low = ROOT.TLine(8.2,0,8.2,hM4b.GetMaximum())
cutEb_low.SetLineColor(ROOT.kRed)
hM4b.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","EnPB","(125,5,12)","CUT(unusedTracks,zProton,chi2,unusedE,cet0103,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
cutEb_low.Draw("same")
cutEb_hi = ROOT.TLine(8.8,0,8.8,hM4b.GetMaximum())
cutEb_hi.SetLineColor(ROOT.kRed)
cutEb_hi.Draw("same")
c1.cd(5)
hM4 = FSModeHistogram.getTH1F(FND,NT,"m101_1","Chi2DOF","(80,0,20)","CUT(unusedTracks,zProton,unusedE,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
hM4.SetTitle("MC: #chi^{2}/dof for -t in (0.1,0.3)")
hM4.SetXTitle("#chi^{2}/dof")
hM4.SetYTitle("Events")
cutChi2 = ROOT.TLine(3.3,0,3.3,hM4.GetMaximum())
cutChi2.SetLineColor(ROOT.kRed)
hM4.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","Chi2DOF","(80,0,20)","CUT(unusedTracks,zProton,unusedE,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
cutChi2.Draw("same")
c1.cd(6)
hM8 = FSModeHistogram.getTH1F(FND,NT,"m101_1","acos(COSINE([eta]a))*180/3.141","(240,0.,60)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,delta,rf,eta,pi0,rejectOmega,protMom)")
hM8.SetTitle("MC: #vartheta_{#gamma} for -t in (0.1,0.3)")
hM8.SetXTitle("#vartheta_{#gamma} [#circ]")
hM8.SetYTitle("Entries")
hM8.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","acos(COSINE([eta]a))*180/3.141","(240,0.,60)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,delta,rf,eta,pi0,rejectOmega,protMom)")
cutFidu1 = ROOT.TLine(2.5,0,2.5,hM8.GetMaximum())
cutFidu1.SetLineColor(ROOT.kRed)
cutFidu1.Draw("same")
cutFidu2 = ROOT.TLine(10.3,0,10.3,hM8.GetMaximum())
cutFidu2.SetLineColor(ROOT.kRed)
cutFidu2.Draw("same")
cutFidu3 = ROOT.TLine(11.9,0,11.9,hM8.GetMaximum())
cutFidu3.SetLineColor(ROOT.kRed)
cutFidu3.Draw("same")
c11 = ROOT.TCanvas("c11","c11",800,800)
c11.Divide(2,2)
c11.cd(1)
hM2 = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([eta])","(60,0.2,0.8)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,rejectOmega,protMom)")
hM2.SetTitle("MC: M(#gamma#gamma) for -t in (0.1,0.3)")
hM2.SetXTitle("M(#gamma#gamma) [GeV/c^{2}]")
hM2.SetYTitle("Events / 10 MeV/c^{2}")
hM2.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([eta])","(60,0.2,0.8)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,rejectOmega,protMom)")
cutEtaSigL = ROOT.TLine(0.548-0.05,0,0.548-0.05,hM2.GetMaximum())
cutEtaSigR = ROOT.TLine(0.548+0.05,0,0.548+0.05,hM2.GetMaximum())
cutEtaSBLowL = ROOT.TLine(0.42,0,0.42,hM2.GetMaximum())
cutEtaSBLowR = ROOT.TLine(0.47,0,0.47,hM2.GetMaximum())
cutEtaSBHiL = ROOT.TLine(0.63,0,0.63,hM2.GetMaximum())
cutEtaSBHiR = ROOT.TLine(0.68,0,0.68,hM2.GetMaximum())
cutEtaSigL.SetLineColor(ROOT.kRed)
cutEtaSigR.SetLineColor(ROOT.kRed)
cutEtaSBLowL.SetLineColor(ROOT.kRed)
cutEtaSBLowR.SetLineColor(ROOT.kRed)
cutEtaSBHiL.SetLineColor(ROOT.kRed)
cutEtaSBHiR.SetLineColor(ROOT.kRed)
cutEtaSigL.Draw("same")
cutEtaSigR.Draw("same")
cutEtaSBLowL.Draw("same")
cutEtaSBLowR.Draw("same")
cutEtaSBHiL.Draw("same")
cutEtaSBHiR.Draw("same")
c11.cd(2)
hM9 = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([pi0])","(60,0.0,0.3)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,rejectOmega,protMom)")
hM9.SetTitle("MC: M(#gamma#gamma) for -t in (0.1,0.3)")
hM9.SetXTitle("M(#gamma#gamma) [GeV/c^{2}]")
hM9.SetYTitle("Events / 10 MeV/c^{2}")
hM9.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([pi0])","(60,0.0,0.3)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,rejectOmega,protMom)")
hM9.SetTitle("MC: M(#gamma#gamma) for -t in (0.1,0.3)")
cutPi0SigL = ROOT.TLine(0.135-0.025,0,0.135-0.025,hM9.GetMaximum())
cutPi0SigR = ROOT.TLine(0.135+0.025,0,0.135+0.025,hM9.GetMaximum())
cutPi0SBLowL = ROOT.TLine(0.08,0,0.08,hM9.GetMaximum())
cutPi0SBLowR = ROOT.TLine(0.105,0,0.105,hM9.GetMaximum())
cutPi0SBHiL = ROOT.TLine(0.165,0,0.165,hM9.GetMaximum())
cutPi0SBHiR = ROOT.TLine(0.19,0,0.19,hM9.GetMaximum())
cutPi0SigL.SetLineColor(ROOT.kRed)
cutPi0SigR.SetLineColor(ROOT.kRed)
cutPi0SBLowL.SetLineColor(ROOT.kRed)
cutPi0SBLowR.SetLineColor(ROOT.kRed)
cutPi0SBHiL.SetLineColor(ROOT.kRed)
cutPi0SBHiR.SetLineColor(ROOT.kRed)
cutPi0SigL.Draw("same")
cutPi0SigR.Draw("same")
cutPi0SBLowL.Draw("same")
cutPi0SBLowR.Draw("same")
cutPi0SBHiL.Draw("same")
cutPi0SBHiR.Draw("same")
c11.cd(3)
hM3 = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([p+],[pi0])","(100,0.8,1.8)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,rf,eta,pi0,rejectOmega,protMom)")
hM3.SetTitle("MC: M(p#pi^{0}) for -t in (0.1,0.3)")
hM3.SetXTitle("M(#p#pi^{0}) [GeV/c^{2}]")
hM3.SetYTitle("Events / 10 MeV/c^{2}")
cutDelta = ROOT.TLine(1.4,0,1.4,hM3.GetMaximum())
cutDelta.SetLineColor(ROOT.kRed)
hM3.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([p+],[pi0])","(100,0.8,1.8)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,rf,eta,pi0,rejectOmega,protMom)")
cutDelta.Draw("same")
c11.cd(4)
hM5 = FSModeHistogram.getTH1F(FND,NT,"m101_1","RFDeltaT","(400,-20,20)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,eta,pi0,rejectOmega,protMom)")
hM5.SetTitle("MC: #Delta t_{RF} for -t in (0.1,0.3)")
hM5.SetXTitle("#Delta t_{RF}")
hM5.SetYTitle("Events")
cutRFSigL = ROOT.TLine(-2.0,0,-2.0,hM5.GetMaximum())
cutRFSigR = ROOT.TLine(2.0,0,2.0,hM5.GetMaximum())
cutRFSigL.SetLineColor(ROOT.kRed)
cutRFSigR.SetLineColor(ROOT.kRed)
hM5.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","RFDeltaT","(400,-20,20)","CUT(unusedTracks,zProton,unusedE,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,eta,pi0,rejectOmega,protMom)")
cutRFSigL.Draw("same")
cutRFSigR.Draw("same")
c1.Print("plots/p001_etapi.pdf")
CREATING HISTOGRAM... FSRootHist:000001 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000002 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000003 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000004 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000005 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000006 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000007 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000008 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000009 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000010 101_1 !!NO_FILE!! (entries = 0)
Info in <TCanvas::Print>: pdf file plots/p001_etapi.pdf has been created
Draw second set of plots#
c2 = ROOT.TCanvas("c2","c2",1200,400)
c2.Divide(3,1)
c2.cd(1)
hMetapi = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([eta],[pi0])","(100,0.5,2.5)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rf,eta,pi0,rejectOmega,protMom)")
hMetapiSig = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([eta],[pi0])","(100,0.5,2.5)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rejectOmega,protMom)*CUTWT(rf,eta,pi0)")
hMetapiBg = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([eta],[pi0])","(100,0.5,2.5)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rejectOmega,protMom)*CUTSBWT(rf,eta,pi0)*(-1.0)")
hMetapiSig.SetFillColor(ROOT.kBlue)
hMetapiBg.SetFillColor(ROOT.kRed)
hMetapi.SetTitle("MC: M(#eta#pi^{0}) for -t in (0.1,0.3)")
hMetapi.SetXTitle("M(#eta#pi^{0}) [GeV/c^{2}]")
hMetapi.SetYTitle("Events")
hMetapi.Draw()
hMetapiSig.Draw("hist,same")
hMetapiBg.Draw("hist,same")
c2.cd(2)
hMetapiSig.SetTitle("MC: M(#eta#pi^{0}) for -t in (0.1,0.3)")
hMetapiSig.SetXTitle("M(#eta#pi^{0}) [GeV/c^{2}]")
hMetapiSig.SetYTitle("Events")
hMetapiSig.DrawClone()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([eta],[pi0])","(100,0.5,2.5)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom)*CUTWT(rf,eta,pi0)")
print(hMetapiSig.Integral())
FSFitUtilities.createFunction(FSFitPOLY("p",1.04,1.56,1,0.0),1600.0,-900)
FSFitUtilities.createFunction(FSFitGAUS("g",1.04,1.56),500.0,1.32,0.05)
FSFitUtilities.createFunction("pg","p+g")
FSFitUtilities.fixParameters("p")
FSFitUtilities.migrad(hMetapiSig,"pg")
fpg = FSFitUtilities.getTF1("pg")
fpg.SetLineColor(ROOT.kRed)
fpg.SetLineStyle(ROOT.kSolid)
fpg.Draw("same")
fg = FSFitUtilities.getTF1("g")
fg.SetLineColor(ROOT.kBlue)
fg.SetLineStyle(ROOT.kSolid)
fg.Draw("same")
print("fg Integral:", fg.Integral(1.04,1.56))
c2.cd(3)
hMetapiVsGJCosTheta = FSModeHistogram.getTH2F(FND,NT,"m101_1","GJCOSTHETA([eta];[pi0];GLUEXBEAM):MASS([eta],[pi0])","(100,0.7,2.7,50,-1.,1.)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,rejectOmega,protMom)*CUTWT(rf,eta,pi0)")
hMetapiVsGJCosTheta.SetXTitle("M(#eta#pi^{0}) [GeV/c^{2}]")
hMetapiVsGJCosTheta.SetYTitle("cos#theta_{GJ}")
hMetapiVsGJCosTheta.Draw("colz")
c2.Print("plots/p002_etapi.pdf")
0.0
fg Integral: 499.999597944485
CREATING HISTOGRAM... FSRootHist:000011 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000012 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000013 101_1 !!NO_FILE!! (entries = 0)
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
1 p:a 1.60000e+03 1.60000e+02 no limits
**********
** 1 **FIX 1
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
2 p:b -9.00000e+02 9.00000e+01 no limits
**********
** 2 **FIX 2
**********
PARAMETER DEFINITIONS:
NO. NAME VALUE STEP SIZE LIMITS
3 g:N 5.00000e+02 5.00000e+01 no limits
4 g:M 1.32000e+00 1.32000e-01 no limits
5 g:W 5.00000e-02 5.00000e-03 no limits
**********
** 3 **SET STRATEGY 1
**********
NOW USING STRATEGY 1: TRY TO BALANCE SPEED AGAINST RELIABILITY
**********
** 4 **MIGRAD 1e+04
**********
FIRST CALL TO USER FUNCTION AT NEW START POINT, WITH IFLAG=4.
START MIGRAD MINIMIZATION. STRATEGY 1. CONVERGENCE WHEN EDM .LT. 1.00e-04
FCN=0 FROM MIGRAD STATUS=INITIATE 97 CALLS 98 TOTAL
EDM= unknown STRATEGY= 1 NO ERROR MATRIX
EXT PARAMETER CURRENT GUESS STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 p:a 1.60000e+03 fixed
2 p:b -9.00000e+02 fixed
3 g:N 5.00000e+02 5.00000e+01 0.00000e+00 0.00000e+00
4 g:M 1.32000e+00 1.32000e-01 0.00000e+00 0.00000e+00
5 g:W 5.00000e-02 5.00000e-03 0.00000e+00 0.00000e+00
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
MINUIT WARNING IN HESSE
============== Second derivative zero for parameter3
MNHESS FAILS AND WILL RETURN DIAGONAL MATRIX.
FCN=0 FROM HESSE STATUS=FAILED 11 CALLS 109 TOTAL
EDM=0 STRATEGY= 1 ERROR MATRIX UNCERTAINTY 100.0 per cent
EXT PARAMETER APPROXIMATE STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 p:a 1.60000e+03 fixed
2 p:b -9.00000e+02 fixed
3 g:N 5.00000e+02 1.41421e+00 0.00000e+00 0.00000e+00
4 g:M 1.32000e+00 1.41421e+00 0.00000e+00 0.00000e+00
5 g:W 5.00000e-02 1.41421e+00 0.00000e+00 0.00000e+00
MIGRAD MINIMIZATION HAS CONVERGED.
MIGRAD WILL VERIFY CONVERGENCE AND ERROR MATRIX.
MINUIT WARNING IN HESSE
============== Second derivative zero for parameter3
MNHESS FAILS AND WILL RETURN DIAGONAL MATRIX.
FCN=0 FROM MIGRAD STATUS=CONVERGED 119 CALLS 120 TOTAL
EDM=0 STRATEGY= 1 ERROR MATRIX UNCERTAINTY 100.0 per cent
EXT PARAMETER APPROXIMATE STEP FIRST
NO. NAME VALUE ERROR SIZE DERIVATIVE
1 p:a 1.60000e+03 fixed
2 p:b -9.00000e+02 fixed
3 g:N 5.00000e+02 1.41421e+00 0.00000e+00 0.00000e+00
4 g:M 1.32000e+00 1.41421e+00 0.00000e+00 0.00000e+00
5 g:W 5.00000e-02 1.41421e+00 0.00000e+00 0.00000e+00
EXTERNAL ERROR MATRIX. NDIM=1001 NPAR= 3 ERR DEF=1
2.000e+00 0.000e+00 0.000e+00
0.000e+00 2.000e+00 0.000e+00
0.000e+00 0.000e+00 2.000e+00
ERR MATRIX APPROXIMATE
PARAMETER CORRELATION COEFFICIENTS
NO. GLOBAL 3 4 5
3 0.00000 1.000 0.000 0.000
4 0.00000 0.000 1.000 0.000
5 0.00000 0.000 0.000 1.000
ERR MATRIX APPROXIMATE
CREATING HISTOGRAM... FSRootHist:000014 101_1 !!NO_FILE!! (entries = 0)
Info in <TCanvas::Print>: pdf file plots/p002_etapi.pdf has been created
Draw third set of plots#
c3 = ROOT.TCanvas("c3","c3",800,800)
c3.Divide(2,2)
c3.cd(1)
hMetapiSig.DrawClone()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([eta],[pi0])","(100,0.5,2.5)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
c3.cd(2)
hMpiproton = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([p+],[pi0])","(100,0.9,3.9)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
hMpiproton.SetTitle("MC: M(p#pi^{0}) for -t in (0.1,0.3)")
hMpiproton.SetXTitle("M(p#pi^{0}) [GeV/c^{2}]")
hMpiproton.SetYTitle("Events")
hMpiproton.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([p+],[pi0])","(100,0.9,3.9)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
c3.cd(3)
hMetaproton = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([p+],[eta])","(100,1.4,3.9)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
hMetaproton.SetTitle("MC: M(p#eta) for -t in (0.1,0.3)")
hMetaproton.SetXTitle("M(p#eta) [GeV/c^{2}]")
hMetaproton.SetYTitle("Events")
hMetaproton.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([p+],[eta])","(100,1.4,3.9)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
c3.cd(4)
hMpi0g = FSModeHistogram.getTH1F(FND,NT,"m101_1","MASS([pi0],[eta]b)","(100,0.,2.)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
hMpi0g.SetTitle("MC: M(#pi^{0}#gamma) for -t in (0.1,0.3)")
hMpi0g.SetXTitle("M(#pi^{0}#gamma) [GeV/c^{2}]")
hMpi0g.SetYTitle("Events")
hMpi0g.Draw()
if bggen: FSModeHistogram.drawMCComponentsSame(FND,NT,"m101_1","MASS([pi0],[eta]b)","(100,0.,2.)","CUT(unusedTracks,unusedE,zProton,chi2,cet0103,e8288,photFiducialA,photFiducialB,photFiducialC,photFiducialD,delta,protMom,rejectOmega)*CUTWT(rf,eta,pi0)")
FSHistogram.dumpHistogramCache()
c3.Print("plots/p003_etapi.pdf")
CREATING HISTOGRAM... FSRootHist:000015 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000016 101_1 !!NO_FILE!! (entries = 0)
CREATING HISTOGRAM... FSRootHist:000017 101_1 !!NO_FILE!! (entries = 0)
Info in <TCanvas::Print>: pdf file plots/p003_etapi.pdf has been created
# cleanup again
os.system("rm -rf plots")
0