]>
Commit | Line | Data |
---|---|---|
a9e2aefa | 1 | TH1F *mass1 = new TH1F("mass1","Invariant Mass",120,0,12); |
2 | TH1F *mass2 = new TH1F("mass2","Invariant Mass",100,9.0,10.5); | |
3 | TH1F *mass3 = new TH1F("mass3","Invariant Mass",100,2.5,3.5); | |
4 | void MUONstraggling (Int_t evNumber1=0,Int_t evNumber2=0) | |
5 | { | |
6 | ///////////////////////////////////////////////////////////////////////// | |
7 | // This macro is a small example of a ROOT macro | |
8 | // illustrating how to read the output of GALICE | |
9 | // and do some analysis. | |
10 | // | |
11 | ///////////////////////////////////////////////////////////////////////// | |
12 | ||
13 | // Dynamically link some shared libs | |
14 | static Float_t xmuon, ymuon; | |
15 | ||
16 | if (gClassTable->GetID("AliRun") < 0) { | |
17 | gROOT->LoadMacro("loadlibs.C"); | |
18 | loadlibs(); | |
19 | } | |
20 | ||
21 | // Connect the Root Galice file containing Geometry, Kine and Hits | |
22 | ||
23 | TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject("galice.root"); | |
24 | ||
25 | ||
26 | if (!file) { | |
27 | printf("\n Creating galice.root \n"); | |
28 | file = new TFile("galice.root"); | |
29 | } else { | |
30 | printf("\n galice.root found in file list"); | |
31 | } | |
32 | file->ls(); | |
33 | ||
34 | // Get AliRun object from file or create it if not on file | |
35 | if (!gAlice) { | |
36 | gAlice = (AliRun*)(file->Get("gAlice")); | |
37 | if (gAlice) printf("AliRun object found on file\n"); | |
38 | if (!gAlice) { | |
39 | printf("\n Create new gAlice object"); | |
40 | gAlice = new AliRun("gAlice","Alice test program"); | |
41 | } | |
42 | } | |
43 | ||
44 | ||
45 | // Create some histograms | |
46 | ||
47 | ||
48 | AliMUONChamber* iChamber; | |
49 | // | |
50 | // Loop over events | |
51 | // | |
52 | Int_t Nh=0; | |
53 | Int_t Nh1=0; | |
54 | for (Int_t nev=0; nev<= evNumber2; nev++) { | |
55 | Int_t nparticles = gAlice->GetEvent(nev); | |
56 | //cout << "nparticles " << nparticles <<endl; | |
57 | if (nev < evNumber1) continue; | |
58 | if (nparticles <= 0) return; | |
59 | ||
60 | AliMUON *MUON = (AliMUON*) gAlice->GetModule("MUON"); | |
61 | ||
62 | ||
63 | TTree *TH = gAlice->TreeH(); | |
64 | Int_t ntracks = TH->GetEntries(); | |
65 | // | |
66 | // Loop over tracks | |
67 | // | |
68 | ||
69 | Float_t pups[4]={0,0,0,0}; | |
70 | for (Int_t track=0; track<ntracks;track++) { | |
71 | Int_t Nelt=0; | |
72 | gAlice->ResetHits(); | |
73 | Int_t nbytes += TH->GetEvent(track); | |
74 | ||
75 | ||
76 | if (MUON) { | |
77 | for(AliMUONHit* mHit=(AliMUONHit*)MUON->FirstHit(-1); | |
78 | mHit; | |
79 | mHit=(AliMUONHit*)MUON->NextHit()) | |
80 | { | |
81 | Int_t nch = mHit->fChamber; // chamber number | |
82 | Float_t x = mHit->fX; // x-pos of hit | |
83 | Float_t y = mHit->fY; // y-pos | |
84 | Float_t z = mHit->fZ; // y-pos | |
85 | Float_t p=mHit->fPTot; | |
86 | Float_t px=mHit->fCxHit; | |
87 | Float_t py=mHit->fCyHit; | |
88 | Float_t pz=mHit->fCzHit; | |
89 | ||
90 | if (nch != 1) continue; | |
91 | ||
92 | Int_t ipart = mHit->fParticle; | |
93 | TClonesArray *fPartArray = gAlice->Particles(); | |
94 | TParticle *Part; | |
95 | Int_t ftrack = mHit->fTrack; | |
96 | Part = (TParticle*) fPartArray->UncheckedAt(ftrack); | |
97 | Int_t ipart = Part->GetPdgCode(); | |
98 | TParticle *Mother; | |
99 | Float_t px0=Part->Px(); | |
100 | Float_t py0=Part->Py(); | |
101 | Float_t pz0=Part->Pz(); | |
102 | Float_t e0=Part->Energy(); | |
103 | ||
104 | if (ipart == kMuonPlus || ipart == kMuonMinus) { | |
105 | // | |
106 | // Branson Correction | |
107 | // | |
108 | Float_t zch=505.; | |
109 | Float_t r=TMath::Sqrt(x*x+y*y); | |
110 | Float_t zb; | |
111 | ||
112 | if (r<26.3611) { | |
113 | zb=466.; | |
114 | } else { | |
115 | zb=441.; | |
116 | } | |
117 | ||
118 | Float_t xb=x-(zch-zb)*px/pz; | |
119 | Float_t yb=y-(zch-zb)*py/pz; | |
120 | Float_t tx=xb/zb; | |
121 | Float_t ty=yb/zb; | |
122 | pz=zb*p/TMath::Sqrt(zb*zb+xb*xb+yb*yb); | |
123 | px=pz*tx; | |
124 | py=pz*ty; | |
125 | ||
126 | // | |
127 | // Energy Correction | |
128 | // | |
129 | // | |
130 | Float_t corr=(p+CorrectP(p,mHit->fTheta))/p; | |
131 | pups[0]+=p*corr; | |
132 | pups[1]+=px*corr; | |
133 | pups[2]+=py*corr; | |
134 | pups[3]+=pz*corr; | |
135 | } | |
136 | } // hits | |
137 | } // if MUON | |
138 | } // tracks | |
139 | Float_t mass=TMath::Sqrt(pups[0]*pups[0]-pups[1]*pups[1]-pups[2]*pups[2]-pups[3]*pups[3]); | |
140 | mass1->Fill(mass, 1.); | |
141 | mass2->Fill(mass, 1.); | |
142 | mass3->Fill(mass, 1.); | |
143 | } // event | |
144 | //Create a canvas, set the view range, show histograms | |
145 | TCanvas *c1 = new TCanvas("c1","Vertices from electrons and positrons",400,10,600,700); | |
146 | c1->Divide(2,2); | |
147 | ||
148 | c1->cd(1); | |
149 | mass1->SetFillColor(42); | |
150 | mass1->SetXTitle("Mass (GeV)"); | |
151 | mass1->Draw(); | |
152 | ||
153 | c1->cd(2); | |
154 | mass2->SetFillColor(42); | |
155 | mass2->SetXTitle("Mass (GeV)"); | |
156 | mass2->Draw(); | |
157 | ||
158 | c1->cd(3); | |
159 | mass3->SetFillColor(42); | |
160 | mass3->SetXTitle("Mass (GeV)"); | |
161 | mass3->Draw(); | |
162 | ||
163 | } | |
164 | ||
165 | ||
166 | ||
167 | Float_t CorrectP(Float_t p, Float_t theta) | |
168 | { | |
169 | if (theta<3.) { | |
170 | //W | |
171 | if (p<15) { | |
172 | return 2.737+0.0494*p-0.001123*p*p; | |
173 | } else { | |
174 | return 3.0643+0.01346*p; | |
175 | } | |
176 | } else { | |
177 | //Pb | |
178 | if (p<15) { | |
179 | return 2.1380+0.0351*p-0.000853*p*p; | |
180 | } else { | |
181 | return 2.407+0.00702*p; | |
182 | } | |
183 | } | |
184 | } | |
185 | ||
186 | ||
187 | ||
188 | ||
189 | ||
190 |