1 #include "AliConversionSelection.h"
2 #include "AliAODHeader.h"
4 #include "AliMultiplicity.h"
8 // Author Daniel Lohner (Daniel.Lohner@cern.ch)
12 ClassImp(AliConversionSelection)
15 //________________________________________________________________________
16 AliConversionSelection::AliConversionSelection(AliConversionCuts *convCut, AliConversionMesonCuts *mesonCut) : TObject(),
19 fConversionCut(convCut),
32 fCurrentEventNumber(-1)
35 fInvMassRange=new Double_t[2];
36 fInvMassRange[0]=0.05;
39 //________________________________________________________________________
40 AliConversionSelection::~AliConversionSelection(){
55 delete fPi0Candidates;
68 //________________________________________________________________________
69 Bool_t AliConversionSelection::ProcessEvent(TClonesArray *photons,AliVEvent *inputEvent,AliMCEvent *mcEvent){
70 fInputEvent=inputEvent;
74 Int_t eventnumber=GetEventNumber(inputEvent);
75 if(eventnumber==fCurrentEventNumber){
76 AliWarning("Event already analyzed! Return.");
80 fCurrentEventNumber=eventnumber;
83 // Initialize and Reset Arrays
84 if(fGoodGammas == NULL){
85 fGoodGammas=new TObjArray(30);
89 if(fPi0Candidates == NULL){
90 fPi0Candidates = new TClonesArray("AliAODConversionMother",100);
92 fPi0Candidates->Delete();
95 fBGPi0s = new TClonesArray("AliAODConversionMother",100);
100 if(!photons||!fInputEvent)return kFALSE;
102 if(!fConversionCut->EventIsSelected(fInputEvent,fMCEvent))return kFALSE;
105 for(Int_t i = 0; i < photons->GetEntriesFast(); i++) {
106 AliAODConversionPhoton* gamma =dynamic_cast<AliAODConversionPhoton*>(photons->At(i));
108 if(!fConversionCut->PhotonIsSelected(gamma,fInputEvent))continue;
109 fGoodGammas->Add(gamma);
113 if(fMesonCut->UseMCPSmearing() && fMCEvent){
114 fUnsmearedPx = new Double_t[fGoodGammas->GetEntries()]; // Store unsmeared Momenta
115 fUnsmearedPy = new Double_t[fGoodGammas->GetEntries()];
116 fUnsmearedPz = new Double_t[fGoodGammas->GetEntries()];
117 fUnsmearedE = new Double_t[fGoodGammas->GetEntries()];
119 for(Int_t gamma=0;gamma<fGoodGammas->GetEntries();gamma++){ // Smear the AODPhotons in MC
120 fUnsmearedPx[gamma] = ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->Px();
121 fUnsmearedPy[gamma] = ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->Py();
122 fUnsmearedPz[gamma] = ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->Pz();
123 fUnsmearedE[gamma] = ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->E();
124 fMesonCut->SmearParticle(dynamic_cast<AliAODConversionPhoton*>(fGoodGammas->At(gamma)));
128 // Reconstruct Pi0 and BG
129 CalculatePi0Candidates();
130 CalculateBackground();
131 if(fBGHandler)fBGHandler->AddEvent(fGoodGammas,fInputEvent);
134 if(fMesonCut->UseMCPSmearing() && fMCEvent){
135 for(Int_t gamma=0;gamma<fGoodGammas->GetEntries();gamma++){ // Smear the AODPhotons in MC
136 ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->SetPx(fUnsmearedPx[gamma]); // Reset Unsmeared Momenta
137 ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->SetPy(fUnsmearedPy[gamma]);
138 ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->SetPz(fUnsmearedPz[gamma]);
139 ((AliAODConversionPhoton*)fGoodGammas->At(gamma))->SetE(fUnsmearedE[gamma]);
141 delete[] fUnsmearedPx; fUnsmearedPx = 0x0;
142 delete[] fUnsmearedPy; fUnsmearedPy = 0x0;
143 delete[] fUnsmearedPz; fUnsmearedPz = 0x0;
144 delete[] fUnsmearedE; fUnsmearedE = 0x0;
150 //________________________________________________________________________
151 AliAODConversionMother* AliConversionSelection::GetPi0(Int_t index){
153 if(index>=0&&index<GetNumberOfPi0s()){
154 return dynamic_cast<AliAODConversionMother*>(fPi0Candidates->At(index));
159 //________________________________________________________________________
160 AliAODConversionMother* AliConversionSelection::GetBG(Int_t index){
162 if(index>=0&&index<GetNumberOfBGs()){
163 return dynamic_cast<AliAODConversionMother*>(fBGPi0s->At(index));
168 //________________________________________________________________________
169 AliAODConversionPhoton* AliConversionSelection::GetPhoton(Int_t index){
171 if(index>=0&&index<GetNumberOfPhotons()){
172 return dynamic_cast<AliAODConversionPhoton*>(fGoodGammas->At(index));
177 //________________________________________________________________________
178 void AliConversionSelection::CalculatePi0Candidates(){
181 if(fGoodGammas->GetEntriesFast()>1){
183 for(Int_t firstGammaIndex=0;firstGammaIndex<fGoodGammas->GetEntriesFast()-1;firstGammaIndex++){
185 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGoodGammas->At(firstGammaIndex));
189 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGoodGammas->GetEntriesFast();secondGammaIndex++){
191 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGoodGammas->At(secondGammaIndex));
193 //Check for same Electron ID
194 if(gamma0->GetTrackLabelPositive()==gamma1->GetTrackLabelPositive()||gamma0->GetTrackLabelNegative()==gamma1->GetTrackLabelNegative()
195 ||gamma0->GetTrackLabelNegative()==gamma1->GetTrackLabelPositive()||gamma0->GetTrackLabelPositive()==gamma1->GetTrackLabelNegative())continue;
197 AliAODConversionMother pi0cand(gamma0,gamma1);
198 pi0cand.SetLabels(firstGammaIndex,secondGammaIndex);
204 TParticle *mcgam0=gamma0->GetMCParticle(fMCEvent->Stack());
205 TParticle *mcgam1=gamma1->GetMCParticle(fMCEvent->Stack());
210 if(mcgam0->GetMother(0)==mcgam1->GetMother(0)){
212 pi0cand.SetMCLabel(mcgam0->GetMother(0));
217 if((fMesonCut->MesonIsSelected(&pi0cand))){
218 if(MesonInMassWindow(&pi0cand)){
221 new((*fPi0Candidates)[fPi0Candidates->GetEntriesFast()]) AliAODConversionMother(pi0cand);
229 //________________________________________________________________________
230 Bool_t AliConversionSelection::MesonInMassWindow(AliAODConversionMother *pi0cand)
232 if (pi0cand->M() > fInvMassRange[0] && pi0cand->M() < fInvMassRange[1] ){
238 //________________________________________________________________________
239 void AliConversionSelection::RotateParticle(AliAODConversionPhoton *gamma,Int_t nDegreesPMBackground){
242 fRandomizer=new TRandom3();
243 fRandomizer->SetSeed(0);
245 Double_t nRadiansPM = nDegreesPMBackground*TMath::Pi()/180;
246 Double_t rotationValue = fRandomizer->Rndm()*2*nRadiansPM + TMath::Pi()-nRadiansPM;
247 gamma->RotateZ(rotationValue);
250 //________________________________________________________________________
252 void AliConversionSelection::CalculateBackground(){
255 if(fMesonCut->UseRotationMethod()){
257 // Correct for the number of rotations
258 // BG is for rotation the same, except for factor NRotations
259 Double_t weight=1./Double_t(fMesonCut->NumberOfBGEvents());
261 for(Int_t firstGammaIndex=0;firstGammaIndex<fGoodGammas->GetEntriesFast();firstGammaIndex++){
263 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fGoodGammas->At(firstGammaIndex));
265 for(Int_t secondGammaIndex=firstGammaIndex+1;secondGammaIndex<fGoodGammas->GetEntriesFast();secondGammaIndex++){
266 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fGoodGammas->At(secondGammaIndex));
267 if(!fConversionCut->PhotonIsSelected(gamma1,fInputEvent))continue;
268 for(Int_t nRandom=0;nRandom<fMesonCut->NumberOfBGEvents();nRandom++){
270 RotateParticle(gamma1,fMesonCut->NDegreesRotation());
272 AliAODConversionMother BGcandidate(gamma0,gamma1);
274 if(fMesonCut->MesonIsSelected(&BGcandidate,kFALSE)){
275 if(MesonInMassWindow(&BGcandidate)){
277 new((*fBGPi0s)[fBGPi0s->GetEntriesFast()]) AliAODConversionMother(BGcandidate);
279 dynamic_cast<AliAODConversionMother*>(fBGPi0s->At(fBGPi0s->GetEntriesFast()-1))->SetWeight(weight);
290 if(fBGHandler==NULL){
291 fBGHandler=new AliConversionAODBGHandlerRP(fConversionCut->IsHeavyIon(),fMesonCut->UseTrackMultiplicity());
294 for(Int_t nEventsInBG=0;nEventsInBG <fBGHandler->GetNBGEvents(fGoodGammas,fInputEvent);nEventsInBG++){
296 AliGammaConversionPhotonVector *previousEventGammas = fBGHandler->GetBGGoodGammas(fGoodGammas,fInputEvent,nEventsInBG);
298 if(previousEventGammas){
300 // test weighted background
302 // Correct for the number of eventmixing:
303 // N gammas -> (N-1) + (N-2) +(N-3) ...+ (N-(N-1)) using sum formula sum(i)=N*(N-1)/2 -> N*(N-1)/2
304 // real combinations (since you cannot combine a photon with its own)
305 // but BG leads to N_{a}*N_{b} combinations
306 weight*=0.5*(Double_t(fGoodGammas->GetEntriesFast()-1))/Double_t(previousEventGammas->size());
308 for(Int_t iCurrent=0;iCurrent<fGoodGammas->GetEntriesFast();iCurrent++){
310 AliAODConversionPhoton *gamma0 = (AliAODConversionPhoton*)(fGoodGammas->At(iCurrent));
312 for(UInt_t iPrevious=0;iPrevious<previousEventGammas->size();iPrevious++){
314 AliAODConversionPhoton *gamma1 = (AliAODConversionPhoton*)(previousEventGammas->at(iPrevious));
316 AliAODConversionMother BGcandidate(gamma0,gamma1);
318 if(fMesonCut->MesonIsSelected(&BGcandidate,kFALSE)){
319 if(MesonInMassWindow(&BGcandidate)){
321 new((*fBGPi0s)[fBGPi0s->GetEntriesFast()]) AliAODConversionMother(BGcandidate);
322 dynamic_cast<AliAODConversionMother*>(fBGPi0s->At(fBGPi0s->GetEntriesFast()-1))->SetWeight(weight);
331 //________________________________________________________________________
332 Double_t AliConversionSelection::GetMultiplicity(AliVEvent *inputEvent){
334 switch(fConversionCut->GetMultiplicityMethod())
337 return Double_t(GetNumberOfPhotons());
339 return Double_t(GetNumberOfChargedTracks(inputEvent));
341 return GetVZEROMult(inputEvent);
343 return GetSPDMult(inputEvent);
345 return 1; // if mult is used as a weight, this number can be used to switch off weighting
351 //________________________________________________________________________
352 Int_t AliConversionSelection::GetNumberOfChargedTracks(AliVEvent *inputEvent){
356 AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(inputEvent);
359 fESDTrackCuts= AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
360 fESDTrackCuts->SetMaxDCAToVertexZ(2);
361 Double_t etamax=fConversionCut->GetEtaCut();
362 fESDTrackCuts->SetEtaRange(-etamax, etamax);
363 fESDTrackCuts->SetPtRange(0.15);
365 for(Int_t iTracks = 0; iTracks < inputEvent->GetNumberOfTracks(); iTracks++){
366 AliESDtrack* currentTrack = esdEvent->GetTrack(iTracks);
367 if(!currentTrack) continue;
368 if(fESDTrackCuts->AcceptTrack(currentTrack))ntracks++;
371 for(Int_t ii=0; ii<inputEvent->GetNumberOfTracks(); ii++) {
372 AliVTrack * track = dynamic_cast<AliVTrack*>(inputEvent->GetTrack(ii));
373 if(TMath::Abs(track->Eta())>fConversionCut->GetEtaCut())continue;
381 //________________________________________________________________________
382 Double_t AliConversionSelection::GetVZEROMult(AliVEvent *inputEvent){
384 AliVVZERO *vzero=inputEvent->GetVZEROData();
385 Double_t multV0A=vzero->GetMTotV0A();
386 Double_t multV0C=vzero->GetMTotV0C();
387 Double_t mult=multV0A+multV0C;
392 //________________________________________________________________________
393 Double_t AliConversionSelection::GetSPDMult(AliVEvent *inputEvent){
395 AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(inputEvent);
397 const AliMultiplicity *esdmult=esdEvent->GetMultiplicity();
398 return esdmult->GetNumberOfITSClusters(1);
400 // AOD implementation
401 AliAODHeader *header=(AliAODHeader*)inputEvent->GetHeader();
402 return header->GetNumberOfITSClusters(1);
407 //________________________________________________________________________
408 Int_t AliConversionSelection::GetEventNumber(AliVEvent *inputEvent){
410 AliESDEvent * esdEvent = dynamic_cast<AliESDEvent*>(inputEvent);
412 return esdEvent->GetEventNumberInFile();
415 AliAODHeader *header=(AliAODHeader*)inputEvent->GetHeader();
416 return header->GetEventNumberESDFile();