]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGLF/FORWARD/analysis2/AliFMDMCEventInspector.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Segregated the Landau+Gaus function from the AliForwardUtil dumping ground. Other...
[u/mrichter/AliRoot.git] / PWGLF / FORWARD / analysis2 / AliFMDMCEventInspector.cxx
CommitLineData
e1f47419 1//
2// This class inspects the event
3//
4// Input:
5// - AliESDFMD object possibly corrected for sharing
6//
7// Output:
8// - A histogram of v_z of events with triggers.
9// - A histogram of v_z of events with vertex and triggers
10// - A histogram of trigger counters
11//
12// Note, that these are added to the master output list
13//
14// Corrections used:
15// - None
16//
17#include "AliFMDMCEventInspector.h"
18#include "AliLog.h"
19#include "AliAODForwardMult.h"
20#include "AliForwardUtil.h"
21#include "AliCentrality.h"
22#include "AliESDEvent.h"
23#include "AliMCEvent.h"
1f480471 24#include "AliStack.h"
e1f47419 25#include "AliGenPythiaEventHeader.h"
26#include "AliGenDPMjetEventHeader.h"
27#include "AliGenHijingEventHeader.h"
28// #include "AliGenHydjetEventHeader.h"
29#include "AliGenEposEventHeader.h"
30#include "AliGenHerwigEventHeader.h"
31#include "AliGenGeVSimEventHeader.h"
2d2513f0 32#include "AliAnalysisManager.h"
33#include "AliMCEventHandler.h"
e1f47419 34#include "AliHeader.h"
35#include <TList.h>
e308a636 36#include <TH2F.h>
1f480471 37#include <TParticle.h>
38#include <TMath.h>
a6682cdc 39#include <TParameter.h>
e308a636 40
e1f47419 41//====================================================================
42AliFMDMCEventInspector::AliFMDMCEventInspector()
43 : AliFMDEventInspector(),
44 fHVertex(0),
45 fHPhiR(0),
e308a636 46 fHB(0),
73b32206 47 fHMcC(0),
e308a636 48 fHBvsPart(0),
49 fHBvsBin(0),
50 fHBvsCent(0),
51 fHVzComp(0),
52 fHCentVsPart(0),
f7cfc454 53 fHCentVsBin(0),
73b32206 54 fHCentVsMcC(0),
f7cfc454 55 fProduction("")
e1f47419 56{
57 //
58 // Constructor
59 //
60}
61
62//____________________________________________________________________
63AliFMDMCEventInspector::AliFMDMCEventInspector(const char* /* name */)
64 : AliFMDEventInspector("fmdEventInspector"),
65 fHVertex(0),
66 fHPhiR(0),
e308a636 67 fHB(0),
73b32206 68 fHMcC(0),
e308a636 69 fHBvsPart(0),
70 fHBvsBin(0),
71 fHBvsCent(0),
72 fHVzComp(0),
73 fHCentVsPart(0),
f7cfc454 74 fHCentVsBin(0),
73b32206 75 fHCentVsMcC(0),
f7cfc454 76 fProduction("")
e1f47419 77{
78 //
79 // Constructor
80 //
81 // Parameters:
82 // name Name of object
83 //
e65b8b56 84 fMC = true;
e1f47419 85}
86
87//____________________________________________________________________
88AliFMDMCEventInspector::AliFMDMCEventInspector(const AliFMDMCEventInspector& o)
89 : AliFMDEventInspector(o),
90 fHVertex(0),
91 fHPhiR(0),
e308a636 92 fHB(0),
73b32206 93 fHMcC(0),
e308a636 94 fHBvsPart(0),
95 fHBvsBin(0),
96 fHBvsCent(0),
97 fHVzComp(0),
98 fHCentVsPart(0),
f7cfc454 99 fHCentVsBin(0),
73b32206 100 fHCentVsMcC(0),
f7cfc454 101 fProduction("")
e1f47419 102{
103 //
104 // Copy constructor
105 //
106 // Parameters:
107 // o Object to copy from
108 //
109}
110
111//____________________________________________________________________
112AliFMDMCEventInspector::~AliFMDMCEventInspector()
113{
114 //
115 // Destructor
116 //
117}
118//____________________________________________________________________
119AliFMDMCEventInspector&
120AliFMDMCEventInspector::operator=(const AliFMDMCEventInspector& o)
121{
122 //
123 // Assignement operator
124 //
125 // Parameters:
126 // o Object to assign from
127 //
128 // Return:
129 // Reference to this object
130 //
131 AliFMDEventInspector::operator=(o);
132 return *this;
133}
134
135//____________________________________________________________________
136void
5934a3e3 137AliFMDMCEventInspector::SetupForData(const TAxis& vtxAxis)
e1f47419 138{
139 //
140 // Initialize the object
141 //
142 // Parameters:
143 // vtxAxis Vertex axis in use
144 //
bfab35d9 145
146 // We temporary disable the displaced vertices so we can initialize
147 // the routine ourselves.
148 Bool_t saveDisplaced = fUseDisplacedVertices;
149 fUseDisplacedVertices = false;
5934a3e3 150 AliFMDEventInspector::SetupForData(vtxAxis);
bfab35d9 151 fUseDisplacedVertices = saveDisplaced;
e1f47419 152
e308a636 153 Int_t maxPart = 450;
154 Int_t maxBin = 225;
155 Int_t maxB = 25;
156 Int_t nVtx = vtxAxis.GetNbins();
157 Double_t lVtx = vtxAxis.GetXmin();
158 Double_t hVtx = vtxAxis.GetXmax();
159 fHVertex = new TH1F("vertex", "True vertex distribution", nVtx, lVtx, hVtx);
e1f47419 160 fHVertex->SetXTitle("v_{z} [cm]");
161 fHVertex->SetYTitle("# of events");
162 fHVertex->SetFillColor(kGreen+1);
163 fHVertex->SetFillStyle(3001);
164 fHVertex->SetDirectory(0);
165 // fHVertex->Sumw2();
166 fList->Add(fHVertex);
167
168 fHPhiR = new TH1F("phiR", "Event plane", 120, 0, 2*TMath::Pi());
169 fHPhiR->SetXTitle("#Phi_{R} [radians]");
170 fHPhiR->SetYTitle("# of events");
171 fHPhiR->SetFillColor(kGreen+1);
172 fHPhiR->SetFillStyle(3001);
173 fHPhiR->SetDirectory(0);
174 fList->Add(fHPhiR);
175
e308a636 176 fHB = new TH1F("b", "Impact parameter", 5*maxB, 0, maxB);
e1f47419 177 fHB->SetXTitle("b [fm]");
178 fHB->SetYTitle("# of events");
179 fHB->SetFillColor(kGreen+1);
180 fHB->SetFillStyle(3001);
181 fHB->SetDirectory(0);
182 fList->Add(fHB);
e308a636 183
73b32206 184 fHMcC = static_cast<TH1F*>(fHCent->Clone("mcC"));
185 fHMcC->SetFillColor(kCyan+2);
186 fHMcC->SetDirectory(0);
187 fList->Add(fHMcC);
188
e308a636 189 fHBvsPart = new TH2F("bVsParticipants", "Impact parameter vs Participants",
190 5*maxB, 0, maxB, maxPart, -.5, maxPart-.5);
191 fHBvsPart->SetXTitle("b [fm]");
192 fHBvsPart->SetYTitle("# of participants");
193 fHBvsPart->SetZTitle("Events");
194 fHBvsPart->SetDirectory(0);
195 fList->Add(fHBvsPart);
196
197 fHBvsBin = new TH2F("bVsBinary", "Impact parameter vs Binary Collisions",
198 5*maxB, 0, maxB, maxBin, -.5, maxBin-.5);
199 fHBvsBin->SetXTitle("b [fm]");
200 fHBvsBin->SetYTitle("# of binary collisions");
201 fHBvsBin->SetZTitle("Events");
202 fHBvsBin->SetDirectory(0);
203 fList->Add(fHBvsBin);
e1f47419 204
e308a636 205 fHBvsCent = new TH2F("bVsCentrality", "Impact parameter vs Centrality",
206 5*maxB, 0, maxB, fCentAxis->GetNbins(),
207 fCentAxis->GetXbins()->GetArray());
208 fHBvsCent->SetXTitle("b [fm]");
209 fHBvsCent->SetYTitle("Centrality [%]");
210 fHBvsCent->SetZTitle("Event");
211 fHBvsCent->SetDirectory(0);
212 fList->Add(fHBvsCent);
213
214
215 fHVzComp = new TH2F("vzComparison", "v_{z} truth vs reconstructed",
216 10*nVtx, lVtx, hVtx, 10*nVtx, lVtx, hVtx);
217 fHVzComp->SetXTitle("True v_{z} [cm]");
218 fHVzComp->SetYTitle("Reconstructed v_{z} [cm]");
219 fHVzComp->SetZTitle("Events");
220 fHVzComp->SetDirectory(0);
221 fList->Add(fHVzComp);
222
223 fHCentVsPart = new TH2F("centralityVsParticipans",
224 "# of participants vs Centrality",
225 maxPart, -.5, maxPart-.5, fCentAxis->GetNbins(),
226 fCentAxis->GetXbins()->GetArray());
227 fHCentVsPart->SetXTitle("Participants");
228 fHCentVsPart->SetYTitle("Centrality [%]");
229 fHCentVsPart->SetZTitle("Event");
230 fHCentVsPart->SetDirectory(0);
231 fList->Add(fHCentVsPart);
232
233 fHCentVsBin = new TH2F("centralityVsBinary",
234 "# of binary collisions vs Centrality",
235 maxBin, -.5, maxBin-.5, fCentAxis->GetNbins(),
236 fCentAxis->GetXbins()->GetArray());
237 fHCentVsBin->SetXTitle("Binary collisions");
238 fHCentVsBin->SetYTitle("Centrality [%]");
239 fHCentVsBin->SetZTitle("Event");
240 fHCentVsBin->SetDirectory(0);
241 fList->Add(fHCentVsBin);
73b32206 242
243 Int_t nC = fHCent->GetNbinsX();
244 Double_t cL = fHCent->GetXaxis()->GetXmin();
245 Double_t cH = fHCent->GetXaxis()->GetXmax();
246 fHCentVsMcC = new TH2F("centralityRecoVsMC",
247 "Centrality from reconstruction vs MC derived",
248 nC, cL, cH, nC, cL, cH);
249 fHCentVsMcC->SetDirectory(0);
250 fHCentVsMcC->SetStats(0);
251 fHCentVsMcC->SetXTitle("Centralty from Reco [%]");
252 fHCentVsMcC->SetYTitle("Centralty derived from Impact Par. [%]");
253 fHCentVsMcC->SetZTitle("Events");
254 fList->Add(fHCentVsMcC);
bfab35d9 255
256 if (fUseDisplacedVertices) fDisplacedVertex.SetupForData(fList, "", true);
e1f47419 257}
258
f7cfc454 259namespace
260{
261 TString& AppendField(TString& s, bool test, const char* f)
262 {
263 if (!test) return s;
264 if (!s.IsNull()) s.Append(",");
265 s.Append(f);
266 return s;
267 }
268}
269
270//____________________________________________________________________
271void
272AliFMDMCEventInspector::ReadProductionDetails(AliMCEvent* event)
273{
274 //Assign MC only triggers : True NSD etc.
275 AliHeader* header = event->Header();
276 AliGenEventHeader* genHeader = header->GenEventHeader();
277 AliCollisionGeometry* colGeometry =
278 dynamic_cast<AliCollisionGeometry*>(genHeader);
279 AliGenPythiaEventHeader* pythiaHeader =
280 dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
281 AliGenDPMjetEventHeader* dpmHeader =
282 dynamic_cast<AliGenDPMjetEventHeader*>(genHeader);
283 AliGenGeVSimEventHeader* gevHeader =
284 dynamic_cast<AliGenGeVSimEventHeader*>(genHeader);
285 AliGenHerwigEventHeader* herwigHeader =
286 dynamic_cast<AliGenHerwigEventHeader*>(genHeader);
287 AliGenHijingEventHeader* hijingHeader =
288 dynamic_cast<AliGenHijingEventHeader*>(genHeader);
289 // AliGenHydjetEventHeader* hydjetHeader =
290 // dynamic_cast<AliGenHydjetEventHeader*>(genHeader);
291 AliGenEposEventHeader* eposHeader =
292 dynamic_cast<AliGenEposEventHeader*>(genHeader);
293
294 AppendField(fProduction, colGeometry, "Geometry");
295 AppendField(fProduction, pythiaHeader, "Pythia");
296 AppendField(fProduction, dpmHeader, "DPM");
297 AppendField(fProduction, gevHeader, "GevSim");
298 AppendField(fProduction, herwigHeader, "Herwig");
299 AppendField(fProduction, hijingHeader, "Hijing");
300 // AppendField(fProduction, hydjetHeader, "Hydjet");
301 AppendField(fProduction, eposHeader, "EPOS");
302}
303
e1f47419 304//____________________________________________________________________
305UInt_t
306AliFMDMCEventInspector::ProcessMC(AliMCEvent* event,
307 UInt_t& triggers,
308 UShort_t& ivz,
309 Double_t& vz,
310 Double_t& b,
241cca4d 311 Double_t& c,
e308a636 312 Int_t& npart,
313 Int_t& nbin,
e1f47419 314 Double_t& phiR)
315{
316 //
317 // Process the event
318 //
319 // Parameters:
320 // event Input event
321 // triggers On return, the triggers fired
322 // ivz On return, the found vertex bin (1-based). A zero
323 // means outside of the defined vertex range
324 // vz On return, the z position of the interaction
325 // b On return, the impact parameter - < 0 if not available
241cca4d 326 // c On return, centrality estimate - < 0 if not available
e1f47419 327 // phiR On return, the event plane angle - < 0 if not available
328 //
329 // Return:
330 // 0 (or kOk) on success, otherwise a bit mask of error codes
331 //
332
333 // Check that we have an event
334 if (!event) {
335 AliWarning("No MC event found for input event");
336 return kNoEvent;
337 }
338
339 //Assign MC only triggers : True NSD etc.
340 AliHeader* header = event->Header();
341 AliGenEventHeader* genHeader = header->GenEventHeader();
e308a636 342 AliCollisionGeometry* colGeometry =
343 dynamic_cast<AliCollisionGeometry*>(genHeader);
e1f47419 344 AliGenPythiaEventHeader* pythiaHeader =
345 dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
346 AliGenDPMjetEventHeader* dpmHeader =
347 dynamic_cast<AliGenDPMjetEventHeader*>(genHeader);
348 AliGenGeVSimEventHeader* gevHeader =
349 dynamic_cast<AliGenGeVSimEventHeader*>(genHeader);
e1f47419 350 AliGenHerwigEventHeader* herwigHeader =
351 dynamic_cast<AliGenHerwigEventHeader*>(genHeader);
e308a636 352 // AliGenHijingEventHeader* hijingHeader =
353 // dynamic_cast<AliGenHijingEventHeader*>(genHeader);
e1f47419 354 // AliGenHydjetEventHeader* hydjetHeader =
355 // dynamic_cast<AliGenHydjetEventHeader*>(genHeader);
e308a636 356 // AliGenEposEventHeader* eposHeader =
357 // dynamic_cast<AliGenEposEventHeader*>(genHeader);
e1f47419 358
359 // Check if this is a single diffractive event
e308a636 360 Bool_t sd = kFALSE;
361 Double_t phi = -1111;
362 npart = 0;
363 nbin = 0;
364 b = -1;
241cca4d 365 c = -1;
ceb6a94f 366 phiR = -1111;
e308a636 367 if (colGeometry) {
368 b = colGeometry->ImpactParameter();
369 phi = colGeometry->ReactionPlaneAngle();
370 npart = (colGeometry->ProjectileParticipants() +
371 colGeometry->TargetParticipants());
372 nbin = colGeometry->NN();
373 }
d4d486f8 374 if (fDebug && !colGeometry) {
375 AliWarningF("Collision header of class %s is not a CollisionHeader",
376 genHeader->ClassName());
377 }
378
e1f47419 379 if(pythiaHeader) {
380 Int_t pythiaType = pythiaHeader->ProcessType();
1f480471 381 // 92 and 93 are SD
382 // 94 is DD
e1f47419 383 if (pythiaType==92 || pythiaType==93) sd = kTRUE;
e308a636 384 b = pythiaHeader->GetImpactParameter();
385 npart = 2; // Always 2 protons
386 nbin = 1; // Always 1 binary collision
e1f47419 387 }
241cca4d 388 if (b >= 0) {
389 if (b < 3.5) c = 2.5; //0-5%
390 else if (b < 4.95) c = 7.5; //5-10%
391 else if (b < 6.98) c = 15; //10-20%
392 else if (b < 8.55) c = 25; //20-30%
393 else if (b < 9.88) c = 35; //30-40%
394 else if (b < 11.04) c = 45; //40-50%
395 else c = 55; //50-60%
396 }
e308a636 397 if(dpmHeader) { // Also an AliCollisionGeometry
e1f47419 398 Int_t processType = dpmHeader->ProcessType();
1f480471 399 // 1 & 4 are ND
400 // 5 & 6 are SD
401 // 7 is DD
e1f47419 402 if (processType == 5 || processType == 6) sd = kTRUE;
e1f47419 403 }
404 if (gevHeader) {
405 phi = gevHeader->GetEventPlane();
406 }
e1f47419 407 if (herwigHeader) {
408 Int_t processType = herwigHeader->ProcessType();
409 // This is a guess
410 if (processType == 5 || processType == 6) sd = kTRUE;
e308a636 411 npart = 2; // Always 2 protons
412 nbin = 1; // Always 1 binary collision
e1f47419 413 }
e308a636 414 // if (hijingHeader) {
415 // b = hijingHeader->ImpactParameter();
416 // phi = hijingHeader->ReactionPlaneAngle();
417 // }
e1f47419 418 // if (hydjetHeader) {
419 // b = hydjetHeader->ImpactParameter();
420 // phi = hydjetHeader->ReactionPlaneAngle();
421 // }
e308a636 422 // if (eposHeader) {
423 // b = eposHeader->ImpactParameter();
424 // phi = eposHeader->ReactionPlaneAngle();
425 // }
e1f47419 426
427 // Normalize event plane angle to [0,2pi]
428 if (phi <= -1111) phiR = -1;
429 else {
430 while (true) {
431 if (phi < 0) phi += 2*TMath::Pi();
432 else if (phi > 2*TMath::Pi()) phi -= 2*TMath::Pi();
433 else break;
434 }
ceb6a94f 435 phiR = phi;
e1f47419 436 }
437
1f480471 438 // Do a check on particles
439 sd = IsSingleDiffractive(event->Stack());
440
e1f47419 441 // Set NSD flag
442 if(!sd) {
443 triggers |= AliAODForwardMult::kMCNSD;
444 fHTriggers->Fill(kMCNSD+0.5);
445 }
446
447 // Get the primary vertex from EG
448 TArrayF vtx;
449 genHeader->PrimaryVertex(vtx);
450 vz = vtx[2];
451
4274f0b5 452 DMSG(fDebug, 2, "vz=%f, phiR=%f, b=%f, npart=%d, nbin=%d",
453 vz, phiR, b, npart, nbin);
454
e1f47419 455 fHVertex->Fill(vz);
456 fHPhiR->Fill(phiR);
457 fHB->Fill(b);
73b32206 458 fHMcC->Fill(c);
e308a636 459 fHBvsPart->Fill(b, npart);
460 fHBvsBin->Fill(b, nbin);
e1f47419 461
241cca4d 462 if(fUseDisplacedVertices) {
bfab35d9 463#if 0
241cca4d 464 // Put the vertex at fixed locations
465 Double_t zvtx = vz;
466 Double_t ratio = zvtx/37.5;
467 if(ratio > 0) ratio = ratio + 0.5;
468 if(ratio < 0) ratio = ratio - 0.5;
469 Int_t ratioInt = Int_t(ratio);
470 zvtx = 37.5*((Double_t)ratioInt);
471 if(TMath::Abs(zvtx) > 999)
472 return kBadVertex;
bfab35d9 473#endif
474 if (!fDisplacedVertex.ProcessMC(event))
475 return kBadVertex;
476 if (fDisplacedVertex.IsSatellite())
477 vz = fDisplacedVertex.GetVertexZ();
241cca4d 478 }
479
e1f47419 480 // Check for the vertex bin
d2286e26 481 ivz = fVtxAxis.FindBin(vz);
e1f47419 482 if (ivz <= 0 || ivz > fHEventsTr->GetXaxis()->GetNbins()) {
483 if (fDebug > 3) {
484 AliWarning(Form("Vertex @ %f outside of range [%f,%f]",
d2286e26 485 vz, fVtxAxis.GetXmin(), fVtxAxis.GetXmax())); }
e1f47419 486 ivz = 0;
487 return kBadVertex;
488 }
489
490
491 return kOk;
492}
73b32206 493//____________________________________________________________________
494Bool_t
495AliFMDMCEventInspector::ReadCentrality(const AliESDEvent& esd,
496 Double_t& cent,
497 UShort_t& qual) const
498{
499 //
500 // Read centrality from event
501 //
502 // Parameters:
503 // esd Event
504 // cent On return, the centrality or negative if not found
505 //
506 // Return:
507 // False on error, true otherwise
508 //
509 Bool_t ret = AliFMDEventInspector::ReadCentrality(esd, cent, qual);
510 if (qual != 0) {
511 AliCentrality* centObj = const_cast<AliESDEvent&>(esd).GetCentrality();
512 if (!centObj) return ret;
513
514 // For MC, we allow `bad' centrality selections
515 cent = centObj->GetCentralityPercentileUnchecked(fCentMethod);
516 }
517 return ret;
518}
1f480471 519
520//____________________________________________________________________
521namespace {
522 Double_t rapidity(TParticle* p, Double_t mass)
523 {
524 Double_t pt = p->Pt();
525 Double_t pz = p->Pz();
526 Double_t ee = TMath::Sqrt(pt*pt+pz*pz+mass*mass);
527 if (TMath::Abs(ee - TMath::Abs(pz)) < 1e-9) return TMath::Sign(1e30, pz);
528 return .5 * TMath::Log((ee + pz) / (ee - pz));
529 }
530}
531
532//____________________________________________________________________
533Bool_t
534AliFMDMCEventInspector::IsSingleDiffractive(AliStack* stack,
535 Double_t xiMin,
536 Double_t xiMax) const
537{
538 // Re-implementation of AliPWG0Helper::IsHadronLevelSingleDiffrative
539 //
540 // This is re-implemented here to be indendent of the PWG0 library.
541 TParticle* p1 = 0; // Particle with least y
542 TParticle* p2 = 0; // Particle with largest y
543 Double_t y1 = 1e10; // y of p1
544 Double_t y2 = -1e10; // y of p2
545
546 // Loop over primaries
547 for (Int_t i = 0; i < stack->GetNprimary(); i++) {
548 TParticle* p = stack->Particle(i);
549 if (!p) continue;
550
551 Int_t pdg = TMath::Abs(p->GetPdgCode());
552 Int_t c1 = p->GetFirstDaughter();
553 Int_t s = p->GetStatusCode();
554 Int_t mfl = Int_t(pdg/TMath::Power(10,Int_t(TMath::Log10(pdg))));
555
556 // Select final state charm and beauty
557 if (c1 > -1 || s != 1) mfl = 0;
558
559 // Check if this is a primary, pi0, Sigma0, ???, or most
560 // significant digit is larger than or equal to 4
561 if (!(stack->IsPhysicalPrimary(i) ||
562 pdg == 111 ||
563 pdg == 3212 ||
564 pdg == 3124 ||
565 mfl >= 4)) continue;
566
567 Int_t m1 = p->GetFirstMother();
568 if (m1 > 0) {
569 TParticle* pm1 = stack->Particle(m1);
570 Int_t mpdg = TMath::Abs(pm1->GetPdgCode());
571 // Check if mother is a p0, Simga0, or ???
572 if (mpdg == 111 || mpdg == 3124 || mpdg == 3212) continue;
573 }
574
575 // Calculate the rapidity of the particle
576 Double_t mm = (pdg != 3124 ? p->GetMass() : 1.5195);
577 Double_t yy = rapidity(p, mm);
578
579 // Check if the rapidity of this particle is further out than any
580 // of the preceding particles
581 if (yy < y1) {
582 y1 = yy;
583 p1 = p;
584 }
585 if (yy > y2) {
586 y2 = yy;
587 p2 = p;
588 }
589 }
590 if (!p1 || !p2) return false;
591
592 // Calculate rapidities assuming protons
593 y1 = TMath::Abs(rapidity(p1, 0.938));
594 y2 = TMath::Abs(rapidity(p2, 0.938));
595
596 // Check if both or just one is a proton
597 Int_t pdg1 = p1->GetPdgCode();
598 Int_t pdg2 = p2->GetPdgCode();
599 Int_t arm = -99999;
600 if (pdg1 == 2212 && pdg2 == 2212)
601 arm = (y1 > y2 ? 0 : 1);
602 else if (pdg1 == 2212)
603 arm = 0;
604 else if (pdg2 == 2212)
605 arm = 1;
606 else
607 return false;
608
5bb5d1f6 609 // Rapidity shift
d4d486f8 610 Double_t m02s = (fEnergy > 0 ? 1 - 2 * p1->Energy() / fEnergy : 0);
611 Double_t m12s = (fEnergy > 0 ? 1 - 2 * p2->Energy() / fEnergy : 0);
1f480471 612
613 if (arm == 0 && m02s > xiMin && m02s < xiMax) return true;
614 if (arm == 1 && m12s > xiMin && m12s < xiMax) return true;
615
616 return false;
617}
e1f47419 618
e308a636 619//____________________________________________________________________
620Bool_t
621AliFMDMCEventInspector::CompareResults(Double_t vz, Double_t trueVz,
73b32206 622 Double_t cent, Double_t mcC,
623 Double_t b,
e308a636 624 Int_t npart, Int_t nbin)
625{
626 fHVzComp->Fill(trueVz, vz);
627 fHBvsCent->Fill(b, cent);
628 fHCentVsPart->Fill(npart, cent);
629 fHCentVsBin->Fill(nbin, cent);
73b32206 630 fHCentVsMcC->Fill(cent, mcC);
e308a636 631
632 return true;
633}
73b32206 634
635
e1f47419 636//
637// EOF
638//
639