]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/GammaConv/AliConversionMesonCuts.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
- changes in Conversion & Meson cuts to print full list of cuts during intialization
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliConversionMesonCuts.cxx
CommitLineData
2bb2434e 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Authors: Svein Lindal, Daniel Lohner *
5 * Version 1.0 *
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////////////////////////////////////////////////
4803eb1f 17//---------------------------------------------
2bb2434e 18// Class handling all kinds of selection cuts for
19// Gamma Conversion analysis
20//---------------------------------------------
21////////////////////////////////////////////////
22
23#include "AliConversionMesonCuts.h"
24
25#include "AliKFVertex.h"
26#include "AliAODTrack.h"
27#include "AliESDtrack.h"
28#include "AliAnalysisManager.h"
29#include "AliInputEventHandler.h"
30#include "AliMCEventHandler.h"
31#include "AliAODHandler.h"
32#include "AliPIDResponse.h"
33#include "TH1.h"
34#include "TH2.h"
35#include "AliStack.h"
36#include "AliAODConversionMother.h"
37#include "TObjString.h"
38#include "AliAODEvent.h"
39#include "AliESDEvent.h"
40#include "AliCentrality.h"
41#include "TList.h"
ae947965 42#include "TPDGCode.h"
43#include "TDatabasePDG.h"
44#include "AliAODMCParticle.h"
45
2bb2434e 46class iostream;
47
48using namespace std;
49
50ClassImp(AliConversionMesonCuts)
51
52
53const char* AliConversionMesonCuts::fgkCutNames[AliConversionMesonCuts::kNCuts] = {
0a2b2b4b 54 "MesonKind", //0
55 "BackgroundScheme", //1
56 "NumberOfBGEvents", //2
57 "DegreesForRotationMethod", //3
58 "RapidityMesonCut", //4
59 "RCut", //5
60 "AlphaMesonCut", //6
61 "Chi2MesonCut", //7
62 "SharedElectronCuts", //8
63 "RejectToCloseV0s", //9
64 "UseMCPSmearing", //10
4803eb1f 65 "DcaGammaGamma", //11
66 "DcaRPrimVtx", //12
67 "DcaZPrimVtx" //13
2bb2434e 68};
69
70
71//________________________________________________________________________
e5b6e8a6 72AliConversionMesonCuts::AliConversionMesonCuts(const char *name,const char *title) :
73 AliAnalysisCuts(name,title),
74 fHistograms(NULL),
75 fMesonKind(0),
76 fMaxR(200),
77 fChi2CutMeson(1000),
78 fAlphaMinCutMeson(0),
79 fAlphaCutMeson(1),
80 fRapidityCutMeson(1),
81 fUseRotationMethodInBG(kFALSE),
82 fDoBG(kTRUE),
83 fdoBGProbability(kFALSE),
84 fUseTrackMultiplicityForBG(kFALSE),
85 fnDegreeRotationPMForBG(0),
86 fNumberOfBGEvents(0),
87 fOpeningAngle(0.005),
88 fDoToCloseV0sCut(kFALSE),
89 fminV0Dist(200.),
90 fDoSharedElecCut(kFALSE),
91 fUseMCPSmearing(kFALSE),
92 fPBremSmearing(0),
93 fPSigSmearing(0),
94 fPSigSmearingCte(0),
95 fBrem(NULL),
96 fRandom(0),
976b1f89 97 fFAlphaCut(0),
98 fAlphaPtDepCut(kFALSE),
a280ac15 99 fElectronLabelArraySize(500),
e5b6e8a6 100 fElectronLabelArray(NULL),
4803eb1f 101 fDCAGammaGammaCut(1000),
102 fDCAZMesonPrimVtxCut(1000),
103 fDCARMesonPrimVtxCut(1000),
e5b6e8a6 104 fBackgroundHandler(0),
105 fCutString(NULL),
106 hMesonCuts(NULL),
4803eb1f 107 hMesonBGCuts(NULL),
108 hDCAGammaGammaMesonBefore(NULL),
109 hDCAZMesonPrimVtxBefore(NULL),
110 hDCARMesonPrimVtxBefore(NULL),
111 hDCAGammaGammaMesonAfter(NULL),
112 hDCAZMesonPrimVtxAfter(NULL),
113 hDCARMesonPrimVtxAfter(NULL)
114
2bb2434e 115{
e5b6e8a6 116 for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=0;}
117 fCutString=new TObjString((GetCutNumber()).Data());
a280ac15 118 fElectronLabelArray = new Int_t[fElectronLabelArraySize];
e5b6e8a6 119 if (fBrem == NULL){
120 fBrem = new TF1("fBrem","pow(-log(x),[0]/log(2.0)-1.0)/TMath::Gamma([0]/log(2.0))",0.00001,0.999999999);
121 // tests done with 1.0e-14
122 fBrem->SetParameter(0,fPBremSmearing);
123 fBrem->SetNpx(100000);
124 }
ca91a3e1 125
2bb2434e 126}
127
1d9e6011 128//________________________________________________________________________
129AliConversionMesonCuts::AliConversionMesonCuts(const AliConversionMesonCuts &ref) :
130 AliAnalysisCuts(ref),
131 fHistograms(NULL),
132 fMesonKind(ref.fMesonKind),
133 fMaxR(ref.fMaxR),
134 fChi2CutMeson(ref.fChi2CutMeson),
135 fAlphaMinCutMeson(ref.fAlphaMinCutMeson),
136 fAlphaCutMeson(ref.fAlphaCutMeson),
137 fRapidityCutMeson(ref.fRapidityCutMeson),
138 fUseRotationMethodInBG(ref.fUseRotationMethodInBG),
139 fDoBG(ref.fDoBG),
140 fdoBGProbability(ref.fdoBGProbability),
141 fUseTrackMultiplicityForBG(ref.fUseTrackMultiplicityForBG),
142 fnDegreeRotationPMForBG(ref.fnDegreeRotationPMForBG),
143 fNumberOfBGEvents(ref. fNumberOfBGEvents),
144 fOpeningAngle(ref.fOpeningAngle),
145 fDoToCloseV0sCut(ref.fDoToCloseV0sCut),
146 fminV0Dist(ref.fminV0Dist),
147 fDoSharedElecCut(ref.fDoSharedElecCut),
148 fUseMCPSmearing(ref.fUseMCPSmearing),
149 fPBremSmearing(ref.fPBremSmearing),
150 fPSigSmearing(ref.fPSigSmearing),
151 fPSigSmearingCte(ref.fPSigSmearingCte),
152 fBrem(NULL),
153 fRandom(ref.fRandom),
976b1f89 154 fFAlphaCut(NULL),
155 fAlphaPtDepCut(ref.fAlphaPtDepCut),
a280ac15 156 fElectronLabelArraySize(ref.fElectronLabelArraySize),
1d9e6011 157 fElectronLabelArray(NULL),
4803eb1f 158 fDCAGammaGammaCut(ref.fDCAGammaGammaCut),
159 fDCAZMesonPrimVtxCut(ref.fDCAZMesonPrimVtxCut),
160 fDCARMesonPrimVtxCut(ref.fDCARMesonPrimVtxCut),
1d9e6011 161 fBackgroundHandler(ref.fBackgroundHandler),
162 fCutString(NULL),
163 hMesonCuts(NULL),
4803eb1f 164 hMesonBGCuts(NULL),
165 hDCAGammaGammaMesonBefore(NULL),
166 hDCAZMesonPrimVtxBefore(NULL),
167 hDCARMesonPrimVtxBefore(NULL),
168 hDCAGammaGammaMesonAfter(NULL),
169 hDCAZMesonPrimVtxAfter(NULL),
170 hDCARMesonPrimVtxAfter(NULL)
1d9e6011 171{
a280ac15 172 // Copy Constructor
173 for(Int_t jj=0;jj<kNCuts;jj++){fCuts[jj]=ref.fCuts[jj];}
174 fCutString=new TObjString((GetCutNumber()).Data());
175 fElectronLabelArray = new Int_t[fElectronLabelArraySize];
176 if (fBrem == NULL)fBrem = (TF1*)ref.fBrem->Clone("fBrem");
177 // Histograms are not copied, if you need them, call InitCutHistograms
1d9e6011 178}
179
a280ac15 180
2bb2434e 181//________________________________________________________________________
182AliConversionMesonCuts::~AliConversionMesonCuts() {
e5b6e8a6 183 // Destructor
184 //Deleting fHistograms leads to seg fault it it's added to output collection of a task
185 // if(fHistograms)
186 // delete fHistograms;
187 // fHistograms = NULL;
2bb2434e 188 if(fCutString != NULL){
189 delete fCutString;
190 fCutString = NULL;
191 }
192 if(fElectronLabelArray){
193 delete fElectronLabelArray;
194 fElectronLabelArray = NULL;
195 }
196
197}
198
199//________________________________________________________________________
4803eb1f 200void AliConversionMesonCuts::InitCutHistograms(TString name, Bool_t additionalHists){
e5b6e8a6 201
202 // Initialize Cut Histograms for QA (only initialized and filled if function is called)
ccfa8c0d 203 TH1::AddDirectory(kFALSE);
4803eb1f 204
e5b6e8a6 205 if(fHistograms != NULL){
206 delete fHistograms;
207 fHistograms=NULL;
208 }
209 if(fHistograms==NULL){
210 fHistograms=new TList();
ccfa8c0d 211 fHistograms->SetOwner(kTRUE);
e5b6e8a6 212 if(name=="")fHistograms->SetName(Form("ConvMesonCuts_%s",GetCutNumber().Data()));
213 else fHistograms->SetName(Form("%s_%s",name.Data(),GetCutNumber().Data()));
214 }
215
216 // Meson Cuts
4803eb1f 217 hMesonCuts=new TH1F(Form("MesonCuts %s",GetCutNumber().Data()),"MesonCuts",13,-0.5,12.5);
e5b6e8a6 218 hMesonCuts->GetXaxis()->SetBinLabel(1,"in");
219 hMesonCuts->GetXaxis()->SetBinLabel(2,"undef rapidity");
220 hMesonCuts->GetXaxis()->SetBinLabel(3,"rapidity cut");
221 hMesonCuts->GetXaxis()->SetBinLabel(4,"opening angle");
222 hMesonCuts->GetXaxis()->SetBinLabel(5,"alpha max");
223 hMesonCuts->GetXaxis()->SetBinLabel(6,"alpha min");
4803eb1f 224 hMesonCuts->GetXaxis()->SetBinLabel(7,"dca gamma gamma");
225 hMesonCuts->GetXaxis()->SetBinLabel(8,"dca R prim Vtx");
226 hMesonCuts->GetXaxis()->SetBinLabel(9,"dca Z prim Vtx");
227 hMesonCuts->GetXaxis()->SetBinLabel(10,"out");
e5b6e8a6 228 fHistograms->Add(hMesonCuts);
229
4803eb1f 230 hMesonBGCuts=new TH1F(Form("MesonBGCuts %s",GetCutNumber().Data()),"MesonBGCuts",13,-0.5,12.5);
e5b6e8a6 231 hMesonBGCuts->GetXaxis()->SetBinLabel(1,"in");
232 hMesonBGCuts->GetXaxis()->SetBinLabel(2,"undef rapidity");
233 hMesonBGCuts->GetXaxis()->SetBinLabel(3,"rapidity cut");
234 hMesonBGCuts->GetXaxis()->SetBinLabel(4,"opening angle");
235 hMesonBGCuts->GetXaxis()->SetBinLabel(5,"alpha max");
236 hMesonBGCuts->GetXaxis()->SetBinLabel(6,"alpha min");
4803eb1f 237 hMesonBGCuts->GetXaxis()->SetBinLabel(7,"dca gamma gamma");
238 hMesonBGCuts->GetXaxis()->SetBinLabel(8,"dca R prim Vtx");
239 hMesonBGCuts->GetXaxis()->SetBinLabel(9,"dca Z prim Vtx");
240 hMesonBGCuts->GetXaxis()->SetBinLabel(10,"out");
241
e5b6e8a6 242 fHistograms->Add(hMesonBGCuts);
4803eb1f 243
244 if (additionalHists){
245 hDCAGammaGammaMesonBefore=new TH1F(Form("DCAGammaGammaMeson Before %s",GetCutNumber().Data()),"DCAGammaGammaMeson Before",200,0,10);
246 fHistograms->Add(hDCAGammaGammaMesonBefore);
247
248 hDCARMesonPrimVtxBefore=new TH1F(Form("DCARMesonPrimVtx Before %s",GetCutNumber().Data()),"DCARMesonPrimVtx Before",200,0,10);
249 fHistograms->Add(hDCARMesonPrimVtxBefore);
250
251 hDCAZMesonPrimVtxBefore=new TH1F(Form("DCAZMesonPrimVtx Before %s",GetCutNumber().Data()),"DCAZMesonPrimVtx Before",401,-10,10);
252 fHistograms->Add(hDCAZMesonPrimVtxBefore);
253
254 }
255
256 hDCAGammaGammaMesonAfter=new TH1F(Form("DCAGammaGammaMeson After %s",GetCutNumber().Data()),"DCAGammaGammaMeson After",200,0,10);
257 fHistograms->Add(hDCAGammaGammaMesonAfter);
258
259 hDCAZMesonPrimVtxAfter=new TH2F(Form("InvMassDCAZMesonPrimVtx After %s",GetCutNumber().Data()),"InvMassDCAZMesonPrimVtx After",800,0,0.8,401,-10,10);
260 fHistograms->Add(hDCAZMesonPrimVtxAfter);
261
262 hDCARMesonPrimVtxAfter=new TH1F(Form("DCARMesonPrimVtx After %s",GetCutNumber().Data()),"DCARMesonPrimVtx After",200,0,10);
263 fHistograms->Add(hDCARMesonPrimVtxAfter);
264
ccfa8c0d 265 TH1::AddDirectory(kTRUE);
2bb2434e 266}
267
2bb2434e 268//________________________________________________________________________
11c1e680 269Bool_t AliConversionMesonCuts::MesonIsSelectedMC(TParticle *fMCMother,AliStack *fMCStack, Double_t fRapidityShift){
e5b6e8a6 270 // Returns true for all pions within acceptance cuts for decay into 2 photons
271 // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
4803eb1f 272
e5b6e8a6 273 if(!fMCStack)return kFALSE;
4803eb1f 274
275 if(fMCMother->GetPdgCode()==111 || fMCMother->GetPdgCode()==221){
e5b6e8a6 276 if(fMCMother->R()>fMaxR) return kFALSE; // cuts on distance from collision point
4803eb1f 277
e5b6e8a6 278 Double_t rapidity = 10.;
279 if(fMCMother->Energy() - fMCMother->Pz() == 0 || fMCMother->Energy() + fMCMother->Pz() == 0){
11c1e680 280 rapidity=8.-fRapidityShift;
e5b6e8a6 281 } else{
11c1e680 282 rapidity = 0.5*(TMath::Log((fMCMother->Energy()+fMCMother->Pz()) / (fMCMother->Energy()-fMCMother->Pz())))-fRapidityShift;
4803eb1f 283 }
284
e5b6e8a6 285 // Rapidity Cut
286 if(abs(rapidity)>fRapidityCutMeson)return kFALSE;
4803eb1f 287
e5b6e8a6 288 // Select only -> 2y decay channel
289 if(fMCMother->GetNDaughters()!=2)return kFALSE;
4803eb1f 290
e5b6e8a6 291 for(Int_t i=0;i<2;i++){
292 TParticle *MDaughter=fMCStack->Particle(fMCMother->GetDaughter(i));
293 // Is Daughter a Photon?
294 if(MDaughter->GetPdgCode()!=22)return kFALSE;
295 // Is Photon in Acceptance?
296 // if(bMCDaughtersInAcceptance){
297 // if(!PhotonIsSelectedMC(MDaughter,fMCStack)){return kFALSE;}
298 // }
299 }
300 return kTRUE;
301 }
302 return kFALSE;
2bb2434e 303}
2bb2434e 304//________________________________________________________________________
ae947965 305Bool_t AliConversionMesonCuts::MesonIsSelectedAODMC(AliAODMCParticle *MCMother,TClonesArray *AODMCArray, Double_t fRapidityShift){
e5b6e8a6 306 // Returns true for all pions within acceptance cuts for decay into 2 photons
307 // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
2bb2434e 308
ae947965 309 if(!AODMCArray)return kFALSE;
310
311 if(MCMother->GetPdgCode()==111 || MCMother->GetPdgCode()==221){
312 Double_t rMeson = sqrt( (MCMother->Xv()*MCMother->Xv()) + (MCMother->Yv()*MCMother->Yv()) ) ;
313 if(rMeson>fMaxR) return kFALSE; // cuts on distance from collision point
e5b6e8a6 314
315 Double_t rapidity = 10.;
ae947965 316 if(MCMother->E() - MCMother->Pz() == 0 || MCMother->E() + MCMother->Pz() == 0){
11c1e680 317 rapidity=8.-fRapidityShift;
ae947965 318 } else{
319 rapidity = 0.5*(TMath::Log((MCMother->E()+MCMother->Pz()) / (MCMother->E()-MCMother->Pz())))-fRapidityShift;
e5b6e8a6 320 }
ae947965 321
e5b6e8a6 322 // Rapidity Cut
323 if(abs(rapidity)>fRapidityCutMeson)return kFALSE;
324
ae947965 325 // Select only -> 2y decay channel
326 if(MCMother->GetNDaughters()!=2)return kFALSE;
e5b6e8a6 327
ae947965 328 for(Int_t i=0;i<2;i++){
329 AliAODMCParticle *MDaughter=static_cast<AliAODMCParticle*>(AODMCArray->At(MCMother->GetDaughter(i)));
330 // Is Daughter a Photon?
331 if(MDaughter->GetPdgCode()!=22)return kFALSE;
332 // Is Photon in Acceptance?
333 // if(bMCDaughtersInAcceptance){
334 // if(!PhotonIsSelectedMC(MDaughter,fMCStack)){return kFALSE;}
335 // }
e5b6e8a6 336 }
ae947965 337 return kTRUE;
e5b6e8a6 338 }
339 return kFALSE;
ae947965 340}
209b710e 341
ae947965 342//________________________________________________________________________
343Bool_t AliConversionMesonCuts::MesonIsSelectedMCDalitz(TParticle *fMCMother,AliStack *fMCStack, Int_t &labelelectron, Int_t &labelpositron, Int_t &labelgamma, Double_t fRapidityShift){
4803eb1f 344
ae947965 345 // Returns true for all pions within acceptance cuts for decay into 2 photons
346 // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
347
348 if( !fMCStack )return kFALSE;
4803eb1f 349
ae947965 350 if( fMCMother->GetPdgCode() != 111 && fMCMother->GetPdgCode() != 221 ) return kFALSE;
4803eb1f 351
ae947965 352 if( fMCMother->R()>fMaxR ) return kFALSE; // cuts on distance from collision point
353
354 Double_t rapidity = 10.;
4803eb1f 355
ae947965 356 if( fMCMother->Energy() - fMCMother->Pz() == 0 || fMCMother->Energy() + fMCMother->Pz() == 0 ){
4803eb1f 357 rapidity=8.-fRapidityShift;
ae947965 358 }
359 else{
4803eb1f 360 rapidity = 0.5*(TMath::Log((fMCMother->Energy()+fMCMother->Pz()) / (fMCMother->Energy()-fMCMother->Pz())))-fRapidityShift;
361 }
362
363 // Rapidity Cut
ae947965 364 if( abs(rapidity) > fRapidityCutMeson )return kFALSE;
365
4803eb1f 366 // Select only -> Dalitz decay channel
ae947965 367 if( fMCMother->GetNDaughters() != 3 )return kFALSE;
368
369 TParticle *positron = 0x0;
370 TParticle *electron = 0x0;
371 TParticle *gamma = 0x0;
4803eb1f 372
ae947965 373 for(Int_t index= fMCMother->GetFirstDaughter();index<= fMCMother->GetLastDaughter();index++){
4803eb1f 374
375 TParticle* temp = (TParticle*)fMCStack->Particle( index );
376
377 switch( temp->GetPdgCode() ) {
378 case ::kPositron:
379 positron = temp;
380 labelpositron = index;
381 break;
382 case ::kElectron:
383 electron = temp;
384 labelelectron = index;
385 break;
386 case ::kGamma:
387 gamma = temp;
388 labelgamma = index;
389 break;
390 }
391 }
392
393 if( positron && electron && gamma) return kTRUE;
394 return kFALSE;
209b710e 395}
396
397//________________________________________________________________________
398Bool_t AliConversionMesonCuts::MesonIsSelectedMCEtaPiPlPiMiGamma(TParticle *fMCMother,AliStack *fMCStack, Int_t &labelNegPion, Int_t &labelPosPion, Int_t &labelGamma, Double_t fRapidityShift){
399
400 // Returns true for all pions within acceptance cuts for decay into 2 photons
401 // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
402
403 if( !fMCStack )return kFALSE;
404
405 if( fMCMother->GetPdgCode() != 221 ) return kFALSE;
406
407 if( fMCMother->R()>fMaxR ) return kFALSE; // cuts on distance from collision point
4803eb1f 408
209b710e 409 Double_t rapidity = 10.;
4803eb1f 410
209b710e 411 if( fMCMother->Energy() - fMCMother->Pz() == 0 || fMCMother->Energy() + fMCMother->Pz() == 0 ){
412 rapidity=8.-fRapidityShift;
413 }
414 else{
415 rapidity = 0.5*(TMath::Log((fMCMother->Energy()+fMCMother->Pz()) / (fMCMother->Energy()-fMCMother->Pz())))-fRapidityShift;
416 }
417
418 // Rapidity Cut
419 if( abs(rapidity) > fRapidityCutMeson )return kFALSE;
420
421 // Select only -> Dalitz decay channel
422 if( fMCMother->GetNDaughters() != 3 )return kFALSE;
423
424 TParticle *posPion = 0x0;
425 TParticle *negPion = 0x0;
426 TParticle *gamma = 0x0;
427
428 for(Int_t index= fMCMother->GetFirstDaughter();index<= fMCMother->GetLastDaughter();index++){
429
430 TParticle* temp = (TParticle*)fMCStack->Particle( index );
431
432 switch( temp->GetPdgCode() ) {
433 case 211:
434 posPion = temp;
435 labelPosPion = index;
436 break;
437 case -211:
438 negPion = temp;
439 labelNegPion = index;
440 break;
441 case ::kGamma:
442 gamma = temp;
443 labelGamma = index;
444 break;
445 }
446 }
447
448 if( posPion && negPion && gamma) return kTRUE;
449 return kFALSE;
2bb2434e 450}
209b710e 451
0f8c33c1 452//________________________________________________________________________
453Bool_t AliConversionMesonCuts::MesonIsSelectedMCChiC(TParticle *fMCMother,AliStack *fMCStack,Int_t & labelelectronChiC, Int_t & labelpositronChiC, Int_t & labelgammaChiC, Double_t fRapidityShift){
454 // Returns true for all ChiC within acceptance cuts for decay into JPsi + gamma -> e+ + e- + gamma
455 // If bMCDaughtersInAcceptance is selected, it requires in addition that both daughter photons are within acceptance cuts
456
457 if(!fMCStack)return kFALSE;
4803eb1f 458 // if(fMCMother->GetPdgCode()==20443 ){
459 // return kFALSE;
460 // }
0f8c33c1 461 if(fMCMother->GetPdgCode()==10441 || fMCMother->GetPdgCode()==10443 || fMCMother->GetPdgCode()==445 ){
4803eb1f 462 if(fMCMother->R()>fMaxR) return kFALSE; // cuts on distance from collision point
463
464 Double_t rapidity = 10.;
465 if(fMCMother->Energy() - fMCMother->Pz() == 0 || fMCMother->Energy() + fMCMother->Pz() == 0){
466 rapidity=8.-fRapidityShift;
467 }
468 else{
469 rapidity = 0.5*(TMath::Log((fMCMother->Energy()+fMCMother->Pz()) / (fMCMother->Energy()-fMCMother->Pz())))-fRapidityShift;
470 }
471
472 // Rapidity Cut
473 if(abs(rapidity)>fRapidityCutMeson)return kFALSE;
474
475 // Select only -> ChiC radiative (JPsi+gamma) decay channel
0f8c33c1 476 if(fMCMother->GetNDaughters()!=2)return kFALSE;
477
4803eb1f 478 TParticle *jpsi = 0x0;
479 TParticle *gamma = 0x0;
480 TParticle *positron = 0x0;
481 TParticle *electron = 0x0;
482
483 Int_t labeljpsiChiC = -1;
484
485 for(Int_t index= fMCMother->GetFirstDaughter();index<= fMCMother->GetLastDaughter();index++){
486
487 TParticle* temp = (TParticle*)fMCStack->Particle( index );
488
489 switch( temp->GetPdgCode() ) {
490 case 443:
491 jpsi = temp;
492 labeljpsiChiC = index;
493 break;
494 case 22:
495 gamma = temp;
496 labelgammaChiC = index;
497 break;
498 }
499 }
0f8c33c1 500
501 if ( !jpsi || ! gamma) return kFALSE;
4803eb1f 502 if(jpsi->GetNDaughters()!=2)return kFALSE;
503
504
505 for(Int_t index= jpsi->GetFirstDaughter();index<= jpsi->GetLastDaughter();index++){
506 TParticle* temp = (TParticle*)fMCStack->Particle( index );
507 switch( temp->GetPdgCode() ) {
508 case -11:
509 electron = temp;
510 labelelectronChiC = index;
511 break;
512 case 11:
513 positron = temp;
514 labelpositronChiC = index;
515 break;
516 }
517 }
518 if( !electron || !positron) return kFALSE;
519 if( positron && electron && gamma) return kTRUE;
0f8c33c1 520 }
521 return kFALSE;
522}
523
2bb2434e 524
525
526///________________________________________________________________________
11c1e680 527Bool_t AliConversionMesonCuts::MesonIsSelected(AliAODConversionMother *pi0,Bool_t IsSignal, Double_t fRapidityShift)
2bb2434e 528{
4803eb1f 529
e5b6e8a6 530 // Selection of reconstructed Meson candidates
531 // Use flag IsSignal in order to fill Fill different
532 // histograms for Signal and Background
533 TH1 *hist=0x0;
534
535 if(IsSignal){hist=hMesonCuts;}
536 else{hist=hMesonBGCuts;}
537
538 Int_t cutIndex=0;
539
540 if(hist)hist->Fill(cutIndex);
541 cutIndex++;
542
543 // Undefined Rapidity -> Floating Point exception
544 if((pi0->E()+pi0->Pz())/(pi0->E()-pi0->Pz())<=0){
545 if(hist)hist->Fill(cutIndex);
546 cutIndex++;
d9d6352b 547// cout << "undefined rapidity" << endl;
e5b6e8a6 548 return kFALSE;
549 }
550 else{
551 // PseudoRapidity Cut --> But we cut on Rapidity !!!
552 cutIndex++;
11c1e680 553 if(abs(pi0->Rapidity()-fRapidityShift)>fRapidityCutMeson){
e5b6e8a6 554 if(hist)hist->Fill(cutIndex);
d9d6352b 555// cout << abs(pi0->Rapidity()-fRapidityShift) << ">" << fRapidityCutMeson << endl;
e5b6e8a6 556 return kFALSE;
557 }
558 }
559 cutIndex++;
560
561 // Opening Angle Cut
562 //fOpeningAngle=2*TMath::ATan(0.134/pi0->P());// physical minimum opening angle
563 if( pi0->GetOpeningAngle() < fOpeningAngle){
d9d6352b 564// cout << pi0->GetOpeningAngle() << "<" << fOpeningAngle << endl;
e5b6e8a6 565 if(hist)hist->Fill(cutIndex);
566 return kFALSE;
567 }
568 cutIndex++;
569
976b1f89 570 if ( fAlphaPtDepCut == kTRUE ) {
571
572 fAlphaCutMeson = fFAlphaCut->Eval( pi0->Pt() );
573 }
574
575
e5b6e8a6 576 // Alpha Max Cut
577 if(pi0->GetAlpha()>fAlphaCutMeson){
d9d6352b 578// cout << pi0->GetAlpha() << ">" << fAlphaCutMeson << endl;
e5b6e8a6 579 if(hist)hist->Fill(cutIndex);
580 return kFALSE;
581 }
582 cutIndex++;
583
584 // Alpha Min Cut
585 if(pi0->GetAlpha()<fAlphaMinCutMeson){
d9d6352b 586// cout << pi0->GetAlpha() << "<" << fAlphaMinCutMeson << endl;
e5b6e8a6 587 if(hist)hist->Fill(cutIndex);
588 return kFALSE;
589 }
590 cutIndex++;
591
4803eb1f 592 if (hDCAGammaGammaMesonBefore)hDCAGammaGammaMesonBefore->Fill(pi0->GetDCABetweenPhotons());
593 if (hDCARMesonPrimVtxBefore)hDCARMesonPrimVtxBefore->Fill(pi0->GetDCARMotherPrimVtx());
594
595 if (pi0->GetDCABetweenPhotons() > fDCAGammaGammaCut){
d9d6352b 596// cout << pi0->GetDCABetweenPhotons() << ">" << fDCAGammaGammaCut << endl;
4803eb1f 597 if(hist)hist->Fill(cutIndex);
598 return kFALSE;
599 }
600 cutIndex++;
601
602 if (pi0->GetDCARMotherPrimVtx() > fDCARMesonPrimVtxCut){
d9d6352b 603// cout << pi0->GetDCARMotherPrimVtx() << ">" << fDCARMesonPrimVtxCut << endl;
4803eb1f 604 if(hist)hist->Fill(cutIndex);
605 return kFALSE;
606 }
607 cutIndex++;
608
609
610 if (hDCAZMesonPrimVtxBefore)hDCAZMesonPrimVtxBefore->Fill(pi0->GetDCAZMotherPrimVtx());
611
612 if (abs(pi0->GetDCAZMotherPrimVtx()) > fDCAZMesonPrimVtxCut){
d9d6352b 613// cout << pi0->GetDCAZMotherPrimVtx() << ">" << fDCAZMesonPrimVtxCut << endl;
4803eb1f 614 if(hist)hist->Fill(cutIndex);
615 return kFALSE;
616 }
617 cutIndex++;
618
619
620 if (hDCAGammaGammaMesonAfter)hDCAGammaGammaMesonAfter->Fill(pi0->GetDCABetweenPhotons());
621 if (hDCARMesonPrimVtxAfter)hDCARMesonPrimVtxAfter->Fill(pi0->GetDCARMotherPrimVtx());
622 if (hDCAZMesonPrimVtxAfter)hDCAZMesonPrimVtxAfter->Fill(pi0->M(),pi0->GetDCAZMotherPrimVtx());
623
e5b6e8a6 624 if(hist)hist->Fill(cutIndex);
625 return kTRUE;
2bb2434e 626}
627
628
629
630///________________________________________________________________________
631///________________________________________________________________________
e5b6e8a6 632Bool_t AliConversionMesonCuts::UpdateCutString() {
633 ///Update the cut string (if it has been created yet)
634
635 if(fCutString && fCutString->GetString().Length() == kNCuts) {
e5b6e8a6 636 fCutString->SetString(GetCutNumber());
637 } else {
e5b6e8a6 638 return kFALSE;
639 }
e5b6e8a6 640 return kTRUE;
2bb2434e 641}
642
643///________________________________________________________________________
644Bool_t AliConversionMesonCuts::InitializeCutsFromCutString(const TString analysisCutSelection ) {
e5b6e8a6 645 // Initialize Cuts from a given Cut string
646 AliInfo(Form("Set Meson Cutnumber: %s",analysisCutSelection.Data()));
647 if(analysisCutSelection.Length()!=kNCuts) {
648 AliError(Form("Cut selection has the wrong length! size is %d, number of cuts is %d", analysisCutSelection.Length(), kNCuts));
649 return kFALSE;
650 }
651 if(!analysisCutSelection.IsDigit()){
652 AliError("Cut selection contains characters");
653 return kFALSE;
654 }
4803eb1f 655
e5b6e8a6 656 const char *cutSelection = analysisCutSelection.Data();
657#define ASSIGNARRAY(i) fCuts[i] = cutSelection[i] - '0'
658 for(Int_t ii=0;ii<kNCuts;ii++){
659 ASSIGNARRAY(ii);
660 }
661
662 // Set Individual Cuts
663 for(Int_t ii=0;ii<kNCuts;ii++){
2bb2434e 664 if(!SetCut(cutIds(ii),fCuts[ii]))return kFALSE;
e5b6e8a6 665 }
2bb2434e 666
e656a821 667 PrintCutsWithValues();
e5b6e8a6 668 return kTRUE;
2bb2434e 669}
670///________________________________________________________________________
671Bool_t AliConversionMesonCuts::SetCut(cutIds cutID, const Int_t value) {
e5b6e8a6 672 ///Set individual cut ID
673
674 //cout << "Updating cut " << fgkCutNames[cutID] << " (" << cutID << ") to " << value << endl;
675 switch (cutID) {
676 case kMesonKind:
677 if( SetMesonKind(value)) {
678 fCuts[kMesonKind] = value;
679 UpdateCutString();
680 return kTRUE;
681 } else return kFALSE;
682 case kchi2MesonCut:
683 if( SetChi2MesonCut(value)) {
684 fCuts[kchi2MesonCut] = value;
685 UpdateCutString();
686 return kTRUE;
687 } else return kFALSE;
688 case kalphaMesonCut:
689 if( SetAlphaMesonCut(value)) {
690 fCuts[kalphaMesonCut] = value;
691 UpdateCutString();
692 return kTRUE;
693 } else return kFALSE;
694 case kRCut:
695 if( SetRCut(value)) {
696 fCuts[kRCut] = value;
697 UpdateCutString();
698 return kTRUE;
699 } else return kFALSE;
700
701 case kRapidityMesonCut:
702 if( SetRapidityMesonCut(value)) {
703 fCuts[kRapidityMesonCut] = value;
704 UpdateCutString();
705 return kTRUE;
706 } else return kFALSE;
707
708 case kBackgroundScheme:
709 if( SetBackgroundScheme(value)) {
710 fCuts[kBackgroundScheme] = value;
711 UpdateCutString();
712 return kTRUE;
713 } else return kFALSE;
714
715 case kDegreesForRotationMethod:
716 if( SetNDegreesForRotationMethod(value)) {
717 fCuts[kDegreesForRotationMethod] = value;
718 UpdateCutString();
719 return kTRUE;
720 } else return kFALSE;
721
722 case kNumberOfBGEvents:
723 if( SetNumberOfBGEvents(value)) {
724 fCuts[kNumberOfBGEvents] = value;
725 UpdateCutString();
726 return kTRUE;
727 } else return kFALSE;
728
729 case kuseMCPSmearing:
730 if( SetMCPSmearing(value)) {
731 fCuts[kuseMCPSmearing] = value;
732 UpdateCutString();
733 return kTRUE;
734 } else return kFALSE;
735 case kElecShare:
736 if( SetSharedElectronCut(value)) {
737 fCuts[kElecShare] = value;
738 UpdateCutString();
739 return kTRUE;
740 } else return kFALSE;
741 case kToCloseV0s:
742 if( SetToCloseV0sCut(value)) {
743 fCuts[kToCloseV0s] = value;
744 UpdateCutString();
745 return kTRUE;
746 } else return kFALSE;
4803eb1f 747 case kDcaGammaGamma:
748 if( SetDCAGammaGammaCut(value)) {
749 fCuts[kDcaGammaGamma] = value;
750 UpdateCutString();
751 return kTRUE;
752 } else return kFALSE;
753 case kDcaZPrimVtx:
754 if( SetDCAZMesonPrimVtxCut(value)) {
755 fCuts[kDcaZPrimVtx] = value;
756 UpdateCutString();
757 return kTRUE;
758 } else return kFALSE;
759 case kDcaRPrimVtx:
760 if( SetDCARMesonPrimVtxCut(value)) {
761 fCuts[kDcaRPrimVtx] = value;
762 UpdateCutString();
763 return kTRUE;
764 } else return kFALSE;
e5b6e8a6 765
766 case kNCuts:
767 cout << "Error:: Cut id out of range"<< endl;
768 return kFALSE;
769 }
4803eb1f 770
e5b6e8a6 771 cout << "Error:: Cut id " << cutID << " not recognized "<< endl;
772 return kFALSE;
4803eb1f 773
2bb2434e 774}
775
776
777///________________________________________________________________________
778void AliConversionMesonCuts::PrintCuts() {
e5b6e8a6 779 // Print out current Cut Selection
780 for(Int_t ic = 0; ic < kNCuts; ic++) {
781 printf("%-30s : %d \n", fgkCutNames[ic], fCuts[ic]);
782 }
2bb2434e 783}
784
e656a821 785///________________________________________________________________________
786void AliConversionMesonCuts::PrintCutsWithValues() {
787 // Print out current Cut Selection with values
788 printf("\nMeson cutnumber \n");
789 for(Int_t ic = 0; ic < kNCuts; ic++) {
790 printf("%d",fCuts[ic]);
791 }
792 printf("\n\n");
793
794 printf("Meson cuts \n");
795 printf("\t |y| < %3.2f \n", fRapidityCutMeson);
796 printf("\t theta_{open} < %3.2f\n", fOpeningAngle);
797 if (!fAlphaPtDepCut) printf("\t %3.2f < alpha < %3.2f\n", fAlphaMinCutMeson, fAlphaCutMeson);
798 printf("\t dca_{gamma,gamma} > %3.2f\n", fDCAGammaGammaCut);
799 printf("\t dca_{R, prim Vtx} > %3.2f\n", fDCARMesonPrimVtxCut);
800 printf("\t dca_{Z, prim Vtx} > %3.2f\n\n", fDCAZMesonPrimVtxCut);
801
802 printf("Meson BG settings \n");
803 if (!fDoBG){
804 if (!fUseRotationMethodInBG & !fUseTrackMultiplicityForBG) printf("\t BG scheme: mixing V0 mult \n");
805 if (!fUseRotationMethodInBG & fUseTrackMultiplicityForBG) printf("\t BG scheme: mixing track mult \n");
806 if (fUseRotationMethodInBG )printf("\t BG scheme: rotation \n");
807 if (fdoBGProbability) printf("\t -> use BG probability \n");
808 if (fBackgroundHandler) printf("\t -> use new BG handler \n");
809 printf("\t depth of pool: %d\n", fNumberOfBGEvents);
810 if (fUseRotationMethodInBG )printf("\t degree's for BG rotation: %d\n", fnDegreeRotationPMForBG);
811 }
812
813}
814
815
2bb2434e 816///________________________________________________________________________
817Bool_t AliConversionMesonCuts::SetMesonKind(Int_t mesonKind){
e5b6e8a6 818 // Set Cut
819 switch(mesonKind){
820 case 0:
821 fMesonKind=0;
822 break;
823 case 1:
824 fMesonKind=1;;
825 break;
826 default:
827 cout<<"Warning: Meson kind not defined"<<mesonKind<<endl;
828 return kFALSE;
829 }
830 return kTRUE;
2bb2434e 831}
832
833///________________________________________________________________________
834Bool_t AliConversionMesonCuts::SetRCut(Int_t rCut){
e5b6e8a6 835 // Set Cut
836 switch(rCut){
837 case 0:
838 fMaxR = 180.;
839 break;
840 case 1:
841 fMaxR = 180.;
842 break;
843 case 2:
844 fMaxR = 180.;
845 break;
846 case 3:
4803eb1f 847 fMaxR = 70.;
e5b6e8a6 848 break;
849 case 4:
850 fMaxR = 70.;
851 break;
852 case 5:
853 fMaxR = 180.;
854 break;
855 // High purity cuts for PbPb
856 case 9:
857 fMaxR = 180.;
858 break;
859 default:
860 cout<<"Warning: rCut not defined"<<rCut<<endl;
861 return kFALSE;
862 }
863 return kTRUE;
2bb2434e 864}
865
866///________________________________________________________________________
867Bool_t AliConversionMesonCuts::SetChi2MesonCut(Int_t chi2MesonCut){
e5b6e8a6 868 // Set Cut
869 switch(chi2MesonCut){
870 case 0: // 100.
871 fChi2CutMeson = 100.;
872 break;
873 case 1: // 50.
874 fChi2CutMeson = 50.;
875 break;
876 case 2: // 30.
877 fChi2CutMeson = 30.;
878 break;
879 case 3:
880 fChi2CutMeson = 200.;
881 break;
882 case 4:
883 fChi2CutMeson = 500.;
884 break;
885 case 5:
886 fChi2CutMeson = 1000.;
887 break;
888 default:
889 cout<<"Warning: Chi2MesonCut not defined "<<chi2MesonCut<<endl;
890 return kFALSE;
891 }
892 return kTRUE;
2bb2434e 893}
894
895
896///________________________________________________________________________
897Bool_t AliConversionMesonCuts::SetAlphaMesonCut(Int_t alphaMesonCut)
898{ // Set Cut
e5b6e8a6 899 switch(alphaMesonCut){
900 case 0: // 0- 0.7
901 fAlphaMinCutMeson = 0.0;
902 fAlphaCutMeson = 0.7;
976b1f89 903 fAlphaPtDepCut = kFALSE;
e5b6e8a6 904 break;
976b1f89 905 case 1: // Updated 31 October 2013 before 0.0 - 0.5
906 if( fFAlphaCut ) delete fFAlphaCut;
907 fFAlphaCut= new TF1("fFAlphaCut","[0]*tanh([1]*x)",0.,100.);
908 fFAlphaCut->SetParameter(0,0.7);
909 fFAlphaCut->SetParameter(1,1.2);
910 fAlphaMinCutMeson = 0.0;
911 fAlphaCutMeson = -1.0;
912 fAlphaPtDepCut = kTRUE;
e5b6e8a6 913 break;
976b1f89 914 case 2: // Updated 31 October 2013 before 0.5-1
915 if( fFAlphaCut ) delete fFAlphaCut;
916 fFAlphaCut= new TF1("fFAlphaCut","[0]*tanh([1]*x)",0.,100.);
917 fFAlphaCut->SetParameter(0,0.8);
918 fFAlphaCut->SetParameter(1,1.2);
919 fAlphaMinCutMeson = 0.0;
920 fAlphaCutMeson = -1.0;
921 fAlphaPtDepCut = kTRUE;
e5b6e8a6 922 break;
923 case 3: // 0.0-1
924 fAlphaMinCutMeson = 0.0;
925 fAlphaCutMeson = 1.;
976b1f89 926 fAlphaPtDepCut = kFALSE;
e5b6e8a6 927 break;
928 case 4: // 0-0.65
929 fAlphaMinCutMeson = 0.0;
930 fAlphaCutMeson = 0.65;
976b1f89 931 fAlphaPtDepCut = kFALSE;
e5b6e8a6 932 break;
933 case 5: // 0-0.75
934 fAlphaMinCutMeson = 0.0;
935 fAlphaCutMeson = 0.75;
976b1f89 936 fAlphaPtDepCut = kFALSE;
e5b6e8a6 937 break;
938 case 6: // 0-0.8
939 fAlphaMinCutMeson = 0.0;
940 fAlphaCutMeson = 0.8;
976b1f89 941 fAlphaPtDepCut = kFALSE;
e5b6e8a6 942 break;
943 case 7: // 0.0-0.85
944 fAlphaMinCutMeson = 0.0;
945 fAlphaCutMeson = 0.85;
976b1f89 946 fAlphaPtDepCut = kFALSE;
e5b6e8a6 947 break;
948 case 8: // 0.0-0.6
949 fAlphaMinCutMeson = 0.0;
950 fAlphaCutMeson = 0.6;
976b1f89 951 fAlphaPtDepCut = kFALSE;
e5b6e8a6 952 break;
49af2ef8 953 case 9: // Updated 11 November 2013 before 0.0 - 0.3
954 if( fFAlphaCut ) delete fFAlphaCut;
955 fFAlphaCut= new TF1("fFAlphaCut","[0]*tanh([1]*x)",0.,100.);
956 fFAlphaCut->SetParameter(0,0.65);
957 fFAlphaCut->SetParameter(1,1.2);
958 fAlphaMinCutMeson = 0.0;
959 fAlphaCutMeson = -1.0;
960 fAlphaPtDepCut = kTRUE;
e5b6e8a6 961 break;
962 default:
963 cout<<"Warning: AlphaMesonCut not defined "<<alphaMesonCut<<endl;
964 return kFALSE;
965 }
966 return kTRUE;
2bb2434e 967}
968
969///________________________________________________________________________
4803eb1f 970Bool_t AliConversionMesonCuts::SetRapidityMesonCut(Int_t RapidityMesonCut){
e5b6e8a6 971 // Set Cut
972 switch(RapidityMesonCut){
621d4c94 973 case 0: // changed from 0.9 to 1.35
974 fRapidityCutMeson = 1.35;
e5b6e8a6 975 break;
976 case 1: //
977 fRapidityCutMeson = 0.8;
978 break;
979 case 2: //
980 fRapidityCutMeson = 0.7;
981 break;
982 case 3: //
983 fRapidityCutMeson = 0.6;
984 break;
985 case 4: //
986 fRapidityCutMeson = 0.5;
987 break;
988 case 5: //
989 fRapidityCutMeson = 0.85;
990 break;
991 case 6: //
992 fRapidityCutMeson = 0.75;
993 break;
11c1e680 994 case 7: //
995 fRapidityCutMeson = 0.3;
996 break;
3c56f740 997 case 8: //changed, before 0.35
998 fRapidityCutMeson = 0.25;
11c1e680 999 break;
1000 case 9: //
1001 fRapidityCutMeson = 0.4;
1002 break;
e5b6e8a6 1003 default:
1004 cout<<"Warning: RapidityMesonCut not defined "<<RapidityMesonCut<<endl;
1005 return kFALSE;
1006 }
1007 return kTRUE;
2bb2434e 1008}
1009
1010
1011///________________________________________________________________________
1012Bool_t AliConversionMesonCuts::SetBackgroundScheme(Int_t BackgroundScheme){
e5b6e8a6 1013 // Set Cut
1014 switch(BackgroundScheme){
1015 case 0: //Rotation
1016 fUseRotationMethodInBG=kTRUE;
1017 fdoBGProbability=kFALSE;
1018 break;
1019 case 1: // mixed event with V0 multiplicity
1020 fUseRotationMethodInBG=kFALSE;
1021 fUseTrackMultiplicityForBG=kFALSE;
1022 fdoBGProbability=kFALSE;
1023 break;
1024 case 2: // mixed event with track multiplicity
1025 fUseRotationMethodInBG=kFALSE;
1026 fUseTrackMultiplicityForBG=kTRUE;
1027 fdoBGProbability=kFALSE;
1028 break;
1029 case 3: //Rotation
1030 fUseRotationMethodInBG=kTRUE;
1031 fdoBGProbability=kTRUE;
1032 break;
1033 case 4: //No BG calculation
1034 cout << "no BG calculation should be done" << endl;
1035 fUseRotationMethodInBG=kFALSE;
1036 fdoBGProbability=kFALSE;
1037 fDoBG=kFALSE;
1038 break;
1039 case 5: //Rotation
1040 fUseRotationMethodInBG=kTRUE;
1041 fdoBGProbability=kFALSE;
1042 fBackgroundHandler = 1;
1043 break;
1044 case 6: // mixed event with V0 multiplicity
1045 fUseRotationMethodInBG=kFALSE;
1046 fUseTrackMultiplicityForBG=kFALSE;
1047 fdoBGProbability=kFALSE;
1048 fBackgroundHandler = 1;
1049 break;
1050 case 7: // mixed event with track multiplicity
1051 fUseRotationMethodInBG=kFALSE;
1052 fUseTrackMultiplicityForBG=kTRUE;
1053 fdoBGProbability=kFALSE;
1054 fBackgroundHandler = 1;
1055 break;
1056 case 8: //Rotation
1057 fUseRotationMethodInBG=kTRUE;
1058 fdoBGProbability=kTRUE;
1059 fBackgroundHandler = 1;
1060 break;
1061 default:
1062 cout<<"Warning: BackgroundScheme not defined "<<BackgroundScheme<<endl;
1063 return kFALSE;
1064 }
1065 return kTRUE;
2bb2434e 1066}
1067
ca91a3e1 1068
2bb2434e 1069///________________________________________________________________________
1070Bool_t AliConversionMesonCuts::SetNDegreesForRotationMethod(Int_t DegreesForRotationMethod){
e5b6e8a6 1071 // Set Cut
1072 switch(DegreesForRotationMethod){
1073 case 0:
1074 fnDegreeRotationPMForBG = 5;
1075 break;
1076 case 1:
1077 fnDegreeRotationPMForBG = 10;
1078 break;
1079 case 2:
1080 fnDegreeRotationPMForBG = 15;
1081 break;
1082 case 3:
1083 fnDegreeRotationPMForBG = 20;
1084 break;
1085 default:
1086 cout<<"Warning: DegreesForRotationMethod not defined "<<DegreesForRotationMethod<<endl;
1087 return kFALSE;
1088 }
1089 fCuts[kDegreesForRotationMethod]=DegreesForRotationMethod;
1090 return kTRUE;
2bb2434e 1091}
1092
1093///________________________________________________________________________
ca91a3e1 1094Bool_t AliConversionMesonCuts::SetNumberOfBGEvents(Int_t NumberOfBGEvents){
e5b6e8a6 1095 // Set Cut
1096 switch(NumberOfBGEvents){
1097 case 0:
1098 fNumberOfBGEvents = 5;
1099 break;
1100 case 1:
1101 fNumberOfBGEvents = 10;
1102 break;
1103 case 2:
1104 fNumberOfBGEvents = 15;
1105 break;
1106 case 3:
1107 fNumberOfBGEvents = 20;
1108 break;
1109 case 4:
1110 fNumberOfBGEvents = 2;
1111 break;
1112 case 5:
1113 fNumberOfBGEvents = 50;
1114 break;
1115 case 6:
1116 fNumberOfBGEvents = 80;
1117 break;
1118 case 7:
1119 fNumberOfBGEvents = 100;
1120 break;
1121 default:
1122 cout<<"Warning: NumberOfBGEvents not defined "<<NumberOfBGEvents<<endl;
1123 return kFALSE;
1124 }
1125 return kTRUE;
2bb2434e 1126}
2bb2434e 1127///________________________________________________________________________
1128Bool_t AliConversionMesonCuts::SetSharedElectronCut(Int_t sharedElec) {
1129
e5b6e8a6 1130 switch(sharedElec){
1131 case 0:
1132 fDoSharedElecCut = kFALSE;
1133 break;
1134 case 1:
1135 fDoSharedElecCut = kTRUE;
1136 break;
1137 default:
1138 cout<<"Warning: Shared Electron Cut not defined "<<sharedElec<<endl;
1139 return kFALSE;
1140 }
4803eb1f 1141
e5b6e8a6 1142 return kTRUE;
2bb2434e 1143}
1144
1145///________________________________________________________________________
1146Bool_t AliConversionMesonCuts::SetToCloseV0sCut(Int_t toClose) {
1147
e5b6e8a6 1148 switch(toClose){
1149 case 0:
1150 fDoToCloseV0sCut = kFALSE;
1151 fminV0Dist = 250;
1152 break;
1153 case 1:
1154 fDoToCloseV0sCut = kTRUE;
1155 fminV0Dist = 1;
1156 break;
1157 case 2:
1158 fDoToCloseV0sCut = kTRUE;
1159 fminV0Dist = 2;
1160 break;
1161 case 3:
1162 fDoToCloseV0sCut = kTRUE;
1163 fminV0Dist = 3;
1164 break;
1165 default:
1166 cout<<"Warning: Shared Electron Cut not defined "<<toClose<<endl;
1167 return kFALSE;
1168 }
1169 return kTRUE;
2bb2434e 1170}
1171
ca91a3e1 1172///________________________________________________________________________
1173Bool_t AliConversionMesonCuts::SetMCPSmearing(Int_t useMCPSmearing)
1174{// Set Cut
e5b6e8a6 1175 switch(useMCPSmearing){
1176 case 0:
1177 fUseMCPSmearing=0;
1178 fPBremSmearing=1.;
1179 fPSigSmearing=0.;
1180 fPSigSmearingCte=0.;
1181 break;
1182 case 1:
1183 fUseMCPSmearing=1;
1184 fPBremSmearing=1.0e-14;
1185 fPSigSmearing=0.;
1186 fPSigSmearingCte=0.;
1187 break;
1188 case 2:
1189 fUseMCPSmearing=1;
1190 fPBremSmearing=1.0e-15;
1191 fPSigSmearing=0.0;
1192 fPSigSmearingCte=0.;
1193 break;
1194 case 3:
1195 fUseMCPSmearing=1;
1196 fPBremSmearing=1.;
1197 fPSigSmearing=0.003;
1198 fPSigSmearingCte=0.002;
1199 break;
1200 case 4:
1201 fUseMCPSmearing=1;
1202 fPBremSmearing=1.;
1203 fPSigSmearing=0.003;
1204 fPSigSmearingCte=0.007;
1205 break;
1206 case 5:
1207 fUseMCPSmearing=1;
1208 fPBremSmearing=1.;
1209 fPSigSmearing=0.003;
1210 fPSigSmearingCte=0.016;
1211 break;
1212 case 6:
1213 fUseMCPSmearing=1;
1214 fPBremSmearing=1.;
1215 fPSigSmearing=0.007;
1216 fPSigSmearingCte=0.016;
1217 break;
1218 case 7:
1219 fUseMCPSmearing=1;
1220 fPBremSmearing=1.0e-16;
1221 fPSigSmearing=0.0;
1222 fPSigSmearingCte=0.;
1223 break;
1224 case 8:
1225 fUseMCPSmearing=1;
1226 fPBremSmearing=1.;
1227 fPSigSmearing=0.007;
1228 fPSigSmearingCte=0.014;
1229 break;
1230 case 9:
1231 fUseMCPSmearing=1;
1232 fPBremSmearing=1.;
1233 fPSigSmearing=0.007;
1234 fPSigSmearingCte=0.011;
1235 break;
1236
1237 default:
1238 AliError("Warning: UseMCPSmearing not defined");
1239 return kFALSE;
1240 }
1241 return kTRUE;
ca91a3e1 1242}
4803eb1f 1243
1244
1245///________________________________________________________________________
1246Bool_t AliConversionMesonCuts::SetDCAGammaGammaCut(Int_t DCAGammaGamma){
1247 // Set Cut
1248 switch(DCAGammaGamma){
1249 case 0: //
1250 fDCAGammaGammaCut = 1000;
1251 break;
1252 case 1: //
1253 fDCAGammaGammaCut = 10;
1254 break;
1255 case 2: //
1256 fDCAGammaGammaCut = 5;
1257 break;
1258 case 3: //
1259 fDCAGammaGammaCut = 4;
1260 break;
1261 case 4: //
1262 fDCAGammaGammaCut = 3;
1263 break;
1264 case 5: //
1265 fDCAGammaGammaCut = 2.5;
1266 break;
1267 case 6: //
1268 fDCAGammaGammaCut = 2;
1269 break;
1270 case 7: //
1271 fDCAGammaGammaCut = 1.5;
1272 break;
1273 case 8: //
1274 fDCAGammaGammaCut = 1;
1275 break;
1276 case 9: //
1277 fDCAGammaGammaCut = 0.5;
1278 break;
1279 default:
1280 cout<<"Warning: DCAGammaGamma not defined "<<DCAGammaGamma<<endl;
1281 return kFALSE;
1282 }
1283 return kTRUE;
1284}
1285
1286///________________________________________________________________________
1287Bool_t AliConversionMesonCuts::SetDCAZMesonPrimVtxCut(Int_t DCAZMesonPrimVtx){
1288 // Set Cut
1289 switch(DCAZMesonPrimVtx){
1290 case 0: //
1291 fDCAZMesonPrimVtxCut = 1000;
1292 break;
1293 case 1: //
1294 fDCAZMesonPrimVtxCut = 10;
1295 break;
1296 case 2: //
1297 fDCAZMesonPrimVtxCut = 5;
1298 break;
1299 case 3: //
1300 fDCAZMesonPrimVtxCut = 4;
1301 break;
1302 case 4: //
1303 fDCAZMesonPrimVtxCut = 3;
1304 break;
1305 case 5: //
1306 fDCAZMesonPrimVtxCut = 2.5;
1307 break;
1308 case 6: //
1309 fDCAZMesonPrimVtxCut = 2;
1310 break;
1311 case 7: //
1312 fDCAZMesonPrimVtxCut = 1.5;
1313 break;
1314 case 8: //
1315 fDCAZMesonPrimVtxCut = 1;
1316 break;
1317 case 9: //
1318 fDCAZMesonPrimVtxCut = 0.5;
1319 break;
1320 default:
1321 cout<<"Warning: DCAZMesonPrimVtx not defined "<<DCAZMesonPrimVtx<<endl;
1322 return kFALSE;
1323 }
1324 return kTRUE;
1325}
1326
1327///________________________________________________________________________
1328Bool_t AliConversionMesonCuts::SetDCARMesonPrimVtxCut(Int_t DCARMesonPrimVtx){
1329 // Set Cut
1330 switch(DCARMesonPrimVtx){
1331 case 0: //
1332 fDCARMesonPrimVtxCut = 1000;
1333 break;
1334 case 1: //
1335 fDCARMesonPrimVtxCut = 10;
1336 break;
1337 case 2: //
1338 fDCARMesonPrimVtxCut = 5;
1339 break;
1340 case 3: //
1341 fDCARMesonPrimVtxCut = 4;
1342 break;
1343 case 4: //
1344 fDCARMesonPrimVtxCut = 3;
1345 break;
1346 case 5: //
1347 fDCARMesonPrimVtxCut = 2.5;
1348 break;
1349 case 6: //
1350 fDCARMesonPrimVtxCut = 2;
1351 break;
1352 case 7: //
1353 fDCARMesonPrimVtxCut = 1.5;
1354 break;
1355 case 8: //
1356 fDCARMesonPrimVtxCut = 1;
1357 break;
1358 case 9: //
1359 fDCARMesonPrimVtxCut = 0.5;
1360 break;
1361 default:
1362 cout<<"Warning: DCARMesonPrimVtx not defined "<<DCARMesonPrimVtx<<endl;
1363 return kFALSE;
1364 }
1365 return kTRUE;
1366}
1367
1368
2bb2434e 1369///________________________________________________________________________
1370TString AliConversionMesonCuts::GetCutNumber(){
e5b6e8a6 1371 // returns TString with current cut number
1372 TString a(kNCuts);
1373 for(Int_t ii=0;ii<kNCuts;ii++){
1374 a.Append(Form("%d",fCuts[ii]));
1375 }
1376 return a;
2bb2434e 1377}
1378
1379///________________________________________________________________________
1380void AliConversionMesonCuts::FillElectonLabelArray(AliAODConversionPhoton* photon, Int_t nV0){
4803eb1f 1381
e5b6e8a6 1382 Int_t posLabel = photon->GetTrackLabelPositive();
1383 Int_t negLabel = photon->GetTrackLabelNegative();
4803eb1f 1384
e5b6e8a6 1385 fElectronLabelArray[nV0*2] = posLabel;
1386 fElectronLabelArray[(nV0*2)+1] = negLabel;
2bb2434e 1387}
1388
1389///________________________________________________________________________
1390Bool_t AliConversionMesonCuts::RejectSharedElectronV0s(AliAODConversionPhoton* photon, Int_t nV0, Int_t nV0s){
1391
e5b6e8a6 1392 Int_t posLabel = photon->GetTrackLabelPositive();
1393 Int_t negLabel = photon->GetTrackLabelNegative();
2bb2434e 1394
e5b6e8a6 1395 for(Int_t i = 0; i<nV0s*2;i++){
1396 if(i==nV0*2) continue;
1397 if(i==(nV0*2)+1) continue;
1398 if(fElectronLabelArray[i] == posLabel){
1399 return kFALSE;}
1400 if(fElectronLabelArray[i] == negLabel){
1401 return kFALSE;}
1402 }
2bb2434e 1403
e5b6e8a6 1404 return kTRUE;
2bb2434e 1405}
1406///________________________________________________________________________
1407Bool_t AliConversionMesonCuts::RejectToCloseV0s(AliAODConversionPhoton* photon, TList *photons, Int_t nV0){
e5b6e8a6 1408 Double_t posX = photon->GetConversionX();
1409 Double_t posY = photon->GetConversionY();
1410 Double_t posZ = photon->GetConversionZ();
1411
1412 for(Int_t i = 0;i<photons->GetEntries();i++){
1413 if(nV0 == i) continue;
1414 AliAODConversionPhoton *photonComp = (AliAODConversionPhoton*) photons->At(i);
1415 Double_t posCompX = photonComp->GetConversionX();
1416 Double_t posCompY = photonComp->GetConversionY();
1417 Double_t posCompZ = photonComp->GetConversionZ();
4803eb1f 1418
e5b6e8a6 1419 Double_t dist = pow((posX - posCompX),2)+pow((posY - posCompY),2)+pow((posZ - posCompZ),2);
2bb2434e 1420
e5b6e8a6 1421 if(dist < fminV0Dist*fminV0Dist){
1422 if(photon->GetChi2perNDF() < photonComp->GetChi2perNDF()) return kTRUE;
1423 else {
1424 return kFALSE;}
1425 }
4803eb1f 1426
e5b6e8a6 1427 }
1428 return kTRUE;
2bb2434e 1429}
ca91a3e1 1430
1431///________________________________________________________________________
1432void AliConversionMesonCuts::SmearParticle(AliAODConversionPhoton* photon)
1433{
4803eb1f 1434
0a2b2b4b 1435 if (photon==NULL) return;
ca91a3e1 1436 Double_t facPBrem = 1.;
1437 Double_t facPSig = 0.;
1438
1439 Double_t phi=0.;
1440 Double_t theta=0.;
1441 Double_t P=0.;
1442
4803eb1f 1443
ca91a3e1 1444 P=photon->P();
1445 phi=photon->Phi();
1446 if( photon->P()!=0){
1447 theta=acos( photon->Pz()/ photon->P());
1448 }
1449
4803eb1f 1450 if( fPSigSmearing != 0. || fPSigSmearingCte!=0. ){
ca91a3e1 1451 facPSig = TMath::Sqrt(fPSigSmearingCte*fPSigSmearingCte+fPSigSmearing*fPSigSmearing*P*P)*fRandom.Gaus(0.,1.);
1452 }
4803eb1f 1453
ca91a3e1 1454 if( fPBremSmearing != 1.){
1455 if(fBrem!=NULL){
1456 facPBrem = fBrem->GetRandom();
1457 }
1458 }
1459
1460 photon->SetPx(facPBrem* (1+facPSig)* P*sin(theta)*cos(phi)) ;
1461 photon->SetPy(facPBrem* (1+facPSig)* P*sin(theta)*sin(phi)) ;
1462 photon->SetPz(facPBrem* (1+facPSig)* P*cos(theta)) ;
1463 photon->SetE(photon->P());
1464}
a1c2f90c 1465///________________________________________________________________________
1466void AliConversionMesonCuts::SmearVirtualPhoton(AliAODConversionPhoton* photon)
1467{
1468
1469 if (photon==NULL) return;
1470 Double_t facPBrem = 1.;
1471 Double_t facPSig = 0.;
1472
1473 Double_t phi=0.;
1474 Double_t theta=0.;
1475 Double_t P=0.;
1476
1477
1478 P=photon->P();
1479 phi=photon->Phi();
1480 if( photon->P()!=0){
1481 theta=acos( photon->Pz()/ photon->P());
1482 }
1483
1484 if( fPSigSmearing != 0. || fPSigSmearingCte!=0. ){
1485 facPSig = TMath::Sqrt(fPSigSmearingCte*fPSigSmearingCte+fPSigSmearing*fPSigSmearing*P*P)*fRandom.Gaus(0.,1.);
1486 }
1487
1488 if( fPBremSmearing != 1.){
1489 if(fBrem!=NULL){
1490 facPBrem = fBrem->GetRandom();
1491 }
1492 }
1493
1494 photon->SetPx(facPBrem* (1+facPSig)* P*sin(theta)*cos(phi)) ;
1495 photon->SetPy(facPBrem* (1+facPSig)* P*sin(theta)*sin(phi)) ;
1496 photon->SetPz(facPBrem* (1+facPSig)* P*cos(theta)) ;
1497
7ab8e97e 1498}
1499///________________________________________________________________________
1500TLorentzVector AliConversionMesonCuts::SmearElectron(TLorentzVector particle)
1501{
1502
1503 //if (particle==0) return;
1504 Double_t facPBrem = 1.;
1505 Double_t facPSig = 0.;
1506
1507 Double_t phi=0.;
1508 Double_t theta=0.;
1509 Double_t P=0.;
1510
1511
1512 P=particle.P();
1513
1514
1515 phi=particle.Phi();
1516 if (phi < 0.) phi += 2. * TMath::Pi();
1517
1518 if( particle.P()!=0){
1519 theta=acos( particle.Pz()/ particle.P());
1520 }
1521
1522
1523 Double_t fPSigSmearingHalf = fPSigSmearing / 2.0; //The parameter was set for gammas with 2 particles and here we have just one electron
1524 Double_t sqrtfPSigSmearingCteHalf = fPSigSmearingCte / 2.0 ; //The parameter was set for gammas with 2 particles and here we have just one electron
1525
1526
1527
1528 if( fPSigSmearingHalf != 0. || sqrtfPSigSmearingCteHalf!=0. ){
1529 facPSig = TMath::Sqrt(sqrtfPSigSmearingCteHalf*sqrtfPSigSmearingCteHalf+fPSigSmearingHalf*fPSigSmearingHalf*P*P)*fRandom.Gaus(0.,1.);
1530 }
1531
1532 if( fPBremSmearing != 1.){
1533 if(fBrem!=NULL){
1534 facPBrem = fBrem->GetRandom();
1535 }
1536 }
1537
1538 TLorentzVector SmearedParticle;
1539
1540 SmearedParticle.SetXYZM( facPBrem* (1+facPSig)* P*sin(theta)*cos(phi) , facPBrem* (1+facPSig)* P*sin(theta)*sin(phi) ,
1541 facPBrem* (1+facPSig)* P*cos(theta) , TDatabasePDG::Instance()->GetParticle( ::kElectron )->Mass()) ;
1542
1543 //particle.SetPx(facPBrem* (1+facPSig)* P*sin(theta)*cos(phi)) ;
1544 //particle.SetPy(facPBrem* (1+facPSig)* P*sin(theta)*sin(phi)) ;
1545 //particle.SetPz(facPBrem* (1+facPSig)* P*cos(theta)) ;
1546
1547 return SmearedParticle;
a1c2f90c 1548
1549}
1550