Useful macros to test the absorber correction for muon reconstruction
[u/mrichter/AliRoot.git] / MUON / MUONmassPlot.C
CommitLineData
c2a319d4 1// macro to make invariant mass plots
2// for combinations of 2 muons with opposite charges,
3// from root file "MUONtrackReco.root" containing the result of track reconstruction,
4// generated by the macro "MUONrecoNtuple.C".
5// Histograms are stored on the "MUONmassPlot.root" file.
6// A model for macros using the Ntuple in the file "MUONtrackReco.root"
7// may be found in "MUONtrackRecoModel.C":
8// it has been obtained by reading with Root a file "MUONtrackReco.root",
9// and executing the command:
10// MUONtrackReco->MakeCode("MUONtrackRecoModel.C")
11
12// Arguments:
13// FirstEvent (default 0)
14// LastEvent (default 0)
15// ResType (default 553)
16// 553 for Upsilon, anything else for J/Psi
17// NSigma (default 3)
18// the number of combinations is counted around the resonance mass
19// within +/- NSigma times the nominal sigma's
20// (0.099 GeV for Upsilon, 0.0615 GeV for J/Psi)
21// Chi2Cut (default 100)
22// to keep only tracks with chi2 per d.o.f. < Chi2Cut
23// PtCut (default 1)
24// to keep only tracks with transverse momentum > PtCut
25
eff57104 26// IMPORTANT NOTICE FOR USERS:
27// under "root" or "root.exe", execute the following commands:
28// 1. "gSystem->SetIncludePath("-I$ALICE_ROOT/MUON -I$ALICE_ROOT/STEER")" to get the right path at compilation time
29// 2. ".x loadlibs.C" to load the shared libraries
30// 3. ".L MUONrecoNtuple.C+"
31// 4. ".x MUONmassPlot.C()" with the right arguments according to the list above
32
c2a319d4 33void MUONmassPlot(Int_t FirstEvent = 0, Int_t LastEvent = 0, Int_t ResType = 553, Float_t Nsig = 3., Float_t Chi2Cut = 100., Float_t PtCut = 1.)
34{
35 cout << "MUONmassPlot" << endl;
36 cout << "FirstEvent" << FirstEvent << endl;
37 cout << "LastEvent" << LastEvent << endl;
38 cout << "ResType" << ResType << endl;
39 cout << "Nsig" << Nsig << endl;
40 cout << "Chi2Cut" << Chi2Cut << endl;
41 cout << "PtCut" << PtCut << endl;
42
43 //////////////////////////////////////////////////////////
44 // This file has been automatically generated
45 // (Thu Sep 21 14:53:11 2000 by ROOT version2.25/02)
46 // from TTree MUONtrackReco/MUONtrackReco
47 // found on file: MUONtrackReco.root
48 //////////////////////////////////////////////////////////
49
50
51 //Reset ROOT and connect tree file
52 gROOT->Reset();
53 TFile *f = (TFile*)gROOT->GetListOfFiles()->FindObject("MUONtrackReco.root");
54 if (!f) {
55 f = new TFile("MUONtrackReco.root");
56 }
57 TTree *MUONtrackReco = (TTree*)gDirectory->Get("MUONtrackReco");
eff57104 58 MUONtrackReco->SetMakeClass(1);
c2a319d4 59
60//Declaration of leaves types
61 Int_t fEvent;
62 UInt_t fUniqueID;
63 UInt_t fBits;
64 Int_t Tracks_;
65 Int_t Tracks_fCharge[5];
66 Float_t Tracks_fPxRec[5];
67 Float_t Tracks_fPyRec[5];
68 Float_t Tracks_fPzRec[5];
69 Float_t Tracks_fZRec[5];
70 Float_t Tracks_fZRec1[5];
71 Int_t Tracks_fNHits[5];
72 Float_t Tracks_fChi2[5];
73 Float_t Tracks_fPxGen[5];
74 Float_t Tracks_fPyGen[5];
75 Float_t Tracks_fPzGen[5];
76 UInt_t Tracks_fUniqueID[5];
77 UInt_t Tracks_fBits[5];
78
79 //Set branch addresses
80 //MUONtrackReco->SetBranchAddress("Header",&Header);
81 MUONtrackReco->SetBranchAddress("fEvent",&fEvent);
82 MUONtrackReco->SetBranchAddress("fUniqueID",&fUniqueID);
83 MUONtrackReco->SetBranchAddress("fBits",&fBits);
eff57104 84// MUONtrackReco->SetBranchAddress("Tracks_",&Tracks_);
85 MUONtrackReco->SetBranchAddress("Tracks",&Tracks_);
c2a319d4 86 MUONtrackReco->SetBranchAddress("Tracks.fCharge",Tracks_fCharge);
87 MUONtrackReco->SetBranchAddress("Tracks.fPxRec",Tracks_fPxRec);
88 MUONtrackReco->SetBranchAddress("Tracks.fPyRec",Tracks_fPyRec);
89 MUONtrackReco->SetBranchAddress("Tracks.fPzRec",Tracks_fPzRec);
90 MUONtrackReco->SetBranchAddress("Tracks.fZRec",Tracks_fZRec);
91 MUONtrackReco->SetBranchAddress("Tracks.fZRec1",Tracks_fZRec1);
92 MUONtrackReco->SetBranchAddress("Tracks.fNHits",Tracks_fNHits);
93 MUONtrackReco->SetBranchAddress("Tracks.fChi2",Tracks_fChi2);
94 MUONtrackReco->SetBranchAddress("Tracks.fPxGen",Tracks_fPxGen);
95 MUONtrackReco->SetBranchAddress("Tracks.fPyGen",Tracks_fPyGen);
96 MUONtrackReco->SetBranchAddress("Tracks.fPzGen",Tracks_fPzGen);
97 MUONtrackReco->SetBranchAddress("Tracks.fUniqueID",Tracks_fUniqueID);
98 MUONtrackReco->SetBranchAddress("Tracks.fBits",Tracks_fBits);
99
100// This is the loop skeleton
101// To read only selected branches, Insert statements like:
102// MUONtrackReco->SetBranchStatus("*",0); // disable all branches
103// TTreePlayer->SetBranchStatus("branchname",1); // activate branchname
104
105 Int_t nentries = MUONtrackReco->GetEntries();
106
107 Int_t nbytes = 0;
108 // for (Int_t i=0; i<nentries;i++) {
109 // nbytes += MUONtrackReco->GetEntry(i);
110 // }
111
112 /////////////////////////////////////////////////////////////////
113 // Here comes the specialized part for MUONmassPlot
114 /////////////////////////////////////////////////////////////////
115
116 // File for histograms and histogram booking
117 TFile *histoFile = new TFile("MUONmassPlot.root", "RECREATE");
118 TH1F *hPtMuon = new TH1F("hPtMuon", "Muon Pt (GeV/c)", 100, 0., 20.);
119 TH1F *hChi2PerDof = new TH1F("hChi2PerDof", "Muon track chi2/d.o.f.", 100, 0., 20.);
120 TH1F *hInvMassAll = new TH1F("hInvMassAll", "Mu+Mu- invariant mass (GeV/c2)", 240, 0., 12.);
121 if (ResType = 553) TH1F *hInvMassRes = new TH1F("hInvMassRes", "Mu+Mu- invariant mass (GeV/c2) around Upsilon", 60, 8., 11.);
122 else TH1F *hInvMassRes = new TH1F("hInvMassRes", "Mu+Mu- invariant mass (GeV/c2) around J/Psi", 80, 1., 5.);
123
124 // Loop over events
125 for (Int_t event = FirstEvent; event <= TMath::Min(LastEvent, nentries - 1); event++) {
126 // get current event
127 nbytes += MUONtrackReco->GetEntry(event);
128 // loop over all reconstructed tracks (also first track of combination)
129 for (Int_t t1 = 0; t1 < Tracks_; t1++) {
130 // transverse momentum
131 Float_t pt1 = TMath::Sqrt(Tracks_fPxRec[t1] * Tracks_fPxRec[t1] + Tracks_fPyRec[t1] * Tracks_fPyRec[t1]);
132 // chi2 per d.o.f.
133 Float_t ch1 = Tracks_fChi2[t1] / (2.0 * Tracks_fNHits[t1] - 5);
134 // condition for good track (Chi2Cut and PtCut)
135 if ((ch1 < Chi2Cut) && (pt1 > PtCut)) {
136 // fill histos hPtMuon and hChi2PerDof
137 hPtMuon->Fill(pt1);
138 hChi2PerDof->Fill(ch1);
139 // loop over second track of combination
140 for (Int_t t2 = t1 + 1; t2 < Tracks_; t2++) {
141 // transverse momentum
142 Float_t pt2 = TMath::Sqrt(Tracks_fPxRec[t2] * Tracks_fPxRec[t2] + Tracks_fPyRec[t2] * Tracks_fPyRec[t2]);
143 // chi2 per d.o.f.
144 Float_t ch2 = Tracks_fChi2[t2] / (2.0 * Tracks_fNHits[t2] - 5);
145 // condition for good track (Chi2Cut and PtCut)
146 if ((ch2 < Chi2Cut) && (pt2 > PtCut)) {
147 // condition for opposite charges
148 if ((Tracks_fCharge[t1] * Tracks_fCharge[t2]) == -1) {
149 // invariant mass
150 Float_t invMass = MuPlusMuMinusMass(Tracks_fPxRec[t1], Tracks_fPyRec[t1], Tracks_fPzRec[t1], Tracks_fPxRec[t2], Tracks_fPyRec[t2], Tracks_fPzRec[t2]);
151 // fill histos hInvMassAll and hInvMassRes
152 hInvMassAll->Fill(invMass);
153 hInvMassRes->Fill(invMass);
154 } //if ((Tracks_fCharge[t1] * Tracks_fCharge[t2]) == -1)
155 } // if ((Tracks_fChi2[t2] < Chi2Cut) && (pt2 > PtCut))
156 } // for (Int_t t2 = t1 + 1; t2 < Tracks_; t2++)
157 } // if ((Tracks_fChi2[t1] < Chi2Cut) && (pt1 > PtCut))
158 } // for (Int_t t1 = 0; t1 < Tracks_; t1++)
159 } // for (Int_t event = FirstEvent;
160
161 histoFile->Write();
162 histoFile->Close();
163}
164
165Float_t MuPlusMuMinusMass(Float_t Px1, Float_t Py1, Float_t Pz1, Float_t Px2, Float_t Py2, Float_t Pz2)
166{
167 Float_t muonMass = 0.10566;
168 Float_t e1 = TMath::Sqrt(muonMass * muonMass + Px1 * Px1 + Py1 * Py1 + Pz1 * Pz1);
169 Float_t e2 = TMath::Sqrt(muonMass * muonMass + Px2 * Px2 + Py2 * Py2 + Pz2 * Pz2);
170 return (TMath::Sqrt(2.0 * (muonMass * muonMass + e1 * e2 - Px1 * Px2 - Py1 * Py2 - Pz1 * Pz2)));
171}