git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Fixed up sNN recognision for Heavy Ion
[u/mrichter/AliRoot.git] / PWG2 / FORWARD / analysis / AliFMDAnaParameters.cxx
CommitLineData
d7346eed 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15//
16//The design of this class is based on the AliFMDParameters class. Its purpose
17//is to hold parameters for the analysis such as background correction and
18//fit functions.
19//
20//Author: Hans Hjersing Dalsgaard, NBI, hans.dalsgaard@cern.ch
21//
22
ffa07380 23// #include "AliFMDDebug.h" // ALILOG_H
d7346eed 24#include "AliFMDAnaParameters.h" // ALIFMDPARAMETERS_H
78f6f750 25//#include <AliCDBManager.h> // ALICDBMANAGER_H
26//#include <AliCDBEntry.h> // ALICDBMANAGER_H
27//#include "AliFMDRing.h"
d7346eed 28#include <AliLog.h>
29#include <Riostream.h>
30#include <sstream>
31#include <TSystem.h>
32#include <TH2D.h>
8dc823cc 33#include <TF1.h>
f58a4769 34#include <TMath.h>
5a79fd59 35#include "AliTriggerAnalysis.h"
36#include "AliPhysicsSelection.h"
da0805e2 37#include "AliBackgroundSelection.h"
38#include "AliMultiplicity.h"
b64db9b1 39#include "AliESDEvent.h"
40#include "AliESDVertex.h"
3d4a1473 41#include "AliFMDAnaCalibBackgroundCorrection.h"
42#include "AliFMDAnaCalibEnergyDistribution.h"
43#include "AliFMDAnaCalibEventSelectionEfficiency.h"
44#include "AliFMDAnaCalibSharingEfficiency.h"
45#include "TFile.h"
507687cd 46
d7346eed 47//====================================================================
48ClassImp(AliFMDAnaParameters)
49#if 0
50 ; // This is here to keep Emacs for indenting the next line
51#endif
52
d05586f1 53//const char* AliFMDAnaParameters::fgkBackgroundCorrection = "FMD/Correction/Background";
54//const char* AliFMDAnaParameters::fgkEnergyDists = "FMD/Correction/EnergyDistribution";
7e2bf482 55const char* AliFMDAnaParameters::fgkBackgroundID = "background";
d05586f1 56const char* AliFMDAnaParameters::fgkEnergyDistributionID = "energydistributions";
b64db9b1 57const char* AliFMDAnaParameters::fgkEventSelectionEffID = "eventselectionefficiency";
7e2bf482 58const char* AliFMDAnaParameters::fgkSharingEffID = "sharingefficiency";
d7346eed 59//____________________________________________________________________
60AliFMDAnaParameters* AliFMDAnaParameters::fgInstance = 0;
61
62//____________________________________________________________________
63
64AliFMDAnaParameters*
65AliFMDAnaParameters::Instance()
66{
67 // Get static instance
68 if (!fgInstance) fgInstance = new AliFMDAnaParameters;
69 return fgInstance;
70}
71
72//____________________________________________________________________
73AliFMDAnaParameters::AliFMDAnaParameters() :
74 fIsInit(kFALSE),
b82e76e0 75 fBackground(0),
78f6f750 76 fEnergyDistribution(0),
b64db9b1 77 fEventSelectionEfficiency(0),
ae26bdd7 78 fSharingEfficiency(0),
78f6f750 79 fCorner1(4.2231, 26.6638),
80 fCorner2(1.8357, 27.9500),
e34665a4 81 fEnergyPath("$ALICE_ROOT/PWG2/FORWARD/corrections/EnergyDistribution"),
82 fBackgroundPath("$ALICE_ROOT/PWG2/FORWARD/corrections/Background"),
83 fEventSelectionEffPath("$ALICE_ROOT/PWG2/FORWARD/corrections/EventSelectionEfficiency"),
84 fSharingEffPath("$ALICE_ROOT/PWG2/FORWARD/corrections/SharingEfficiency"),
b64db9b1 85 fProcessPrimary(kFALSE),
86 fProcessHits(kFALSE),
0b0a4ae5 87 fTrigger(kMB1),
11e702b3 88 fEnergy(k900),
01622a51 89 fMagField(k5G),
02a4ce29 90 fSpecies(kPP),
5a79fd59 91 fPhysicsSelection(0),
92 fRealData(kFALSE),
93 fSPDlowLimit(0),
94 fSPDhighLimit(999999999),
70d74659 95 fCentralSelection(kFALSE),
ff293341 96 fSharingObjectPresent(kTRUE),
1a7dad8a 97 fNumberOfEtaBinsToCut(1),
98 fEtaLowBinLimits(),
059c7c6b 99 fEtaHighBinLimits(),
100 fTriggerInel(kFALSE),
101 fTriggerNSD(kFALSE),
541c19ed 102 fTriggerEmpty(kFALSE),
507687cd 103 fUseBuiltInNSD(kFALSE),
cbfdb0cc 104 fInelGtZero(kFALSE),
105 fRunDndeta(kTRUE),
106 fRunBFCorrelation(kTRUE),
107 fRunMultiplicity(kTRUE)
d7346eed 108{
aa303f0c 109 // Default constructor
507687cd 110 // fPhysicsSelection = new AliPhysicsSelection;
111 // fPhysicsSelection->SetAnalyzeMC(kTRUE); //For the background correction. This is reset in Init if relevant
4b8bdb60 112 // fPhysicsSelection->SetUseBXNumbers(kFALSE);
113
059c7c6b 114
507687cd 115 // AliBackgroundSelection* backgroundSelection = new AliBackgroundSelection("bg","bg");
116 // backgroundSelection->Init();
117 // fPhysicsSelection->AddBackgroundIdentification(backgroundSelection);
5a79fd59 118 //fPhysicsSelection->Initialize(104792);
85da855f 119 // Do not use this - it is only for IO
120 fgInstance = this;
aa303f0c 121
d7346eed 122}
123//____________________________________________________________________
879c0f42 124const char* AliFMDAnaParameters::GetPath(const char* species) const
125{
aa303f0c 126 //Get path of object
ca2331d7 127 static TString* path = new TString();
128
0b0a4ae5 129 if(species == fgkBackgroundID)
ca2331d7 130 path->Form("%s/%s_%d_%d_%d_%d_%d_%d.root",
131 fBackgroundPath.Data(),
132 fgkBackgroundID,
133 fEnergy,
134 fTrigger,
135 fMagField,
136 fSpecies,
137 fInelGtZero,
138 0);
0b0a4ae5 139 if(species == fgkEnergyDistributionID)
ca2331d7 140 path->Form("%s/%s_%d_%d_%d_%d_%d_%d.root",
141 fEnergyPath.Data(),
142 fgkEnergyDistributionID,
143 fEnergy,
144 fTrigger,
145 fMagField,
146 fSpecies,
147 fRealData,
148 0);
0b0a4ae5 149 if(species == fgkEventSelectionEffID)
ca2331d7 150 path->Form("%s/%s_%d_%d_%d_%d_%d_%d.root",
151 fEventSelectionEffPath.Data(),
152 fgkEventSelectionEffID,
153 fEnergy,
154 fTrigger,
155 fMagField,
156 fSpecies,
157 fInelGtZero,
158 0);
7e2bf482 159 if(species == fgkSharingEffID)
ca2331d7 160 path->Form("%s/%s_%d_%d_%d_%d_%d_%d.root",
161 fSharingEffPath.Data(),
162 fgkSharingEffID,
163 fEnergy,
164 fTrigger,
165 fMagField,
166 fSpecies,
167 0,
168 0);
169 return path->Data();
0b0a4ae5 170}
171//____________________________________________________________________
d7346eed 172void AliFMDAnaParameters::Init(Bool_t forceReInit, UInt_t what)
173{
aa303f0c 174 // Initialize the parameters manager. We need to get stuff from files here.
507687cd 175
176 /* AliPhysicsSelection* test = (AliPhysicsSelection*)((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetEventSelection();
177
178 if(fPhysicsSelection) {
179 if(!fRealData) {
180 fPhysicsSelection->SetAnalyzeMC(kTRUE);
181 }
182 else fPhysicsSelection->SetAnalyzeMC(kFALSE);
183 }
184 */
d7346eed 185 if (forceReInit) fIsInit = kFALSE;
186 if (fIsInit) return;
187 if (what & kBackgroundCorrection) InitBackground();
188 if (what & kEnergyDistributions) InitEnergyDists();
b64db9b1 189 if (what & kEventSelectionEfficiency) InitEventSelectionEff();
7e2bf482 190 if (what & kSharingEfficiency) InitSharingEff();
5a79fd59 191
d7346eed 192 fIsInit = kTRUE;
ff293341 193
194 if(fBackground)
195 FindEtaLimits();
196
507687cd 197
059c7c6b 198
d7346eed 199}
200//____________________________________________________________________
201
202void AliFMDAnaParameters::InitBackground() {
aa303f0c 203 //Init background correction objects.
0b0a4ae5 204
205 TFile* fin = TFile::Open(GetPath(fgkBackgroundID));
78f6f750 206
207 if (!fin) return;
d7346eed 208
d05586f1 209 fBackground = dynamic_cast<AliFMDAnaCalibBackgroundCorrection*>(fin->Get(fgkBackgroundID));
b82e76e0 210 if (!fBackground) AliFatal("Invalid background object from CDB");
d7346eed 211
ff293341 212
213
d7346eed 214}
b64db9b1 215
d7346eed 216//____________________________________________________________________
217
218void AliFMDAnaParameters::InitEnergyDists() {
aa303f0c 219 //Init energy distributions
ca2331d7 220
0b0a4ae5 221 TFile* fin = TFile::Open(GetPath(fgkEnergyDistributionID));
aa303f0c 222
78f6f750 223 if (!fin) return;
d7346eed 224
d05586f1 225 fEnergyDistribution = dynamic_cast<AliFMDAnaCalibEnergyDistribution*>(fin->Get(fgkEnergyDistributionID));
d7346eed 226
b82e76e0 227 if (!fEnergyDistribution) AliFatal("Invalid background object from CDB");
d7346eed 228
b64db9b1 229}
230
231//____________________________________________________________________
232
233void AliFMDAnaParameters::InitEventSelectionEff() {
aa303f0c 234 //Init event selection objects
235
0b0a4ae5 236 TFile* fin = TFile::Open(GetPath(fgkEventSelectionEffID));
237
b64db9b1 238 if (!fin) return;
239
240 fEventSelectionEfficiency = dynamic_cast<AliFMDAnaCalibEventSelectionEfficiency*>(fin->Get(fgkEventSelectionEffID));
241 if (!fEventSelectionEfficiency) AliFatal("Invalid background object from CDB");
242
aa303f0c 243}
244
245//____________________________________________________________________
246
247void AliFMDAnaParameters::InitSharingEff() {
248
70d74659 249 fSharingObjectPresent = kTRUE;
aa303f0c 250 TFile* fin = TFile::Open(GetPath(fgkSharingEffID));
251
70d74659 252 if (!fin) {
253 fSharingObjectPresent = kFALSE;
254 return;
255 }
aa303f0c 256
257 fSharingEfficiency = dynamic_cast<AliFMDAnaCalibSharingEfficiency*>(fin->Get(fgkSharingEffID));
70d74659 258 if (!fSharingEfficiency) {
259 fSharingObjectPresent = kFALSE;
260 return;
261 }
aa303f0c 262
ff293341 263}
264//____________________________________________________________________
265void AliFMDAnaParameters::FindEtaLimits() {
266
2f323ab5 267 fEtaLowBinLimits.SetBins(4,0,4,2,0,2,GetNvtxBins(),0,GetNvtxBins());
268 fEtaHighBinLimits.SetBins(4,0,4,2,0,2,GetNvtxBins(),0,GetNvtxBins());
269 for(Int_t det=0; det<=3;det++) {
270 Int_t nRings = (det<=1 ? 1 : 2);
ff293341 271 for (UShort_t ir = 0; ir < nRings; ir++) {
272 Char_t ringChar = (ir == 0 ? 'I' : 'O');
273 for(Int_t v =0; v<GetNvtxBins(); v++) {
274 fEtaLowBinLimits.SetBinContent(det,ir,v,GetFirstEtaBinFromMap(v, det, ringChar));
275 fEtaHighBinLimits.SetBinContent(det,ir,v,GetLastEtaBinFromMap(v, det, ringChar));
2f323ab5 276 //std::cout<<det<<" "<<ringChar<<" "<<fEtaLowBinLimits.GetBinContent(det,ir,v)<<" "<<fEtaHighBinLimits.GetBinContent(det,ir,v)<<std::endl;
ff293341 277 }
278 }
279 }
280
7e2bf482 281}
cbfdb0cc 282//____________________________________________________________________
879c0f42 283void AliFMDAnaParameters::SetEnergy(Float_t cmsNNGeV)
284{
285 if (TMath::Abs(cmsNNGeV - 900.) < 10) fEnergy = k900;
286 if (TMath::Abs(cmsNNGeV - 2400.) < 10) fEnergy = k2400;
287 if (TMath::Abs(cmsNNGeV - 2750.) < 10) fEnergy = k2750;
288 if (TMath::Abs(cmsNNGeV - 5500.) < 40) fEnergy = k5500;
289 if (TMath::Abs(cmsNNGeV - 7000.) < 10) fEnergy = k7000;
290 if (TMath::Abs(cmsNNGeV - 10000.) < 10) fEnergy = k10000;
291 if (TMath::Abs(cmsNNGeV - 14000.) < 10) fEnergy = k14000;
292}
293//____________________________________________________________________
294void AliFMDAnaParameters::SetMagField(Float_t bkG)
295{
296 if (TMath::Abs(bkG - 5.) < 1 ) fMagField = k5G;
297 if (TMath::Abs(bkG + 5.) < 1 ) fMagField = k5Gnegative;
298 if (TMath::Abs(bkG) < 1) fMagField = k0G;
299}
300//____________________________________________________________________
301void AliFMDAnaParameters::SetCollisionSystem(const TString& sys)
302{
6afac6e9 303 TString s(sys);
304 s.ToLower();
305 if (s.Contains("p-p") || s.Contains("pp")) fSpecies = kPP;
306 else if (s.Contains("pb-pb") || s.Contains("pbpb")) fSpecies = kPbPb;
307 else if (s.Contains("a-a") || s.Contains("aa")) fSpecies = kPbPb;
879c0f42 308}
309
310//____________________________________________________________________
311void AliFMDAnaParameters::SetParametersFromESD(AliESDEvent* esd)
312{
cbfdb0cc 313
6afac6e9 314
879c0f42 315 SetCollisionSystem(esd->GetBeamType());
6afac6e9 316
317 Float_t energy = esd->GetBeamEnergy();
318 // Correct to center of mass per nucleon (cmsNN) - LHC gives it as
319 // cmsNN * Z / 2
320 if (fSpecies == kPbPb) energy = energy / 208 * 82;
321 SetEnergy(2*energy);
322 SetMagField(esd->GetMagneticField());
323
324 // Float_t magfield = esd->GetCurrentL3();
325 // if (TMath::Abs(magfield - 30000.) < 10 ) fMagField = k5G;
326 // if (TMath::Abs(magfield + 30000.) < 10 ) fMagField = k5Gnegative;
327 // if (TMath::Abs(magfield) < 10 ) fMagField = k0G;
cbfdb0cc 328
cbfdb0cc 329
cbfdb0cc 330
331 Init(kTRUE);
332
333}
334
7e2bf482 335//____________________________________________________________________
336
879c0f42 337void AliFMDAnaParameters::PrintStatus(Bool_t showpaths) const
85da855f 338{
aa303f0c 339 //Print current status
f6b21230 340 TString energystring;
341 switch(fEnergy) {
342 case k900:
343 energystring.Form("900 GeV"); break;
879c0f42 344 case k2400:
345 energystring.Form("2400 GeV"); break;
346 case k2750:
347 energystring.Form("2750 GeV"); break;
348 case k5500:
349 energystring.Form("5500 GeV"); break;
f6b21230 350 case k7000:
351 energystring.Form("7000 GeV"); break;
352 case k10000:
353 energystring.Form("10000 GeV"); break;
354 case k14000:
355 energystring.Form("14000 GeV"); break;
356 default:
357 energystring.Form("invalid energy"); break;
358 }
359 TString triggerstring;
360 switch(fTrigger) {
361 case kMB1:
362 triggerstring.Form("Minimum bias 1"); break;
363 case kMB2:
364 triggerstring.Form("Minimum bias 2"); break;
365 case kSPDFASTOR:
366 triggerstring.Form("SPD FAST OR"); break;
367 case kNOCTP:
368 triggerstring.Form("NO TRIGGER TEST"); break;
369 default:
370 energystring.Form("invalid trigger"); break;
371 }
372 TString magstring;
373 switch(fMagField) {
374 case k5G:
6afac6e9 375 magstring.Form("+5 kGaus"); break;
f6b21230 376 case k0G:
377 magstring.Form("0 kGaus"); break;
f7356393 378 case k5Gnegative:
379 magstring.Form("-5 kGaus"); break;
f6b21230 380 default:
f7356393 381 magstring.Form("invalid mag field %d", fMagField); break;
f6b21230 382 }
383 TString collsystemstring;
384 switch(fSpecies) {
385 case kPP:
6afac6e9 386 collsystemstring.Form("p-p"); break;
f6b21230 387 case kPbPb:
6afac6e9 388 collsystemstring.Form("Pb-Pb"); break;
f6b21230 389 default:
390 collsystemstring.Form("invalid collision system"); break;
391 }
392
059c7c6b 393 TString datastring;
394 switch(fRealData) {
395 case kTRUE:
396 datastring.Form("Nature"); break;
397 case kFALSE:
398 datastring.Form("MC"); break;
399 default:
400 datastring.Form("Unknown"); break ;
401 }
f6b21230 402
507687cd 403 TString InelString;
404 if(fInelGtZero) InelString = "INEL > 0";
405 else InelString = "INEL";
406
407 std::cout<<"Energy = "<<energystring.Data()<<std::endl;
408 std::cout<<"Trigger = "<<triggerstring.Data()<<std::endl;
409 std::cout<<"Mag Field = "<<magstring.Data()<<std::endl;
410 std::cout<<"Coll System = "<<collsystemstring.Data()<<std::endl;
411 std::cout<<"Data origin = "<<datastring.Data()<<std::endl;
412 std::cout<<"Basic trigger: "<<InelString.Data()<<std::endl;
879c0f42 413
414 if (showpaths) {
415 TString bg = GetPath(fgkBackgroundID);
416 TString es = GetPath(fgkEventSelectionEffID);
417 TString ed = GetPath(fgkEnergyDistributionID);
418 TString me = GetPath(fgkSharingEffID);
419 std::cout << "2nd maps: " << bg << "\n"
420 << "Event sel.: " << es << "\n"
421 << "Energy dist.: " << ed << "\n"
422 << "Merge eff.: " << me << std::endl;
423 }
f6b21230 424
425}
426
d7346eed 427//____________________________________________________________________
428Float_t AliFMDAnaParameters::GetVtxCutZ() {
aa303f0c 429 //Get the z vtx cut in analysis
d7346eed 430 if(!fIsInit) {
431 AliWarning("Not initialized yet. Call Init() to remedy");
432 return -1;
433 }
434
b82e76e0 435 return fBackground->GetVtxCutZ();
d7346eed 436}
437
438//____________________________________________________________________
439Int_t AliFMDAnaParameters::GetNvtxBins() {
aa303f0c 440 //Get number of vtx bins
d7346eed 441 if(!fIsInit) {
442 AliWarning("Not initialized yet. Call Init() to remedy");
443 return -1;
444 }
445
b82e76e0 446 return fBackground->GetNvtxBins();
d7346eed 447}
448//____________________________________________________________________
5754671c 449TH1F* AliFMDAnaParameters::GetEnergyDistribution(Int_t det, Char_t ring, Float_t eta) {
da0805e2 450
5754671c 451 return fEnergyDistribution->GetEnergyDistribution(det, ring, eta);
8dc823cc 452}
453//____________________________________________________________________
5a79fd59 454TH1F* AliFMDAnaParameters::GetEmptyEnergyDistribution(Int_t det, Char_t ring) {
da0805e2 455
5a79fd59 456 return fEnergyDistribution->GetEmptyEnergyDistribution(det, ring);
457}
458//____________________________________________________________________
459TH1F* AliFMDAnaParameters::GetRingEnergyDistribution(Int_t det, Char_t ring) {
da0805e2 460
5a79fd59 461 return fEnergyDistribution->GetRingEnergyDistribution(det, ring);
462}
463//____________________________________________________________________
5754671c 464Float_t AliFMDAnaParameters::GetSigma(Int_t det, Char_t ring, Float_t eta) {
aa303f0c 465 //Get sigma of Landau fits
d7346eed 466 if(!fIsInit) {
467 AliWarning("Not initialized yet. Call Init() to remedy");
468 return 0;
469 }
470
5754671c 471 TH1F* hEnergyDist = GetEnergyDistribution(det,ring, eta);
472 TF1* fitFunc = hEnergyDist->GetFunction("FMDfitFunc");
473 if(!fitFunc) {
5a79fd59 474 //AliWarning(Form("No function for FMD%d%c, eta %f",det,ring,eta));
5754671c 475 return 1024;
476 }
477 Float_t sigma = fitFunc->GetParameter(2);
8dc823cc 478 return sigma;
479}
480
481
482//____________________________________________________________________
5754671c 483Float_t AliFMDAnaParameters::GetMPV(Int_t det, Char_t ring, Float_t eta) {
aa303f0c 484 //Get MPV of landau fits
8dc823cc 485 if(!fIsInit) {
486 AliWarning("Not initialized yet. Call Init() to remedy");
487 return 0;
488 }
5a79fd59 489 //AliFMDAnaParameters* pars = AliFMDAnaParameters::Instance();
5754671c 490 TH1F* hEnergyDist = GetEnergyDistribution(det,ring,eta);
491 TF1* fitFunc = hEnergyDist->GetFunction("FMDfitFunc");
daedf077 492
5754671c 493 if(!fitFunc) {
879c0f42 494 AliWarning(Form("No function for FMD%d%c, eta %f (%d)",det,ring,eta,
495 GetEtaBin(eta)));
5754671c 496 return 1024;
497 }
daedf077 498
d05586f1 499 Float_t mpv = fitFunc->GetParameter(1);
500 return mpv;
d7346eed 501}
502//____________________________________________________________________
5a79fd59 503Float_t AliFMDAnaParameters::GetConstant(Int_t det, Char_t ring, Float_t eta) {
aa303f0c 504 //Get constant parameter of Landau fits
5a79fd59 505 if(!fIsInit) {
506 AliWarning("Not initialized yet. Call Init() to remedy");
507 return 0;
508 }
509
510 TH1F* hEnergyDist = GetEnergyDistribution(det,ring,eta);
511 TF1* fitFunc = hEnergyDist->GetFunction("FMDfitFunc");
512 if(!fitFunc) {
513 AliWarning(Form("No function for FMD%d%c, eta %f",det,ring,eta));
514 return 0;
515 }
516
517 Float_t mpv = fitFunc->GetParameter(0);
518 return mpv;
519}
520//____________________________________________________________________
1b418b63 521Float_t
522AliFMDAnaParameters::Get2MIPWeight(Int_t det, Char_t ring, Float_t eta)
523{
aa303f0c 524 //Get 2 MIP weights of convoluted Landau fits
5754671c 525 if(!fIsInit) {
526 AliWarning("Not initialized yet. Call Init() to remedy");
527 return 0;
528 }
529
530 TH1F* hEnergyDist = GetEnergyDistribution(det,ring,eta);
531 TF1* fitFunc = hEnergyDist->GetFunction("FMDfitFunc");
532 if(!fitFunc) return 0;
d05586f1 533 Float_t twoMIPweight = fitFunc->GetParameter(3);
5754671c 534
535
536
d05586f1 537 if(twoMIPweight < 1e-05)
538 twoMIPweight = 0;
5754671c 539
d05586f1 540 return twoMIPweight;
5754671c 541}
542//____________________________________________________________________
1b418b63 543Float_t
544AliFMDAnaParameters::Get3MIPWeight(Int_t det, Char_t ring, Float_t eta)
545{
aa303f0c 546 //Get 3 MIP weights of convoluted Landau fits
5754671c 547 if(!fIsInit) {
548 AliWarning("Not initialized yet. Call Init() to remedy");
549 return 0;
550 }
551
552 TH1F* hEnergyDist = GetEnergyDistribution(det,ring,eta);
553 TF1* fitFunc = hEnergyDist->GetFunction("FMDfitFunc");
554 if(!fitFunc) return 0;
d05586f1 555 Float_t threeMIPweight = fitFunc->GetParameter(4);
5754671c 556
d05586f1 557 if(threeMIPweight < 1e-05)
558 threeMIPweight = 0;
5754671c 559
d05586f1 560 Float_t twoMIPweight = fitFunc->GetParameter(3);
5754671c 561
d05586f1 562 if(twoMIPweight < 1e-05)
563 threeMIPweight = 0;
5754671c 564
d05586f1 565 return threeMIPweight;
5754671c 566}
567//____________________________________________________________________
1b418b63 568Int_t AliFMDAnaParameters::GetNetaBins()
569{
570 return GetBackgroundCorrection(1,'I',5)->GetNbinsX();
5754671c 571}
572//____________________________________________________________________
1b418b63 573Float_t AliFMDAnaParameters::GetEtaMin()
574{
daedf077 575 return GetBackgroundCorrection(1,'I',5)->GetXaxis()->GetXmin();
5754671c 576}
577//____________________________________________________________________
1b418b63 578Float_t AliFMDAnaParameters::GetEtaMax()
579{
580 return GetBackgroundCorrection(1,'I',5)->GetXaxis()->GetXmax();
5a79fd59 581}
582//____________________________________________________________________
1b418b63 583Int_t AliFMDAnaParameters::GetEtaBin(Float_t eta)
584{
5a79fd59 585 TAxis testaxis(GetNetaBins(),GetEtaMin(),GetEtaMax());
586 Int_t binnumber = testaxis.FindBin(eta) ;
587
588 return binnumber;
5754671c 589}
590//____________________________________________________________________
8dc823cc 591
d7346eed 592TH2F* AliFMDAnaParameters::GetBackgroundCorrection(Int_t det,
593 Char_t ring,
594 Int_t vtxbin) {
aa303f0c 595 //Get background correction histogram
d7346eed 596 if(!fIsInit) {
597 AliWarning("Not initialized yet. Call Init() to remedy");
598 return 0;
599 }
b82e76e0 600
601
602
603 if(vtxbin > fBackground->GetNvtxBins()) {
d7346eed 604 AliWarning(Form("No background object for vertex bin %d", vtxbin));
605 return 0;
606 }
607
b82e76e0 608 return fBackground->GetBgCorrection(det,ring,vtxbin);
d7346eed 609}
4fb49e43 610//____________________________________________________________________
507687cd 611TH2F* AliFMDAnaParameters::GetBackgroundCorrectionNSD(Int_t det,
612 Char_t ring,
613 Int_t vtxbin) {
614 //Get background correction histogram for NSD event class
615 if(!fIsInit) {
616 AliWarning("Not initialized yet. Call Init() to remedy");
617 return 0;
618 }
619
620 if(vtxbin > fBackground->GetNvtxBins()) {
621 AliWarning(Form("No background object for vertex bin %d", vtxbin));
622 return 0;
623 }
624
625 if(fBackground->GetNSDBgCorrection(det,ring,vtxbin))
626 return fBackground->GetNSDBgCorrection(det,ring,vtxbin);
627 else
879c0f42 628 AliWarning("No NSD background map. You get usual one. "
629 "Difference is probably negligible");
507687cd 630
631 return fBackground->GetBgCorrection(det,ring,vtxbin);
632}
633//____________________________________________________________________
4fb49e43 634
635TH1F* AliFMDAnaParameters::GetDoubleHitCorrection(Int_t det,
636 Char_t ring) {
aa303f0c 637 //Get correction for several hits in strips for p+p
4fb49e43 638 if(!fIsInit) {
639 AliWarning("Not initialized yet. Call Init() to remedy");
640 return 0;
641 }
642
643 return fBackground->GetDoubleHitCorrection(det,ring);
644}
f58a4769 645//_____________________________________________________________________
04f1ff3d 646TH1F* AliFMDAnaParameters::GetSPDDeadCorrection(Int_t vtxbin) {
647
648 //Get correction for several hits in strips for p+p
649 if(!fIsInit) {
650 AliWarning("Not initialized yet. Call Init() to remedy");
651 return 0;
652 }
653
654 return fBackground->GetSPDDeadCorrection(vtxbin);
655}
656//_____________________________________________________________________
f7356393 657TH1F* AliFMDAnaParameters::GetFMDDeadCorrection(Int_t vtxbin) {
658
659 //Get correction for several hits in strips for p+p
660 if(!fIsInit) {
661 AliWarning("Not initialized yet. Call Init() to remedy");
662 return 0;
663 }
664
665 return fBackground->GetFMDDeadCorrection(vtxbin);
666}
667//_____________________________________________________________________
b64db9b1 668Float_t AliFMDAnaParameters::GetEventSelectionEfficiency(Int_t vtxbin) {
aa303f0c 669 //Get event selection efficiency object
b64db9b1 670 if(!fIsInit) {
671 AliWarning("Not initialized yet. Call Init() to remedy");
672 return 0;
673 }
674 return fEventSelectionEfficiency->GetCorrection(vtxbin);
675
daedf077 676}
677//_____________________________________________________________________
678Float_t AliFMDAnaParameters::GetVtxSelectionEffFromMC() {
679
680 if(!fIsInit) {
681 AliWarning("Not initialized yet. Call Init() to remedy");
682 return 0;
683 }
684 return fEventSelectionEfficiency->GetVtxToTriggerRatio();
685
686
ab3e0abc 687}
688//_____________________________________________________________________
df2a9c32 689TH2F* AliFMDAnaParameters::GetEventSelectionEfficiency(TString trig, Int_t vtxbin, Char_t ring) {
ab3e0abc 690 //Get event selection efficiency object
059c7c6b 691
df2a9c32 692 //TString test = trig;
693 if(!trig.Contains("NSD") && !trig.Contains("INEL")) {
059c7c6b 694 AliWarning("Event selection efficiency only available for INEL and NSD");
695 return 0;
696 }
ab3e0abc 697 if(!fIsInit) {
698 AliWarning("Not initialized yet. Call Init() to remedy");
699 return 0;
700 }
059c7c6b 701 return fEventSelectionEfficiency->GetCorrection(trig,vtxbin,ring);
ab3e0abc 702
7e2bf482 703}
704//_____________________________________________________________________
705TH1F* AliFMDAnaParameters::GetSharingEfficiency(Int_t det, Char_t ring, Int_t vtxbin) {
aa303f0c 706 //Get sharing efficiency object
7e2bf482 707 if(!fIsInit) {
708 AliWarning("Not initialized yet. Call Init() to remedy");
709 return 0;
710 }
711
712 return fSharingEfficiency->GetSharingEff(det,ring,vtxbin);
713
25f47050 714}
715//_____________________________________________________________________
716TH1F* AliFMDAnaParameters::GetSharingEfficiencyTrVtx(Int_t det, Char_t ring, Int_t vtxbin) {
aa303f0c 717 //Get sharing efficiency object TrVtx
25f47050 718 if(!fIsInit) {
719 AliWarning("Not initialized yet. Call Init() to remedy");
720 return 0;
721 }
722
723 return fSharingEfficiency->GetSharingEffTrVtx(det,ring,vtxbin);
724
b64db9b1 725}
726//_____________________________________________________________________
aa303f0c 727Float_t AliFMDAnaParameters::GetMaxR(Char_t ring) const {
728 //Get max R of ring
78f6f750 729 Float_t radius = 0;
730 if(ring == 'I')
731 radius = 17.2;
732 else if(ring == 'O')
733 radius = 28.0;
734 else
735 AliWarning("Unknown ring - must be I or O!");
736
737 return radius;
738}
739//_____________________________________________________________________
d05586f1 740Float_t AliFMDAnaParameters::GetMinR(Char_t ring) const{
aa303f0c 741 //Get min R of ring
78f6f750 742 Float_t radius = 0;
743 if(ring == 'I')
744 radius = 4.5213;
745 else if(ring == 'O')
746 radius = 15.4;
747 else
748 AliWarning("Unknown ring - must be I or O!");
749
750 return radius;
751
752}
753//_____________________________________________________________________
754void AliFMDAnaParameters::SetCorners(Char_t ring) {
aa303f0c 755 //Set corners (taken from nominal geometry)
78f6f750 756 if(ring == 'I') {
757 fCorner1.Set(4.9895, 15.3560);
758 fCorner2.Set(1.8007, 17.2000);
759 }
760 else {
761 fCorner1.Set(4.2231, 26.6638);
762 fCorner2.Set(1.8357, 27.9500);
763 }
764
765}
766//_____________________________________________________________________
d05586f1 767Float_t AliFMDAnaParameters::GetPhiFromSector(UShort_t det, Char_t ring, UShort_t sec) const
f58a4769 768{
aa303f0c 769 //Get phi from sector
f58a4769 770 Int_t nsec = (ring == 'I' ? 20 : 40);
771 Float_t basephi = 0;
772 if(det == 1)
773 basephi = 1.72787594;
774 if(det == 2 && ring == 'I')
775 basephi = 0.15707963;
776 if(det == 2 && ring == 'O')
777 basephi = 0.078539818;
778 if(det == 3 && ring == 'I')
779 basephi = 2.984513044;
780 if(det == 3 && ring == 'O')
781 basephi = 3.06305289;
782
783 Float_t step = 2*TMath::Pi() / nsec;
784 Float_t phi = 0;
785 if(det == 3)
786 phi = basephi - sec*step;
787 else
788 phi = basephi + sec*step;
789
790 if(phi < 0)
791 phi = phi +2*TMath::Pi();
792 if(phi > 2*TMath::Pi() )
793 phi = phi - 2*TMath::Pi();
794
795 return phi;
796}
797//_____________________________________________________________________
d05586f1 798Float_t AliFMDAnaParameters::GetEtaFromStrip(UShort_t det, Char_t ring, UShort_t sec, UShort_t strip, Float_t zvtx) const
f58a4769 799{
aa303f0c 800 //Calculate eta from strip with vertex (redundant with AliESDFMD::Eta)
78f6f750 801 Float_t rad = GetMaxR(ring)-GetMinR(ring);
d05586f1 802 Float_t nStrips = (ring == 'I' ? 512 : 256);
803 Float_t segment = rad / nStrips;
78f6f750 804 Float_t r = GetMinR(ring) + segment*strip;
f58a4769 805 Float_t z = 0;
806 Int_t hybrid = sec / 2;
807
808 if(det == 1) {
809 if(!(hybrid%2)) z = 320.266; else z = 319.766;
810 }
811 if(det == 2 && ring == 'I' ) {
812 if(!(hybrid%2)) z = 83.666; else z = 83.166;
813 }
814 if(det == 2 && ring == 'O' ) {
815 if(!(hybrid%2)) z = 74.966; else z = 75.466;
816 }
817 if(det == 3 && ring == 'I' ) {
818 if(!(hybrid%2)) z = -63.066; else z = -62.566;
819 }
820 if(det == 3 && ring == 'O' ) {
821 if(!(hybrid%2)) z = -74.966; else z = -75.466;
822 }
823
824 //std::cout<<det<<" "<<ring<<" "<<sec<<" "<<hybrid<<" "<<z<<std::endl;
825
826 // Float_t r = TMath::Sqrt(TMath::Power(x,2)+TMath::Power(y,2));
827 Float_t theta = TMath::ATan2(r,z-zvtx);
828 Float_t eta = -1*TMath::Log(TMath::Tan(0.5*theta));
829
830 return eta;
b64db9b1 831}
832
833//_____________________________________________________________________
834
a91555b9 835Bool_t AliFMDAnaParameters::GetVertex(const AliESDEvent* esd, Double_t* vertexXYZ)
d7346eed 836{
aa303f0c 837 //Get the vertex from the ESD
1b418b63 838 const AliESDVertex* vertex = esd->GetPrimaryVertexSPD();
5a79fd59 839
1b418b63 840 if (!vertex) return kFALSE;
841
842 vertex->GetXYZ(vertexXYZ);
059c7c6b 843
5a79fd59 844 //if(vertexXYZ[0] == 0 || vertexXYZ[1] == 0 )
845 // return kFALSE;
70d74659 846
5a79fd59 847 if(vertex->GetNContributors() <= 0)
848 return kFALSE;
849
850 if(vertex->GetZRes() > 0.1 )
851 return kFALSE;
b64db9b1 852
f55d559b 853 return vertex->GetStatus();
b64db9b1 854
855}
5a79fd59 856//____________________________________________________________________
059c7c6b 857void AliFMDAnaParameters::SetTriggerStatus(const AliESDEvent *esd) {
9f55be54 858
059c7c6b 859 //ULong64_t triggerMask = esd->GetTriggerMask();
9f55be54 860
507687cd 861 AliPhysicsSelection* centralPhysicsSelection = (AliPhysicsSelection*)((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetEventSelection();
5a79fd59 862
507687cd 863 if(!centralPhysicsSelection && !fPhysicsSelection) {
864
865 std::cout<<"Creating AliPhysicsSelection object due to absence of central object"<<std::endl;
866 fPhysicsSelection = new AliPhysicsSelection;
867 fPhysicsSelection->SetAnalyzeMC(!fRealData);
868 // fPhysicsSelection->SetUseBXNumbers(kFALSE);
869
870
871 AliBackgroundSelection* backgroundSelection = new AliBackgroundSelection("bg","bg");
872 backgroundSelection->Init();
873 fPhysicsSelection->AddBackgroundIdentification(backgroundSelection);
874
875 }
5a79fd59 876
507687cd 877 TString triggers = esd->GetFiredTriggerClasses();
5a79fd59 878
059c7c6b 879 AliTriggerAnalysis tAna;
059c7c6b 880
881 fTriggerInel = kFALSE;
882 fTriggerNSD = kFALSE;
883 fTriggerEmpty = kFALSE;
884
507687cd 885 UInt_t inel = kFALSE;
886
887 if(centralPhysicsSelection) {
888 inel= ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
889 }
890 else
891 inel = fPhysicsSelection->IsCollisionCandidate(esd);
892
893
894
895 if(fInelGtZero && inel) {
896 const AliMultiplicity* spdmult = esd->GetMultiplicity();
897 Int_t nCentralTracklets = 0;
898 Int_t j = 0;
899 while( nCentralTracklets < 1 && j< spdmult->GetNumberOfTracklets() ) {
900 if(TMath::Abs(spdmult->GetEta(j)) < 1) nCentralTracklets++;
901 j++;
902 }
903 if(nCentralTracklets < 1) inel = kFALSE;
904 }
905
906 //std::cout<<fTriggerInel<<std::endl;
907 if(inel) {
059c7c6b 908 fTriggerInel = kTRUE;
70d74659 909 }
541c19ed 910
911 Bool_t nsd = kFALSE;
912 if(fUseBuiltInNSD) {
913 if ((tAna.FMDTrigger(esd, AliTriggerAnalysis::kASide) || tAna.V0Trigger(esd, AliTriggerAnalysis::kASide, kFALSE) == AliTriggerAnalysis::kV0BB) && (tAna.FMDTrigger(esd, AliTriggerAnalysis::kCSide) || tAna.V0Trigger(esd, AliTriggerAnalysis::kCSide, kFALSE) == AliTriggerAnalysis::kV0BB))
914 nsd = kTRUE;
915 }
916 else nsd = tAna.IsOfflineTriggerFired(esd,AliTriggerAnalysis::kNSD1);
917
918 if(fTriggerInel && nsd) {
059c7c6b 919 fTriggerNSD = kTRUE;
920 }
921 if(triggers.Contains("CBEAMB-ABCE-NOPF-ALL")) {
922 fTriggerEmpty = kTRUE;
923 }
924
925
926 /*switch (fTrigger) {
9f55be54 927 case kMB1: {
059c7c6b 928 if( fPhysicsSelection->IsCollisionCandidate(esd)) {
929 fTriggerInel = kTRUE;
930 }
931 break;
932
9f55be54 933 }
5a79fd59 934 case kMB2: {
9f55be54 935 if (triggerMask & spdFO && ((triggerMask & v0left) || (triggerMask & v0right)))
936 return kTRUE;
937 break;
938 }
939 case kSPDFASTOR: {
940 if (triggerMask & spdFO)
941 return kTRUE;
942 break;
943 }
50a0adf8 944 case kNOCTP: {
93d093c0 945 return kTRUE;
946 break;
947 }
5a79fd59 948 case kEMPTY: {
059c7c6b 949 if(triggers.Contains("CBEAMB-ABCE-NOPF-ALL"))
950 return kTRUE;
951 break;
952 }
953 case kNSD: {
954 if(fPhysicsSelection->IsCollisionCandidate(esd) && tAna.IsOfflineTriggerFired(esd,AliTriggerAnalysis::kNSD1))
5a79fd59 955 return kTRUE;
956 break;
957 }
da0805e2 958
9f55be54 959 }//switch
059c7c6b 960 */
961
962}
963/*
964//____________________________________________________________________
965Bool_t AliFMDAnaParameters::IsEventTriggered(const AliESDEvent *esd, Trigger trig) {
966 //Did we have trig trigger ?
967 Trigger old = fTrigger;
968 fTrigger = trig;
969 Bool_t retval = IsEventTriggered(esd);
970 fTrigger = old;
971 return retval;
972
973}
974*/
975//____________________________________________________________________
976Bool_t AliFMDAnaParameters::IsEventTriggered(Trigger trigger) {
977 // check if the event was triggered
978
979 if (fCentralSelection) return kTRUE;
980 switch (trigger) {
5a79fd59 981
059c7c6b 982 case kMB1:
983 return fTriggerInel;
984 break;
985 case kNSD:
986 return fTriggerNSD;
987 break;
988 case kEMPTY:
989 return fTriggerEmpty;
990 break;
991 case kNOCTP:
992 return kTRUE;
993 break;
994 default:
995 AliWarning("Trigger not implemented!!!");
996 break;
997
998
999 }
9f55be54 1000 return kFALSE;
059c7c6b 1001
d7346eed 1002}
da0805e2 1003
78f6f750 1004//____________________________________________________________________
1005Float_t
1006AliFMDAnaParameters::GetStripLength(Char_t ring, UShort_t strip)
1007{
aa303f0c 1008 //Get length of a strip
78f6f750 1009
1010 Float_t rad = GetMaxR(ring)-GetMinR(ring);
d05586f1 1011 Float_t nStrips = (ring == 'I' ? 512 : 256);
1012 Float_t segment = rad / nStrips;
78f6f750 1013
1014 //TVector2* corner1 = fmdring.GetVertex(2);
1015 // TVector2* corner2 = fmdring.GetVertex(3);
1016
1017 SetCorners(ring);
1018 /*
1019 std::cout<<GetMaxR(ring)<<" "<<fmdring.GetMaxR()<<std::endl;
1020 std::cout<<GetMinR(ring)<<" "<<fmdring.GetMinR()<<std::endl;
1021 std::cout<<corner1->X()<<" "<<fCorner1.X()<<std::endl;
1022 std::cout<<corner2->X()<<" "<<fCorner2.X()<<std::endl;
1023 std::cout<<corner1->Y()<<" "<<fCorner1.Y()<<std::endl;
1024 std::cout<<corner2->Y()<<" "<<fCorner2.Y()<<std::endl;*/
1025 Float_t slope = (fCorner1.Y() - fCorner2.Y()) / (fCorner1.X() - fCorner2.X());
1026 Float_t constant = (fCorner2.Y()*fCorner1.X()-(fCorner2.X()*fCorner1.Y())) / (fCorner1.X() - fCorner2.X());
1027 Float_t radius = GetMinR(ring) + strip*segment;
1028
1029 Float_t d = TMath::Power(TMath::Abs(radius*slope),2) + TMath::Power(radius,2) - TMath::Power(constant,2);
1030
1031 Float_t arclength = GetBaseStripLength(ring,strip);
1032 if(d>0) {
1033
1034 Float_t x = (-1*TMath::Sqrt(d) -slope*constant) / (1+TMath::Power(slope,2));
1035 Float_t y = slope*x + constant;
1036 Float_t theta = TMath::ATan2(x,y);
1037
1038 if(x < fCorner1.X() && y > fCorner1.Y()) {
1039 arclength = radius*theta; //One sector since theta is by definition half-hybrid
1040
1041 }
1042
1043 }
1044
1045 return arclength;
1046
1047
1048}
1049//____________________________________________________________________
1050Float_t
1051AliFMDAnaParameters::GetBaseStripLength(Char_t ring, UShort_t strip)
1052{
aa303f0c 1053 //Get length of strip assuming that corners are not cut away
78f6f750 1054 Float_t rad = GetMaxR(ring)-GetMinR(ring);
d05586f1 1055 Float_t nStrips = (ring == 'I' ? 512 : 256);
1056 Float_t nSec = (ring == 'I' ? 20 : 40);
1057 Float_t segment = rad / nStrips;
1058 Float_t basearc = 2*TMath::Pi() / (0.5*nSec); // One hybrid: 36 degrees inner, 18 outer
78f6f750 1059 Float_t radius = GetMinR(ring) + strip*segment;
1060 Float_t basearclength = 0.5*basearc * radius; // One sector
1061
1062 return basearclength;
1063}
d7346eed 1064//____________________________________________________________________
ff293341 1065Int_t AliFMDAnaParameters::GetFirstEtaBinFromMap(Int_t vtxbin, Int_t det, Char_t ring)
1066{
1067 TH2F* hBg = GetBackgroundCorrection(det,ring,vtxbin);
2f323ab5 1068
1069 if(det == 0) return hBg->GetXaxis()->FindBin(-1.95);
1070
ff293341 1071 Int_t firstbin = -1;
1072 Int_t nNonZeroFirst = 0;
1073
1074 for(Int_t i=1;i<=hBg->GetNbinsX();i++) {
1075 if(nNonZeroFirst == fNumberOfEtaBinsToCut && firstbin==-1) firstbin = i;
1076
1077 for(Int_t j=1;j<=hBg->GetNbinsY();j++) {
1078
1079 Float_t value = hBg->GetBinContent(i,j);
1080
1081 if(value > 0.001 && nNonZeroFirst<fNumberOfEtaBinsToCut)
1082 {nNonZeroFirst++; break;}
1083
1084
1085 }
1086 }
1087
1088 return firstbin;
1089
1090}
1091//____________________________________________________________________
1092Int_t AliFMDAnaParameters::GetLastEtaBinFromMap(Int_t vtxbin, Int_t det, Char_t ring)
1093{
1094 TH2F* hBg = GetBackgroundCorrection(det,ring,vtxbin);
1095 Int_t lastbin=-1;
1096 Int_t nNonZeroLast = 0;
2f323ab5 1097
1098 if(det == 0) return hBg->GetXaxis()->FindBin(1.95);
1099
ff293341 1100 for(Int_t i=hBg->GetNbinsX();i>0;i--) {
1101 if(nNonZeroLast == fNumberOfEtaBinsToCut && lastbin==-1) lastbin = i;
1102
1103 for(Int_t j=1;j<=hBg->GetNbinsY();j++) {
1104
1105 Float_t value = hBg->GetBinContent(i,j);
1106
1107 if(value > 0.001 && nNonZeroLast<fNumberOfEtaBinsToCut)
1108 {nNonZeroLast++; break; }
1109
1110
1111 }
1112 }
1113
1114 return lastbin;
1115}
1116
1117//____________________________________________________________________
1118Int_t AliFMDAnaParameters::GetFirstEtaBinToInclude(Int_t vtxbin, Int_t det, Char_t ring)
1119{
1120 Int_t ringNumber = (ring == 'I' ? 0 : 1);
1121 return fEtaLowBinLimits.GetBinContent(det,ringNumber,vtxbin);
1122
1123}
1124
1125//____________________________________________________________________
1126Int_t AliFMDAnaParameters::GetLastEtaBinToInclude(Int_t vtxbin, Int_t det, Char_t ring)
1127{
1128 Int_t ringNumber = (ring == 'I' ? 0 : 1);
1129 return fEtaHighBinLimits.GetBinContent(det,ringNumber,vtxbin);
1130
1131}
1132//____________________________________________________________________
d7346eed 1133//
1134// EOF
1135//