+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnAnalysis
-// Reconstruction and analysis of a binary resonance
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#include <Riostream.h>
-
-#include <TH1.h>
-#include <TTree.h>
-#include <TRefArray.h>
-#include <TParticle.h>
-#include <TObjectTable.h>
-#include <TDatabasePDG.h>
-
-#include "AliRsnDaughter.h"
-#include "AliRsnDaughterCut.h"
-#include "AliRsnEvent.h"
-#include "AliRsnAnalysis.h"
-
-ClassImp(AliRsnAnalysis)
-
-//--------------------------------------------------------------------------------------------------------
-AliRsnAnalysis::AliRsnAnalysis()
-//
-// Constructor
-// Initializes all pointers and collections to NULL.
-//
-{
- fMixHistograms = 0;
- fHistograms = 0;
- fEventsTree = 0;
- fMixPairDefs = 0;
- fPairDefs = 0;
- fPairCuts = 0;
-
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) fCuts[i] = 0;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnAnalysis::AddCutPair(AliRsnDaughterCut *cut)
-//
-// Add a cut on pairs.
-// This cut is global for all pairs.
-//
-{
- if (!cut->IsPairCut()) {
- Warning("AddCutPair", "This is a single cut, cannot be added");
- return;
- }
-
- if (!fPairCuts) fPairCuts = new TObjArray(0);
-
- fPairCuts->AddLast(cut);
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnAnalysis::AddCutSingle(AliPID::EParticleType type, AliRsnDaughterCut *cut)
-//
-// Add a cut on single particles.
-// This cut must be specified for each particle type.
-//
-{
- if (cut->IsPairCut()) {
- Warning("AddCutSingle", "This is a pair cut, cannot be added");
- return;
- }
-
- if (type >= AliPID::kElectron && type <= AliPID::kProton) {
- if (!fCuts[type]) fCuts[type] = new TObjArray(0);
- fCuts[type]->AddLast(cut);
- }
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnAnalysis::AddMixPairDef
-(AliPID::EParticleType p1, Char_t sign1, AliPID::EParticleType p2, Char_t sign2)
-//
-// Adds a new pair definition to create a new histogram,
-// for the event mixing step.
-// If the pair defs array is NULL, it is initialized here.
-//
-{
- if (!fMixPairDefs) fMixPairDefs = new TObjArray(0);
- fMixPairDefs->AddLast( new AliPairDef(p1, sign1, p2, sign2, 0, kFALSE) );
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnAnalysis::AddPairDef
-(AliPID::EParticleType p1, Char_t sign1, AliPID::EParticleType p2, Char_t sign2, Bool_t onlyTrue)
-//
-// Adds a new pair definition to create a new histogram,
-// for the signal evaluation (same event) step.
-// If the pair defs array is NULL, it is initialized here.
-//
-{
- if (!fPairDefs) fPairDefs = new TObjArray(0);
- fPairDefs->AddLast( new AliPairDef(p1, sign1, p2, sign2, fTrueMotherPDG, onlyTrue) );
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnAnalysis::Clear(Option_t* /* option */)
-//
-// Clear heap
-//
-{
- fHistograms->Clear("C");
- fPairDefs->Clear("C");
-
- fMixHistograms->Clear("C");
- fMixPairDefs->Clear("C");
-
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) fCuts[i]->Clear("C");
- fPairCuts->Clear("C");
-}
-//--------------------------------------------------------------------------------------------------------
-Stat_t AliRsnAnalysis::EventMix(Int_t nmix, Int_t multDiffMax, Double_t vzDiffMax, Bool_t compareTotMult)
-//
-// Performs event mixing.
-// It takes the array of fMixPairDefs and stores results in fMixHistograms.
-// First parameter defines how many events must be mixed with each one.
-// Events to be mixed are chosen using the other parameters:
-//
-// - multDiffMax defines the maximum allowed difference in particle multiplicity,
-// a) when last argument is kFALSE (default) the multiplicity comparison is made with respect
-// to the particles of 'type 2' in the pair
-// b) when last argument is kTRUE, the multiplicity comparison is made with respect of total
-// particle multiplicity in the events
-//
-// - vzDiffMax defines maximum allowed difference in Z coordinate of prim. vertex.
-//
-// If one wants to exchange the particle types, one has to add both combinations of particles
-// as different pair defs.
-//
-// EXAMPLE:
-// analysis->AddMixPairDef(AliRsnEvent::kPion, '+', AliRsnEvent::kKaon, '-');
-// analysis->AddMixPairDef(AliRsnEvent::kKaon, '-', AliRsnEvent::kPion, '+');
-//
-{
- // allocate the histograms array
- Int_t i, npairdefs = (Int_t)fMixPairDefs->GetEntries();
- fMixHistograms = new TObjArray(npairdefs);
- TObjArray &histos = *fMixHistograms;
- AliPairDef *pd = 0;
- for (i = 0; i < npairdefs; i++) {
- pd = (AliPairDef*)fMixPairDefs->At(i);
- histos[i] = new TH1D(Form("hmix_%s", pd->GetName()), pd->GetTitle(), fNBins, fHistoMin, fHistoMax);
- }
-
- // Link events branch
- AliRsnEvent *event2 = new AliRsnEvent;
- fEventsTree->SetBranchAddress("events", &event2);
-
- // define variables to store data about particles
- Stat_t nPairs = 0;
- Int_t iev, ievmix, nEvents = (Int_t)fEventsTree->GetEntries();
-
- // loop on events
- Int_t nmixed, mult1, mult2, diffMult;
- Double_t vz1, vz2, diffVz;
- for (iev = 0; iev < nEvents; iev++) {
-
- // get event
- event2->Clear("DELETE");
- fEventsTree->GetEntry(iev);
-
- // copy this event into 'event 1'
- AliRsnEvent *event1 = new AliRsnEvent(*event2);
-
- // get Z position of primary vertex
- vz1 = event1->GetPrimaryVertexZ();
-
- // if total multiplicities must be used
- // it is computed here
- if (compareTotMult) {
- mult1 = event1->GetMultiplicity(kTRUE);
- }
- else {
- mult1 = 0;
- }
-
- // message
- if (iev % 10 == 0) cout << "\rMixing event " << iev << flush;
-
- // loop on pair definitions
- for (i = 0; i < npairdefs; i++) {
-
- // get pair definition
- pd = (AliPairDef*)fMixPairDefs->At(i);
-
- // get histogram reference
- TH1D *histogram = (TH1D*)histos[i];
-
- // get multiplicity of particle type '2' in event '1'
- if (!mult1) {
- mult1 = (Int_t)event1->GetTracks(pd->GetSign2(), pd->GetParticle2())->GetEntries();
- }
-
- // loop on events for mixing
- nmixed = 0;
- ievmix = iev;
- while (nmixed < nmix) {
-
- // get other event (it starts from iev + 1, and
- // loops again to 0 if reachs end of tree)
- ievmix++;
- if (ievmix >= nEvents) ievmix -= nEvents;
-
- // if event-2-mix is equal to 'iev', then
- // a complete loop has been done, and no more
- // mixing can be done, even if they are too few
- // then the loop breaks anyway
- if (ievmix == iev) break;
-
- // get other event
- event2->Clear("DELETE");
- fEventsTree->GetEntry(ievmix);
-
- // get event stats related to event 2
- vz2 = event2->GetPrimaryVertexZ();
- if (compareTotMult) {
- mult2 = event2->GetMultiplicity(kTRUE);
- }
- else {
- mult2 = event2->GetTracks(pd->GetSign2(), pd->GetParticle2())->GetEntries();
- }
-
- // fill histogram if checks are passed
- diffMult = TMath::Abs(mult1 - mult2);
- diffVz = TMath::Abs(vz1 - vz2);
- if (diffVz <= vzDiffMax && diffMult <= multDiffMax) {
- //cout << ievmix << " " << flush;
- nPairs += Compute(pd, histogram, event1, event2);
- nmixed++;
- }
- }
- }
-
- event1->Clear("DELETE");
- delete event1;
- event1 = 0;
-
- } // end loop on events
- cout << endl;
-
- return nPairs;
-}
-//--------------------------------------------------------------------------------------------------------
-Stat_t AliRsnAnalysis::Process()
-//
-// Reads the list 'fPairDefs', and builds an inv-mass histogram for each definition.
-// For each event, particle of 'type 1' are combined with particles of 'type 2' as
-// defined in each pair definition specified in the list, taking the list of both types
-// from the same event.
-// This method is supposed to evaluate the 'signal' of the resonance, containing the peak.
-// It can also be used when one wants to evaluate the 'background' with the 'wrong sign'
-// particle pairs.
-//
-{
- // allocate the histograms array in order to contain
- // as many objects as the number of pair definitionss
- Int_t i, npairdefs = (Int_t)fPairDefs->GetEntries();
- fHistograms = new TObjArray(npairdefs);
- TObjArray &histos = *fHistograms;
-
- // loop over pair defs in order to retrieve the particle species chosen
- // which are used to set an unique name for the output histogram
- // There will be a direct correspondence between the inder of pairdef in its array
- // and the corresponding histogram in the 'fHistograms' list.
- AliPairDef *pd = 0;
- for (i = 0; i < npairdefs; i++) {
- pd = (AliPairDef*)fPairDefs->At(i);
- histos[i] = new TH1D(Form("h_%s", pd->GetName()), pd->GetTitle(), fNBins, fHistoMin, fHistoMax);
- }
-
- // Link events branch
- AliRsnEvent *event = new AliRsnEvent;
- fEventsTree->SetBranchAddress("events", &event);
-
- // define variables to store data about particles
- Stat_t nPairs = 0;
- Int_t iev, nEvents = (Int_t)fEventsTree->GetEntries();
-
- // loop on events
- for (iev = 0; iev < nEvents; iev++) {
-
- // get event
- event->Clear("DELETE");
- fEventsTree->GetEntry(iev);
- if (iev % 1000 == 0) cout << "\rProcessing event " << iev << flush;
-
- // loop over the collection of pair defs
- for (i = 0; i < npairdefs; i++) {
- pd = (AliPairDef*)fPairDefs->At(i);
- TH1D *histogram = (TH1D*)histos[i];
- // invoke AliPairDef method to fill histogram
- nPairs += Compute(pd, histogram, event, event);
- }
-
- } // end loop on events
- cout << endl;
-
- return nPairs;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnAnalysis::WriteHistograms()
-//
-// Writes histograms in current directory
-//
-{
- TH1D *histogram;
- TObjArrayIter iter(fHistograms);
- while ( (histogram = (TH1D*)iter.Next()) ) {
- histogram->Write();
- }
-
- if (fMixHistograms) {
- TObjArrayIter iterMix(fMixHistograms);
- while ( (histogram = (TH1D*)iterMix.Next()) ) {
- histogram->Write();
- }
- }
-}
-//--------------------------------------------------------------------------------------------------------
-Stat_t AliRsnAnalysis::Compute
-(AliPairDef *pd, TH1D* &histogram, AliRsnEvent *event1, AliRsnEvent *event2)
-//
-// Adds to the specified histogram the invariant mass spectrum calculated taking
-// particles of type 1 from event 1, and particles of type 2 from event 2.
-// Events can be equal (signal) or different (background with event mixing).
-//
-{
- // define two 'cursor' objects
- AliRsnDaughter *track1 = 0, *track2 = 0;
-
- // define iterators for the two collections
- if (!event1 || !event1->GetTracks(pd->GetSign1(), pd->GetParticle1())) {
- Error("Compute", "Null pointer for particle 1 collection");
- return 0.0;
- }
- if (!event2 || !event2->GetTracks(pd->GetSign2(), pd->GetParticle2())) {
- Error("Compute", "Null pointer for particle 2 collection");
- return 0.0;
- }
- TObjArrayIter iter1(event1->GetTracks(pd->GetSign1(), pd->GetParticle1()));
- TObjArrayIter iter2(event2->GetTracks(pd->GetSign2(), pd->GetParticle2()));
-
- // define temporary variables for better code readability
- Stat_t nPairs = 0;
-
- // loop on particle of type 1 (in event 1)
- while ( (track1 = (AliRsnDaughter*)iter1.Next()) ) {
-
- if (fRejectFakes && (track1->GetLabel() < 0)) continue;
- if (!SingleCutCheck(pd->GetParticle1(), track1)) continue;
-
- iter2.Reset();
-
- // loop on particles of type 2 (in event 2)
- while ( (track2 = (AliRsnDaughter*)iter2.Next()) ) {
-
- if (fRejectFakes && (track2->GetLabel() < 0)) continue;
- if (event1 == event2 && track1->GetIndex() == track2->GetIndex()) continue;
- if (!SingleCutCheck(pd->GetParticle2(), track2)) continue;
- if (!PairCutCheck(track1, track2)) continue;
-
- /*
- // check
- Char_t sign1 = (track1->GetSign() > 0) ? '+' : '-';
- Char_t sign2 = (track2->GetSign() > 0) ? '+' : '-';
- Info("Compute", "Particle 1: PDG = %d, sign = %c --- Particle 2: PDG = %d, sign = %c", track1->GetTruePDG(), sign1, track2->GetTruePDG(), sign2);
- */
-
- // total 4-momentum
- track1->SetMass(pd->GetMass1());
- track2->SetMass(pd->GetMass2());
- AliRsnDaughter sum = (*track1) + (*track2);
-
- // if the choice to get ONLY true pairs is selected, a check is made on the mothers
- if (pd->GetOnlyTrue()) {
- // the 'sum' AliRsnDaughter object is initialized with
- // the PDG code of the common mother (if it is the case)
- if (TMath::Abs(sum.GetMotherPDG()) == TMath::Abs(fTrueMotherPDG)) {
- histogram->Fill(sum.GetMass());
- nPairs++;
- }
- }
- else {
- histogram->Fill(sum.GetMass());
- nPairs++;
- }
-
- } // end loop 2
-
- } // end loop 1
-
- return nPairs;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnAnalysis::SingleCutCheck(Int_t itype, AliRsnDaughter *track)
-//
-// Checks a track against single particle cuts (if defined)
-//
-{
- if (!fCuts[itype]) return kTRUE;
-
- TObjArrayIter iter(fCuts[itype]);
- AliRsnDaughterCut *cut = 0;
- while ( (cut = (AliRsnDaughterCut*)iter.Next()) ) {
- if (!cut->Pass(track)) return kFALSE;
- }
-
- return kTRUE;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnAnalysis::PairCutCheck(AliRsnDaughter *track1, AliRsnDaughter *track2)
-//
-// Checks a pair against pair cuts (if defined)
-//
-{
- if (!fPairCuts) return kTRUE;
-
- TObjArrayIter iter(fPairCuts);
- AliRsnDaughterCut *cut = 0;
- while ( (cut = (AliRsnDaughterCut*)iter.Next()) ) {
- if (!cut->Pass(track1, track2)) return kFALSE;
- }
-
- return kTRUE;
-}
-//--------------------------------------------------------------------------------------------------------
-AliRsnAnalysis::AliPairDef::AliPairDef
-(AliPID::EParticleType p1, Char_t sign1,
- AliPID::EParticleType p2, Char_t sign2, Int_t pdgMother, Bool_t onlyTrue)
-//
-// Constructor for nested class
-//
-{
- fOnlyTrue = onlyTrue;
- fTrueMotherPDG = 0;
- if (fOnlyTrue) fTrueMotherPDG = pdgMother;
-
- fSign1 = sign1;
- fParticle1 = p1;
-
- fSign2 = sign2;
- fParticle2 = p2;
-
- // instance a PDG database to recovery true masses of particles
- TDatabasePDG *db = TDatabasePDG::Instance();
- Int_t pdg1 = AliPID::ParticleCode((Int_t)p1);
- Int_t pdg2 = AliPID::ParticleCode((Int_t)p2);
- fMass1 = db->GetParticle(pdg1)->Mass();
- fMass2 = db->GetParticle(pdg2)->Mass();
-
- // set name according to the chosen particles
- fName.Append(Form("%s(%c)_%s(%c)", ParticleName(p1), sign1, ParticleName(p2), sign2));
- fTitle.Append(Form("Inv. mass for %s(%c) and %s(%c)", ParticleName(p1), sign1, ParticleName(p2), sign2));
- if (onlyTrue) {
- fName.Append("_true");
- fTitle.Append(" (true pairs)");
- }
-}
-//--------------------------------------------------------------------------------------------------------
-Text_t* AliRsnAnalysis::AliPairDef::ParticleName(AliPID::EParticleType part)
-//
-// [PRIVATE]
-// Returns the name of the particle in text format
-//
-{
- if (part == AliPID::kElectron) return ("E");
- else if (part == AliPID::kMuon) return ("Mu");
- else if (part == AliPID::kPion) return ("Pi");
- else if (part == AliPID::kKaon) return ("K");
- else if (part == AliPID::kProton) return ("P");
- else {
- Warning("ParticleName", "Unrecognized value of EParticle argument");
- return ("???");
- }
-}
-//--------------------------------------------------------------------------------------------------------
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnAnalysis
-// Reconstruction and analysis of K* Rsn
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#ifndef AliRsnANALYSIS_H
-#define AliRsnANALYSIS_H
-
-#include <TNamed.h>
-
-#include "AliPID.h"
-
-class TH1D;
-class TObjArray;
-class AliRsnEvent;
-class AliRsnDaughter;
-class AliRsnDaughterCut;
-
-class AliRsnAnalysis : public TObject
-{
-public:
-
- AliRsnAnalysis();
- virtual ~AliRsnAnalysis() {Clear();}
-
- void AddCutPair(AliRsnDaughterCut *cut);
- void AddCutSingle(AliPID::EParticleType type, AliRsnDaughterCut *cut);
- void AddMixPairDef(AliPID::EParticleType p1, Char_t s1, AliPID::EParticleType p2, Char_t s2);
- void AddPairDef(AliPID::EParticleType p1, Char_t s1, AliPID::EParticleType p2, Char_t s2, Bool_t onlyTrue = kFALSE);
- void Clear(Option_t *option = "");
- Stat_t EventMix(Int_t nmix = 5, Int_t multDiffMax = 10, Double_t vzDiffMax = 0.01, Bool_t compareTotalMult = kFALSE);
- Stat_t Process();
- void SetBins(Int_t nbins, Double_t min, Double_t max) {fNBins=nbins;fHistoMin=min;fHistoMax=max;}
- void SetEventsTree(TTree *tree) {fEventsTree = tree;}
- void SetRejectFakes(Bool_t doit=kTRUE) {fRejectFakes = doit;}
- void SetTrueMotherPDG(Int_t pdg) {fTrueMotherPDG = pdg;}
- void WriteHistograms();
-
-private:
-
- class AliPairDef : public TNamed
- {
- public:
-
- AliPairDef(AliPID::EParticleType p1, Char_t s1, AliPID::EParticleType p2, Char_t s2, Int_t pdg, Bool_t onlyTrue = kFALSE);
-
- virtual ~AliPairDef() { }
-
- Char_t GetSign1() {return fSign1;}
- Char_t GetSign2() {return fSign2;}
- Bool_t GetOnlyTrue() {return fOnlyTrue;}
- AliPID::EParticleType GetParticle1() {return fParticle1;}
- AliPID::EParticleType GetParticle2() {return fParticle2;}
- Double_t GetMass1() {return fMass1;}
- Double_t GetMass2() {return fMass2;}
-
- void SetSign1(Char_t value) {fSign1 = value;}
- void SetSign2(Char_t value) {fSign2 = value;}
- void SetOnlyTrue(Bool_t value = kTRUE) {fOnlyTrue = value;}
- void SetTrueMotherPDG(Int_t pdg) {fTrueMotherPDG = pdg;}
- void SetParticle1(AliPID::EParticleType p) {fParticle1 = p;}
- void SetParticle2(AliPID::EParticleType p) {fParticle2 = p;}
-
- Text_t* ParticleName(AliPID::EParticleType part);
-
- private:
-
- Bool_t fOnlyTrue; // flag to be used for spectra of true pairs
- Int_t fTrueMotherPDG; // PDG code of true mother (if requested)
-
- Double_t fMass1; //
- Char_t fSign1; // info about particle 1
- AliPID::EParticleType fParticle1; //
-
- Double_t fMass2; //
- Char_t fSign2; // info about particle 2
- AliPID::EParticleType fParticle2; //
- };
-
- Stat_t Compute(AliPairDef *pd, TH1D* &h, AliRsnEvent *ev1, AliRsnEvent *ev2);
- Bool_t SingleCutCheck(Int_t itype, AliRsnDaughter *track);
- Bool_t PairCutCheck(AliRsnDaughter *track1, AliRsnDaughter *track2);
-
- Bool_t fRejectFakes; // reject particles labeled as fake
-
- Int_t fNBins; // number of histogram bins
- Double_t fHistoMin; // minimum of the histograms
- Double_t fHistoMax; // maximum of the histograms
-
- Int_t fTrueMotherPDG; // PDG code of true mother (used to create 'true' histos)
-
- TObjArray *fMixPairDefs; // list of pair definitions for histograms (event mixing)
- TObjArray *fMixHistograms; //! list of invmass histograms created (event mixing)
-
- TObjArray *fPairDefs; // list of pair definitions for histograms
- TObjArray *fHistograms; //! list of invmass histograms created
-
- TObjArray *fCuts[AliPID::kSPECIES]; //! list of single particle cuts for each particle type
- TObjArray *fPairCuts; //! list of pair cuts
- TTree *fEventsTree; //! TTree of events (can not be created here, must be passed)
-
- // Rsn analysis implementation
- ClassDef(AliRsnAnalysis,1)
-};
-
-#endif
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnDaughter
-//
-// A simple object which describes a reconstructed track
-// with some references to its related generated particle
-// and some facilities which could help in composing its
-// 4-momentum, for resonance study.
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#include <Riostream.h>
-
-#include <TParticle.h>
-#include <TParticlePDG.h>
-
-#include "AliESDtrack.h"
-#include "AliRsnDaughter.h"
-
-ClassImp(AliRsnDaughter)
-
-//--------------------------------------------------------------------------------------------------------
-AliRsnDaughter::AliRsnDaughter()
-//
-// Default constructor.
-// Its unique argument defines how many PID weights are allowed.
-// Actually, it should be the same as AliESDtrack::kSPECIES (=5).
-//
-{
- fSign = (Char_t)0;
- fPDG = (UShort_t)0;
- fIndex = (UShort_t)0;
-
- fP[0] = fP[1] = fP[2] = 0.0;
- fV[0] = fV[1] = fV[2] = 0.0;
- fMass = 0.0;
-
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- fPIDwgt[i] = 0.0;
- }
-
- fLabel = -1;
- fTruePDG = 0;
- fMother = -1;
- fMotherPDG = 0;
-}
-//--------------------------------------------------------------------------------------------------------
-AliRsnDaughter::AliRsnDaughter(const AliRsnDaughter ©) : TObject(copy)
-//
-// Copy constructor
-//
-{
- fSign = copy.fSign;
- fPDG = copy.fPDG;
- fIndex = copy.fIndex;
-
- Int_t i;
- for (i = 0; i < 3; i++) {
- fP[i] = copy.fP[i];
- fV[i] = copy.fV[i];
- }
- fMass = copy.fMass;
-
- for (i = 0; i < AliPID::kSPECIES; i++) {
- fPIDwgt[i] = copy.fPIDwgt[i];
- }
-
- fLabel = copy.fLabel;
- fTruePDG = copy.fTruePDG;
- fMother = copy.fMother;
- fMotherPDG = copy.fMotherPDG;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnDaughter::Adopt(const AliESDtrack* esdTrack, Bool_t checkITSRefit)
-//
-// Copies reconstructed data from an AliESDtrack:
-//
-// - charge sign
-// - momentum
-// - point of closest approach to primary vertex
-// - ESD pid weights
-// - track label (AliESDtrack::GetLabel())
-//
-// Makes the following checks:
-//
-// - if argument 'checkITSRefit' is TRUE and the "ITS refit" flag is FALSE
-// in the track, the "fIsOK" flag of (this) is set to "false" (track should be rejected)
-// - if the passed label is negative, track is considered as "fake", and
-// this info is kept to allow fake tracks exclusion when doing analysis
-//
-{
- // check for refit in the ITS (if requested)
- if (checkITSRefit) {
- if ( !(esdTrack->GetStatus() & AliESDtrack::kITSrefit) ) {
- return kFALSE;
- }
- }
-
- // get sign and number of species allowed for PID
- fSign = (Char_t)esdTrack->GetSign();
-
- // get (and check) momentum
- esdTrack->GetPxPyPz(fP);
- if (fP[0] == 0.0 || fP[1] == 0.0 || fP[2] == 0.0) {
- return kFALSE;
- }
-
- // get (and check) vertex
- esdTrack->GetXYZ(fV);
- if (fV[0] == 0.0 || fV[1] == 0.0 || fV[2] == 0.0) {
- return kFALSE;
- }
-
- // get label
- // (other kinematics informations are set to default and meaningless values)
- fLabel = esdTrack->GetLabel();
-
- // get PID weights
- esdTrack->GetESDpid(fPIDwgt);
-
- return kTRUE;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnDaughter::Adopt(TParticle* particle)
-//
-// Copies data from a generated particle:
-//
-// - PDG code
-// - charge sign
-// - momentum
-// - production vertex
-// - GEANT label of mother track
-//
-// When an AliRsnDaughter is copied from a TParticle, it is
-// considered always good for analysis and never fake.
-//
-{
- // get particle sign form the sign of PDG code
- Int_t pdg = particle->GetPdgCode();
- if (TMath::Abs(pdg) < 20) {
- if (pdg > 0) fSign = -1; else fSign = 1;
- }
- else if (TMath::Abs(pdg) < 3000) {
- if (pdg > 0) fSign = 1; else fSign = -1;
- }
-
- // get momentum
- fP[0] = particle->Px();
- fP[1] = particle->Py();
- fP[2] = particle->Pz();
-
- // get vertex
- fV[0] = particle->Vx();
- fV[1] = particle->Vy();
- fV[2] = particle->Vz();
-
- // set simulation data
- fPDG = fTruePDG = (Short_t)particle->GetPdgCode();
- fMother = particle->GetFirstMother();
- fMotherPDG = 0;
-
- return kTRUE;
-}
-//--------------------------------------------------------------------------------------------------------
-AliRsnDaughter operator+(AliRsnDaughter t1, AliRsnDaughter t2)
-//
-// Sum operator overloading.
-// Builds a new AliRsnDaughter object with the sum of momenta of two particles.
-//
-{
- // create new AliRsnDaughter with default useless values for datamembers
- AliRsnDaughter out;
-
- // if the summed particles are daughters of the same resonance
- // their common mother label becomes the label of the sum
- Int_t mum1 = t1.GetMother();
- Int_t mum2 = t2.GetMother();
- if (mum1 == mum2) {
- out.SetLabel(mum1);
- out.SetMotherPDG(t1.GetMotherPDG());
- }
- else {
- out.SetLabel(-1);
- out.SetMotherPDG(0);
- }
-
- // compute total 4-momentum
- Double_t Etot = t1.GetEnergy() + t2.GetEnergy();
- Double_t pxTot = t1.GetPx() + t2.GetPx();
- Double_t pyTot = t1.GetPy() + t2.GetPy();
- Double_t pzTot = t1.GetPz() + t2.GetPz();
- Double_t mass = TMath::Sqrt(Etot*Etot - pxTot*pxTot - pyTot*pyTot - pzTot*pzTot);
-
- //TLorentzVector v1 = track1.Get4Momentum();
- //TLorentzVector v2 = track2.Get4Momentum();
- //TLorentzVector sm = v1 + v2;
- //Double_t pxTot = sum.X();
- //Double_t pyTot = sum.Y();
- //Double_t pzTot = sum.Z();
- //Double_t mass = sm.M();
-
- out.SetPxPyPz(pxTot, pyTot, pzTot);
- out.SetMass(mass);
-
- return out;
-}
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnDaughter
-//
-// A simple object which describes a reconstructed track
-// with some references to its related generated particle
-// and some facilities which could help in composing its
-// 4-momentum, for resonance study.
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#ifndef ALIRSNDAUGHTER_H
-#define ALIRSNDAUGHTER_H
-
-#include "AliPID.h"
-#include <TVector3.h>
-#include <TParticle.h>
-#include <TLorentzVector.h>
-
-class AliESDtrack;
-
-class AliRsnDaughter : public TObject
-{
-public:
- AliRsnDaughter();
- AliRsnDaughter(const AliRsnDaughter ©);
-
- virtual ~AliRsnDaughter() { }
-
- Bool_t Adopt(TParticle* particle);
- Bool_t Adopt(const AliESDtrack* track, Bool_t checkRefit = kTRUE);
- TVector3 Get3Momentum() const {TVector3 v(fP[0],fP[1],fP[2]); return v;}
- TLorentzVector Get4Momentum() const {TLorentzVector v(fP[0],fP[1],fP[2],GetEnergy()); return v;}
- Double_t GetEnergy() const {return TMath::Sqrt(fMass*fMass + GetP2());}
- UShort_t GetIndex() const {return fIndex;}
- Int_t GetLabel() const {return fLabel;}
- Double_t GetMass() const {return fMass;}
- Int_t GetMother() const {return fMother;}
- Short_t GetMotherPDG() const {return fMotherPDG;}
- UShort_t GetPDG() const {return fPDG;}
- Double_t GetPIDweight(Int_t i) const {return ( (i>=0&&i<AliPID::kSPECIES)?fPIDwgt[i]:-1.0 );}
- Char_t GetSign() const {return fSign;}
- Double_t GetP2() const {return fP[0]*fP[0] + fP[1]*fP[1] + fP[2]*fP[2];}
- Double_t GetP() const {return TMath::Sqrt(GetP2());}
- Double_t GetPx() const {return fP[0];}
- Double_t GetPy() const {return fP[1];}
- Double_t GetPz() const {return fP[2];}
- Double_t GetPt() const {return TMath::Sqrt(fP[0]*fP[0] + fP[1]*fP[1]);}
- Short_t GetTruePDG() const {return fTruePDG;}
- TVector3 GetVertex() const {TVector3 v(fV[0],fV[1],fV[2]); return v;}
- Double_t GetVx() const {return fV[0];}
- Double_t GetVy() const {return fV[1];}
- Double_t GetVz() const {return fV[2];}
- Double_t GetVt() const {return TMath::Sqrt(fV[0]*fV[0] + fV[1]*fV[1]);}
- void SetIndex(UShort_t value) {fIndex = value;}
- void SetIndex(Int_t value) {fIndex = (UShort_t)value;}
- void SetLabel(Int_t l) {fLabel = l;}
- void SetMass(Double_t m) {fMass = m;}
- void SetMother(Int_t l) {fMother = l;}
- void SetMotherPDG(Short_t pdg) {fMotherPDG = pdg;}
- void SetPDG(UShort_t pdg) {fPDG = TMath::Abs(pdg);}
- void SetPDG(Int_t pdg) {fPDG = (UShort_t)TMath::Abs(pdg);}
- void SetPIDweights(const Double_t *pid) {Int_t i;for(i=0;i<AliPID::kSPECIES;i++)fPIDwgt[i]=pid[i];}
- void SetPxPyPz(Double_t px,Double_t py,Double_t pz) {fP[0]=px;fP[1]=py;fP[2]=pz;}
- void SetSign(Char_t value) {fSign = value;}
- void SetSign(Int_t value) {fSign = (Char_t)value;}
- void SetTruePDG(Short_t pdg) {fTruePDG = pdg;}
- void SetVxVyVz(Double_t vx,Double_t vy,Double_t vz) {fV[0]=vx;fV[1]=vy;fV[2]=vz;}
-
- friend AliRsnDaughter operator+(AliRsnDaughter t1, AliRsnDaughter t2);
-
-private:
-
- Char_t fSign; // charge sign
- UShort_t fPDG; // assigned PDG code from PID (0=undefined)
- UShort_t fIndex; // reference index in AliESD container
-
- Double_t fP[3]; // vector momentum
- Double_t fV[3]; // production vertex
- Double_t fMass; // mass
-
- Double_t fPIDwgt[AliPID::kSPECIES]; // particle PID weights
-
- // The following data are useful for simulated events only
- // and can be left blank when using real or 'realistic' data
-
- Int_t fLabel; // GEANT label of corresponding particle
- Short_t fTruePDG; // PDG code of corresponding particle
- Int_t fMother; // GEANT label of mother particle
- Short_t fMotherPDG; // PDG code of mother particle
-
- ClassDef(AliRsnDaughter,1)
-};
-
-#endif
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnDaughterCut
-//
-// Implementation of track cuts for analysis.
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#include <Riostream.h>
-
-#include "AliRsnDaughter.h"
-#include "AliRsnDaughterCut.h"
-
-ClassImp(AliRsnDaughterCut)
-
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnDaughterCut::Pass(AliRsnDaughter* /*track1*/, AliRsnDaughter* /*track2*/)
-//
-// Virtual method for cut passing.
-// This function must be overridden and return kTRUE when cut is passed.
-//
-{
- TObject::Error("Pass", "This method must be overridden!");
- return kFALSE;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnDaughterCutPtSingle::Pass(AliRsnDaughter *track1, AliRsnDaughter* /*track2*/)
-//
-// Cut on single track momentum.
-//
-{
- if (!track1) return kFALSE;
- if (track1->GetPt() < fPtMin) return kFALSE;
- if (track1->GetPt() > fPtMax) return kFALSE;
-
- return kTRUE;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnDaughterCutPtPair::Pass(AliRsnDaughter *track1, AliRsnDaughter *track2)
-//
-// Cut on single track momentum.
-//
-{
- if (!track1 || !track2) return kFALSE;
-
- AliRsnDaughter sum = (*track1) + (*track2);
-
- if (sum.GetPt() < fPtMin) return kFALSE;
- if (sum.GetPt() > fPtMax) return kFALSE;
-
- return kTRUE;
-}
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnDaughterCut
-//
-// Implementation of various cut which can be applied
-// during resonance analysis.
-// First is the virtual base class.
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#ifndef ALIRSNDAUGHTERCUT_H
-#define ALIRSNDAUGHTERCUT_H
-
-class AliRsnDaughter;
-
-class AliRsnDaughterCut : public TObject
-{
-public:
- AliRsnDaughterCut() { }
- virtual ~AliRsnDaughterCut() { }
-
- Bool_t IsPairCut() {return fPairCut;}
- virtual Bool_t Pass(AliRsnDaughter *track1, AliRsnDaughter *track2 = 0);
-
-protected:
-
- Bool_t fPairCut;
-
- ClassDef(AliRsnDaughterCut,1)
-};
-
-//-------------------------------------------------------------------------
-
-class AliRsnDaughterCutPtSingle : public AliRsnDaughterCut
-{
-public:
- AliRsnDaughterCutPtSingle(Double_t min, Double_t max) {fPairCut=kFALSE;fPtMin=min;fPtMax=max;}
- virtual ~AliRsnDaughterCutPtSingle() { }
-
- virtual Bool_t Pass(AliRsnDaughter *track1, AliRsnDaughter *track2 = 0);
-
-protected:
-
- Double_t fPtMin; // smallest allowed Pt
- Double_t fPtMax; // largest allowed Pt
-
- ClassDef(AliRsnDaughterCutPtSingle,1)
-};
-
-//-------------------------------------------------------------------------
-
-class AliRsnDaughterCutPtPair : public AliRsnDaughterCut
-{
-public:
- AliRsnDaughterCutPtPair(Double_t min, Double_t max) {fPairCut=kTRUE;fPtMin=min;fPtMax=max;}
- virtual ~AliRsnDaughterCutPtPair() { }
-
- virtual Bool_t Pass(AliRsnDaughter *track1, AliRsnDaughter *track2);
-
-protected:
-
- Double_t fPtMin; // smallest allowed Pt
- Double_t fPtMax; // largest allowed Pt
-
- ClassDef(AliRsnDaughterCutPtPair,1)
-};
-
-//-------------------------------------------------------------------------
-
-#endif
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnEvent
-// Simple collection of reconstructed tracks, selected from an ESD event
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#include <Riostream.h>
-
-#include <TTree.h>
-#include <TString.h>
-#include <TParticle.h>
-#include <TObjArray.h>
-#include <TRefArray.h>
-#include <TClonesArray.h>
-
-#include "AliESD.h"
-#include "AliStack.h"
-#include "AliESDtrack.h"
-#include "AliRsnDaughter.h"
-#include "AliRsnEvent.h"
-
-ClassImp(AliRsnEvent)
-
-//--------------------------------------------------------------------------------------------------------
-AliRsnEvent::AliRsnEvent() :
- fIsESD(kTRUE),
- fPath(""),
- fPVx(0.0),
- fPVy(0.0),
- fPVz(0.0),
- fMultiplicity(-1)
-//
-// Default constructor
-//
-{
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- fPos[i] = NULL;
- fNeg[i] = NULL;
- }
-}
-//--------------------------------------------------------------------------------------------------------
-AliRsnEvent::AliRsnEvent(const AliRsnEvent &event) :
- TObject((TObject)event),
- fIsESD(event.fIsESD),
- fPath(event.fPath),
- fPVx(event.fPVx),
- fPVy(event.fPVy),
- fPVz(event.fPVz),
- fMultiplicity(event.fMultiplicity)
-//
-// Copy constructor.
-// Creates new instances of all collections to store a copy of all objects.
-//
-{
- // clone tracks collections
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- fPos[i] = 0;
- fNeg[i] = 0;
- if (event.fPos[i]) fPos[i] = (TClonesArray*)event.fPos[i]->Clone();
- if (event.fNeg[i]) fNeg[i] = (TClonesArray*)event.fNeg[i]->Clone();
- }
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnEvent::AddTrack(AliRsnDaughter track)
-//
-// Stores a track into the correct array
-//
-{
- Int_t pdg, sign, ifirst, ilast;
-
- // if sign is zero, track is not stored
- sign = (Int_t)track.GetSign();
- if (!sign) return;
-
- // if PDG code is assigned, track is stored in the corresponding collection
- // otherwise, a copy of track is is stored in each collection (undefined track)
- pdg = track.GetPDG();
- if (pdg != 0) {
- ifirst = PDG2Enum(pdg);
- ilast = ifirst;
- if (ifirst < AliPID::kElectron || ifirst > AliPID::kProton) return;
- }
- else {
- ifirst = AliPID::kElectron;
- ilast = AliPID::kProton;
- }
-
- // track is stored
- Int_t i, index;
- for (i = ifirst; i <= ilast; i++) {
- if (sign > 0) {
- index = fPos[i]->GetEntries();
- TClonesArray &array = *fPos[i];
- new(array[index]) AliRsnDaughter(track);
- }
- else {
- index = fNeg[i]->GetEntries();
- TClonesArray &array = *fNeg[i];
- new(array[index]) AliRsnDaughter(track);
- }
- }
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnEvent::Clear(Option_t *option)
-//
-// Clears list of tracks and references.
-// If the string "DELETE" is specified, the collection classes
-// are also cleared from heap.
-//
-{
- // evaluate option
- TString opt(option);
- Bool_t deleteCollections = opt.Contains("DELETE", TString::kIgnoreCase);
-
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- if (fPos[i]) fPos[i]->Delete();
- if (fNeg[i]) fNeg[i]->Delete();
- if (deleteCollections) {
- delete fPos[i];
- delete fNeg[i];
- fPos[i] = 0;
- fNeg[i] = 0;
- }
- }
-}
-//--------------------------------------------------------------------------------------------------------
-Int_t AliRsnEvent::GetMultiplicity(Bool_t recalc)
-//
-// Computes multiplicity.
-// If it is already computed (fMultiplicity > -1), it returns that value,
-// unless one sets the argument to kTRUE.
-//
-{
- if (fMultiplicity < 0) recalc = kTRUE;
- if (recalc) {
- fMultiplicity = 0;
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- fMultiplicity += (Int_t)fPos[i]->GetEntries();
- fMultiplicity += (Int_t)fNeg[i]->GetEntries();
- }
- }
-
- return fMultiplicity;
-}
-//--------------------------------------------------------------------------------------------------------
-const char * AliRsnEvent::GetOriginFileName()
-//
-// Returns the path where input file was stored
-//
-{
- TString str(fPath);
- if (fIsESD) {
- str.Append("/AliESDs.root");
- }
- else {
- str.Append("/galice.root");
- }
-
- return str.Data();
-}
-//--------------------------------------------------------------------------------------------------------
-TClonesArray* AliRsnEvent::GetTracks(Char_t sign, AliPID::EParticleType type)
-//
-// Returns the particle collection specified in argument
-//
-{
- Int_t itype = (Int_t)type;
- if (itype >= 0 && itype < AliPID::kSPECIES) {
- if (sign == '+') return fPos[type]; else return fNeg[type];
- }
- else {
- return NULL;
- }
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnEvent::Init()
-//
-// Action 1: define default values for some data members (including pointers).
-// Action 2: if 'ntracks' > 0 allocates memory to store tracks.
-//
-{
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- fPos[i] = new TClonesArray("AliRsnDaughter", 0);
- fNeg[i] = new TClonesArray("AliRsnDaughter", 0);
- fPos[i]->BypassStreamer(kFALSE);
- fNeg[i]->BypassStreamer(kFALSE);
- }
-}
-//--------------------------------------------------------------------------------------------------------
-Int_t AliRsnEvent::PDG2Enum(Int_t pdgcode)
-//
-// Converts a PDG code into the correct slot in the EParticleType enumeration in AliPID
-//
-{
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- if (AliPID::ParticleCode((AliPID::EParticleType)i) == TMath::Abs(pdgcode)) {
- return i;
- }
- }
-
- return -1;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnEvent::PrintTracks()
-//
-// Print data for particles stored in this event
-//
-{
- cout << endl;
-
- AliRsnDaughter *track = 0;
-
- Int_t type;
- for (type = 0; type < AliPID::kSPECIES; type++) {
- TObjArrayIter iterPos(fPos[type]);
- cout << "Positive " << AliPID::ParticleName(type) << "s" << endl;
- if (!fPos[type]) cout << "NOT INITIALIZED" << endl;
- while ( (track = (AliRsnDaughter*)iterPos.Next()) ) {
- cout << "Index, Sign, PDG = ";
- cout << track->GetIndex() << " ";
- cout << (Int_t)track->GetSign() << " ";
- cout << track->GetPDG() << endl;
- }
- TObjArrayIter iterNeg(fNeg[type]);
- cout << "Negative " << AliPID::ParticleName(type) << "s" << endl;
- if (!fNeg[type]) cout << "NOT INITIALIZED" << endl;
- while ( (track = (AliRsnDaughter*)iterNeg.Next()) ) {
- cout << "Index, Sign, PDG = ";
- cout << track->GetIndex() << " ";
- cout << (Int_t)track->GetSign() << " ";
- cout << track->GetPDG() << endl;
- }
- }
-}
-//--------------------------------------------------------------------------------------------------------
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnEvent
-// Simple collection of reconstructed tracks, selected from an ESD event
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#ifndef ALIRSNEVENT_H
-#define ALIRSNEVENT_H
-
-#include <TMath.h>
-#include "AliPID.h"
-
-class AliRsnDaughter;
-
-class AliRsnEvent : public TObject
-{
-public:
- AliRsnEvent();
- AliRsnEvent(const AliRsnEvent& copy);
-
- virtual ~AliRsnEvent() {Clear("DELETE");}
-
- void AddTrack(AliRsnDaughter track);
- void Clear(Option_t *option = "");
- Int_t GetMultiplicity(Bool_t recalc = kFALSE);
- const char* GetOriginFileName();
- Double_t GetPrimaryVertexX() {return fPVx;}
- Double_t GetPrimaryVertexY() {return fPVy;}
- Double_t GetPrimaryVertexZ() {return fPVz;}
- void GetPrimaryVertex(Double_t &x, Double_t &y, Double_t &z) {x=fPVx;y=fPVy;z=fPVz;}
- TClonesArray* GetTracks(Char_t sign, AliPID::EParticleType type);
- void Init();
- Int_t PDG2Enum(Int_t pdgcode);
- void PrintTracks();
- void SetESD(Bool_t yesno = kTRUE) {fIsESD=yesno;}
- void SetPath(TString path) {fPath=path;}
- void SetPrimaryVertex(Double_t x, Double_t y, Double_t z) {fPVx=x;fPVy=y;fPVz=z;}
-
-private:
-
- Bool_t fIsESD; // if true, it is ESD event, otherwise it comes from Kine
- TString fPath; // complete path where input event file is stored
-
- Double_t fPVx; //
- Double_t fPVy; // primary vertex
- Double_t fPVz; //
-
- Int_t fMultiplicity; // global event multiplicity
-
- TClonesArray *fPos[AliPID::kSPECIES]; // collections of positive particles
- TClonesArray *fNeg[AliPID::kSPECIES]; // collections of negative particles
-
- ClassDef(AliRsnEvent,1);
-};
-
-#endif
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * *
- * Author: The ALICE Off-line Project. *
- * Contributors are mentioned in the code where appropriate. *
- * *
- * Permission to use, copy, modify and distribute this software and its *
- * documentation strictly for non-commercial purposes is hereby granted *
- * without fee, provided that the above copyright notice appears in all *
- * copies and that both the copyright notice and this permission notice *
- * appear in the supporting documentation. The authors make no claims *
- * about the suitability of this software for any purpose. It is *
- * provided "as is" without express or implied warranty. *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnReader
-//
-// Reader for conversion of ESD or Kinematics output into AliRsnEvent
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#include <Riostream.h>
-
-#include <TH1.h>
-#include <TH3.h>
-#include <TFile.h>
-#include <TChain.h>
-#include <TArrayF.h>
-#include <TParticle.h>
-#include <TRandom.h>
-#include <TObjString.h>
-#include <TObjectTable.h>
-#include <TOrdCollection.h>
-
-#include "AliRun.h"
-#include "AliESD.h"
-#include "AliStack.h"
-#include "AliHeader.h"
-#include "AliESDtrack.h"
-#include "AliRunLoader.h"
-#include "AliGenEventHeader.h"
-#include "AliGenPythiaEventHeader.h"
-
-#include "AliRsnDaughter.h"
-#include "AliRsnEvent.h"
-#include "AliRsnReader.h"
-
-ClassImp(AliRsnReader)
-
-//--------------------------------------------------------------------------------------------------------
-AliRsnReader::AliRsnReader()
-//
-// Constructor.
-// Initializes all working parameters to default values:
-// - ESD particle identification
-// - rejection of non-ITS-refitted tracks
-// - maximum distance allowed from primary vertex = 3.0 cm (beam pipe)
-//
-{
- fPIDMethod = kESDPID;
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) fPrior[i] = 1.0;
- fPtLimit4PID = 4.0;
- fProbThreshold = 0.0;
- fMaxRadius = 3.0; // the beam pipe
-
- fEvents = 0;
- for (i = AliPID::kElectron; i <= AliPID::kProton; i++) {
- fEffPos[i] = 0;
- fEffNeg[i] = 0;
- }
- fResPt = 0;
- fResLambda = 0;
- fResPhi = 0;
-}
-//--------------------------------------------------------------------------------------------------------
-AliRsnReader::AliRsnReader(const AliRsnReader ©) : TObject(copy)
-//
-// Copy constructor.
-// Initializes all working parameters to same values of another simlar object.
-// TTree data member is not created.
-//
-{
- fPIDMethod = copy.fPIDMethod;
- Int_t i;
- for (i = 0; i < AliPID::kSPECIES; i++) fPrior[i] = copy.fPrior[i];
- fPtLimit4PID = copy.fPtLimit4PID;
- fProbThreshold = copy.fProbThreshold;
- fMaxRadius = copy.fMaxRadius;
-
- fEvents = 0;
- for (i = AliPID::kElectron; i <= AliPID::kProton; i++) {
- fEffPos[i] = copy.fEffPos[i];
- fEffNeg[i] = copy.fEffNeg[i];
- }
-
- fResPt = copy.fResPt;
- fResLambda = copy.fResLambda;
- fResPhi = copy.fResPhi;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnReader::Clear(Option_t *option)
-//
-// Clear collection of filenames.
-// If requested with the option "DELETELIST",
-// the collection object is also deleted.
-//
-{
- TString opt(option);
-
- if (!opt.CompareTo("TREE", TString::kIgnoreCase)) {
- fEvents->Reset();
- if (!opt.CompareTo("DELTREE", TString::kIgnoreCase)) {
- delete fEvents;
- fEvents = 0;
- }
- }
-
- Int_t i;
- for (i = AliPID::kElectron; i <= AliPID::kProton; i++) {
- fEffPos[i] = 0;
- fEffNeg[i] = 0;
- }
- fResPt = 0;
- fResLambda = 0;
- fResPhi = 0;
-}
-//--------------------------------------------------------------------------------------------------------
-Bool_t AliRsnReader::EffSim(Int_t pdg, Double_t pt, Double_t eta, Double_t z)
-//
-// If efficiency histogram are present, they are used to simulate efficiency.
-// Pt, Eta and Z are used to find reconstruction efficiency value to be used
-// and PDG is used to select efficiency histogram for a given particle.
-// An extraction is done, and it is supposed that particle must be accepted
-// only when this function retunrs kTRUE (= successful extraction).
-//
-{
- // find particle sign from PDG code
- Int_t sign;
- if (TMath::Abs(pdg) >= 211) {
- if (pdg > 0) sign = 1; else sign = -1;
- }
- else {
- if (pdg > 0) sign = -1; else sign = 1;
- }
-
- // convert PDG code into one value in AliPID::kSPECIES
- // (if returned value is not a charged particle, procedure is skipped)
- Int_t index = FindType(pdg);
- TH3D *eff = 0;
- if (index >= AliPID::kElectron && index <= AliPID::kProton) {
- if (sign > 0) eff = fEffPos[index]; else eff = fEffNeg[index];
- }
-
- // if no efficiency histogram is defined, method returns a 'fail' result
- if (!eff) return kFALSE;
-
- // otherwise, a random extraction is made
- Int_t ibin = eff->FindBin(z, eta, pt);
- Double_t ref = (Double_t)eff->GetBinContent(ibin);
- Double_t ran = gRandom->Rndm();
- return (ran <= ref);
-}
-//--------------------------------------------------------------------------------------------------------
-Double_t* AliRsnReader::GetPIDprobabilities(AliRsnDaughter track)
-//
-// Computes PID probabilites starting from priors and weights
-//
-{
- Double_t *prob = new Double_t[AliPID::kSPECIES];
-
- Int_t i;
-
- // step 1 - compute the normalization factor
- Double_t sum = 0.0;
- for (i = 0; i < AliPID::kSPECIES; i++) {
- prob[i] = fPrior[i] * track.GetPIDweight(i);
- sum += prob[i];
- }
- if (sum <= 0.0) return 0;
-
- // step 2 - normalize PID weights by the relative prior probability
- for (Int_t i = 0; i < AliPID::kSPECIES; i++) {
- prob[i] /= sum;
- }
-
- return prob;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnReader::Identify(AliRsnDaughter &track)
-//
-// Identifies a track according to method selected
-//
-{
- Bool_t doESDPID = (fPIDMethod == kESDPID);
- Bool_t keepRecSign = (fPIDMethod == kPerfectPIDwithRecSign);
-
- if (doESDPID) {
- // when doing ESD PID it is supposed that charge sign
- // comes from ESD track and it is not modified
- Double_t pt = track.GetPt();
- if (pt <= fPtLimit4PID) {
- Double_t *prob = GetPIDprobabilities(track);
- if (!prob) track.SetPDG(0);
- Int_t idx[AliPID::kSPECIES];
- TMath::Sort(AliPID::kSPECIES, prob, idx);
- Int_t imax = idx[0];
- Double_t maxprob = prob[imax];
- if (maxprob >= fProbThreshold) {
- track.SetPDG((UShort_t)AliPID::ParticleCode(imax));
- }
- delete [] prob;
- }
- else {
- track.SetPDG(0);
- }
- }
- else {
- Short_t truePDG = track.GetTruePDG();
- track.SetPDG((UShort_t)TMath::Abs(truePDG));
- if (!keepRecSign) {
- if (TMath::Abs(truePDG) <= 13) {
- if (truePDG > 0) track.SetSign(-1); else track.SetSign(1);
- }
- else {
- if (truePDG > 0) track.SetSign(1); else track.SetSign(-1);
- }
- }
- }
-}
-//--------------------------------------------------------------------------------------------------------
-TTree* AliRsnReader::ReadParticles(const char *path, Option_t *option)
-//
-// Opens the file "galice.root" and kinematics in the path specified,
-// loads Kinematics informations and fills and AliRsnEvent object
-// with data coming from simulation.
-// If required, an efficiency simulation can be done which cuts off some tracks
-// depending on a MonteCarlo extraction weighted on the efficiency histograms
-// passed to the class.
-// Moreover, a smearing can be done if requested, using the resolution histograms
-// passed to the class.
-// Allowed options:
-// - "E" --> do efficiency simulation
-// - "P" --> do momentum smearing
-//
-{
- fPIDMethod = kPerfectPID;
-
- TTree *events = new TTree("selection", "AliRsnEvents");
- TTree::SetBranchStyle(1);
- AliRsnEvent *event = new AliRsnEvent;
- TBranch *branch = events->Branch("events", "AliRsnEvent", &event, 32000, 1);
- branch->SetAutoDelete(kFALSE);
-
- // get path
- TString strPath(path);
- if (strPath.Length() < 1) return 0;
- strPath += '/';
-
- // evaluate options
- TString opt(option);
- opt.ToUpper();
- Bool_t doEffSim = (opt.Contains("E"));
- Bool_t doSmearing = (opt.Contains("P"));
-
- // initialize run loader
- AliRunLoader *runLoader = OpenRunLoader(path);
- if (!runLoader) return 0;
- Int_t nEvents = runLoader->GetNumberOfEvents();
- TArrayF vertex(3);
-
- // loop on events
- Int_t i, iev, index = 0, imum, nParticles, nStoredParticles = 0;
- Double_t vtot, px, py, pz, arg, eta;
- TParticle *particle = 0, *mum = 0;
- for (iev = 0; iev < nEvents; iev++) {
-
- // load given event
- cout << "\rEvent " << iev << " of " << nEvents << flush;
- runLoader->GetEvent(iev);
-
- // primary vertex
- runLoader->GetHeader()->GenEventHeader()->PrimaryVertex(vertex);
- //cout << "Process type: " << ((AliGenPythiaEventHeader*)runLoader->GetHeader()->GenEventHeader())->ProcessType() << endl;
-
- // kinematics
- AliStack *stack = runLoader->Stack();
- if (!stack) continue;
-
- // get number of particles to collect
- nParticles = stack->GetNtrack();
- if (!nParticles) return 0;
-
- // create new AliRsnEvent object
- event->Clear("DELETE");
- event->SetESD(kFALSE);
- event->SetPath(strPath);
- event->Init();
-
- // loop on tracks
- index = 0;
- for (i = 0; i < nParticles; i++) {
-
- // get particle
- particle = stack->Particle(i);
- if (!particle) continue;
-
- // check against maximum allowed distance from primary vertex
- vtot = (particle->Vx() - vertex[0])*(particle->Vx() - vertex[0]);
- vtot += (particle->Vy() - vertex[1])*(particle->Vy() - vertex[1]);
- vtot += (particle->Vz() - vertex[2])*(particle->Vz() - vertex[2]);
- vtot = TMath::Sqrt(vtot);
- if (vtot > fMaxRadius) continue;
-
- // efficiency selection
- if (doEffSim) {
- arg = 0.5 * TMath::ATan2(particle->Pz(), particle->Pt());
- if (arg > 0.) eta = -TMath::Log(arg); else continue;
- Bool_t result = EffSim(particle->GetPdgCode(), (Double_t)vertex[2], eta, (Double_t)particle->Pt());
- if (result == kFALSE) continue;
- }
-
- // smear kinematic parameters (if requested)
- px = particle->Px();
- py = particle->Py();
- pz = particle->Pz();
- if (doSmearing) SmearMomentum(px, py, pz);
-
- // add to collection of tracks
- nStoredParticles++;
- AliRsnDaughter *track = new AliRsnDaughter;
- if (!track->Adopt(particle)) {
- delete track;
- continue;
- }
- track->SetIndex(index);
- track->SetLabel(i);
- track->SetPxPyPz(px, py, pz);
- imum = particle->GetFirstMother();
- if (imum >= 0) {
- mum = stack->Particle(imum);
- track->SetMotherPDG( (Short_t)mum->GetPdgCode() );
- }
- Identify(*track);
-
- event->AddTrack(*track);
- delete track;
-
- } // end of loop over tracks
-
- // compute total multiplicity
- event->GetMultiplicity();
-
- // link to events tree and fill
- events->Fill();
- }
-
- runLoader->UnloadKinematics();
- delete runLoader;
-
- return events;
-}
-//--------------------------------------------------------------------------------------------------------
-TTree* AliRsnReader::ReadTracks(const char *path, Option_t *option)
-//
-// Reads AliESD in a given path, with and "experimental" method.
-// When using this method, no Kinematics information is assumed
-// and particle identification is always done with the bayesian method.
-// No Kinematics informations are stored.
-// Allowed options (actually):
-// - "R" : reject tracks which do not have the flag "kITSRefit" set to TRUE
-// - "F" : reject 'fake' tracks (negative label)
-//
-{
- fPIDMethod = kESDPID;
-
- TTree *events = new TTree("selection", "AliRsnEvents");
- TTree::SetBranchStyle(1);
- AliRsnEvent *event = new AliRsnEvent;
- TBranch *branch = events->Branch("events", "AliRsnEvent", &event, 32000, 1);
- branch->SetAutoDelete(kFALSE);
-
- // get path
- TString strPath(path);
- if (strPath.Length() < 1) return 0;
- strPath += '/';
-
- // evaluate options
- TString opt(option);
- opt.ToUpper();
- Bool_t checkITSrefit = (opt.Contains("R"));
- Bool_t rejectFakes = (opt.Contains("F"));
-
- // opens ESD file
- TFile *fileESD = TFile::Open(Form("%s/AliESDs.root", strPath.Data()));
- if (!fileESD) return 0;
- if (fileESD->IsOpen() == kFALSE) return 0;
- TTree* treeESD = (TTree*)fileESD->Get("esdTree");
- AliESD *esd = 0;
- treeESD->SetBranchAddress("ESD", &esd);
- Int_t nev = (Int_t)treeESD->GetEntries();
-
- // loop on events
- Int_t i, nSelTracks = 0;
- Double_t vertex[3];
- for (i = 0; i < nev; i++) {
-
- // message
- cout << "\rEvent " << i << flush;
- treeESD->GetEntry(i);
-
- // primary vertex
- vertex[0] = esd->GetVertex()->GetXv();
- vertex[1] = esd->GetVertex()->GetYv();
- vertex[2] = esd->GetVertex()->GetZv();
-
- // create new AliRsnEvent object
- event->Clear("DELETE");
- event->Init();
- event->SetPath(strPath);
- event->SetESD();
-
- // get number of tracks
- Int_t ntracks = esd->GetNumberOfTracks();
- if (!ntracks) continue;
-
- // store tracks from ESD
- Int_t index, label;
- Double_t vtot, v[3];
- for (index = 0; index < ntracks; index++) {
-
- // get track
- AliESDtrack *esdTrack = esd->GetTrack(index);
-
- // check against vertex constraint
- esdTrack->GetXYZ(v);
- vtot = (v[0] - vertex[0])*(v[0] - vertex[0]);
- vtot += (v[1] - vertex[1])*(v[1] - vertex[1]);
- vtot += (v[2] - vertex[2])*(v[2] - vertex[2]);
- vtot = TMath::Sqrt(vtot);
- if (vtot > fMaxRadius) continue;
-
- // check for ITS refit
- if (checkITSrefit) {
- if (!(esdTrack->GetStatus() & AliESDtrack::kITSrefit)) continue;
- }
-
- // check for fakes
- label = esdTrack->GetLabel();
- if (rejectFakes && (label < 0)) continue;
-
- // create AliRsnDaughter (and make Bayesian PID)
- AliRsnDaughter track;
- if (!track.Adopt(esdTrack, checkITSrefit)) continue;
- track.SetIndex(index);
- track.SetLabel(label);
- Identify(track);
-
- // store in TClonesArray
- event->AddTrack(track);
- nSelTracks++;
- }
-
- // compute total multiplicity
- event->GetMultiplicity();
-
- // link to events tree and fill
- events->Fill();
- }
-
- fileESD->Close();
- return events;
-}
-//--------------------------------------------------------------------------------------------------------
-TTree* AliRsnReader::ReadTracksAndParticles(const char *path, Option_t *option)
-//
-// Reads AliESD in a given path, getting also informations from Kinematics.
-// In this case, the PID method used is the one selected with apposite setter.
-// Allowed options (actually):
-// - "R" : reject tracks which do not have the flag "kITSRefit" set to TRUE
-// - "F" : reject 'fake' tracks (negative label)
-// - "M" : use 'true' momentum instead of reconstructed one
-//
-{
- TTree *events = new TTree("selection", "AliRsnEvents");
- TTree::SetBranchStyle(1);
- AliRsnEvent *event = new AliRsnEvent;
- TBranch *branch = events->Branch("events", "AliRsnEvent", &event, 32000, 1);
- branch->SetAutoDelete(kFALSE);
-
- // get path
- TString strPath(path);
- if (strPath.Length() < 1) return 0;
-
- // evaluate options
- TString opt(option);
- opt.ToUpper();
- Bool_t checkITSrefit = (opt.Contains("R"));
- Bool_t rejectFakes = (opt.Contains("F"));
- Bool_t copyMomentum = (opt.Contains("M"));
-
- // opens ESD file
- TFile *fileESD = TFile::Open(Form("%s/AliESDs.root", strPath.Data()));
- if (!fileESD) return 0;
- if (fileESD->IsOpen() == kFALSE) return 0;
- TTree* treeESD = (TTree*)fileESD->Get("esdTree");
- AliESD *esd = 0;
- treeESD->SetBranchAddress("ESD", &esd);
- Int_t nevRec = (Int_t)treeESD->GetEntries();
-
- // initialize run loader
- AliRunLoader *runLoader = OpenRunLoader(path);
- if (!runLoader) return 0;
- Int_t nevSim = runLoader->GetNumberOfEvents();
-
- // check number of reconstructed and simulated events
- if ( (nevSim != 0 && nevRec != 0) && (nevSim != nevRec) ) {
- cerr << "Count mismatch: sim = " << nevSim << " -- rec = " << nevRec << endl;
- return 0;
- }
- else if (nevSim == 0 && nevRec == 0) {
- cerr << "Count error: sim = rec = 0" << endl;
- return 0;
- }
-
- // loop on events
- Int_t i, procType, ntracks, nSelTracks = 0;
- Double_t vertex[3];
- for (i = 0; i < nevRec; i++) {
-
- // get event
- cout << "\rEvent " << i << " " << flush;
- treeESD->GetEntry(i);
- runLoader->GetEvent(i);
-
- // reject event if it is diffractive
- procType = ((AliGenPythiaEventHeader*)runLoader->GetHeader()->GenEventHeader())->ProcessType();
- if (procType == 92 || procType == 93 || procType == 94) {
- cout << "Skipping diffractive event" << endl;
- continue;
- }
-
- // get particle stack
- AliStack *stack = runLoader->Stack();
-
- // primary vertex
- vertex[0] = esd->GetVertex()->GetXv();
- vertex[1] = esd->GetVertex()->GetYv();
- vertex[2] = esd->GetVertex()->GetZv();
-
- // multiplicity
- ntracks = esd->GetNumberOfTracks();
- if (!ntracks) {
- Warning("ReadTracksAndParticles", "Event %d has no tracks!", i);
- continue;
- }
-
- // create new AliRsnEvent object
- event->Clear("DELETE");
- event->Init();
- event->SetPath(strPath);
- event->SetESD();
-
- // store tracks from ESD
- Int_t index, label;
- Double_t vtot, v[3];
- for (index = 0; index < ntracks; index++) {
-
- // get track
- AliESDtrack *esdTrack = esd->GetTrack(index);
-
- // check against vertex constraint
- esdTrack->GetXYZ(v);
- vtot = (v[0] - vertex[0])*(v[0] - vertex[0]);
- vtot += (v[1] - vertex[1])*(v[1] - vertex[1]);
- vtot += (v[2] - vertex[2])*(v[2] - vertex[2]);
- vtot = TMath::Sqrt(vtot);
- if (vtot > fMaxRadius) continue;
-
- // check for fakes
- label = esdTrack->GetLabel();
- if (rejectFakes) {
- if (label < 0) continue;
- }
-
- // create AliRsnDaughter (and make Bayesian PID)
- AliRsnDaughter track;
- if (!track.Adopt(esdTrack, checkITSrefit)) continue;
- track.SetIndex(index);
-
- // retrieve particle and get Kine info
- TParticle *part = stack->Particle(TMath::Abs(label));
- track.SetTruePDG(part->GetPdgCode());
- Int_t mother = part->GetFirstMother();
- track.SetMother(mother);
- if (mother >= 0) {
- TParticle *mum = stack->Particle(mother);
- track.SetMotherPDG(mum->GetPdgCode());
- }
- if (copyMomentum) track.SetPxPyPz(part->Px(), part->Py(), part->Pz());
-
- // identification
- Identify(track);
-
- // store in TClonesArray
- event->AddTrack(track);
- nSelTracks++;
- }
-
- // compute total multiplicity
- event->GetMultiplicity();
-
- // link to events tree and fill
- events->Fill();
- }
-
- runLoader->UnloadKinematics();
- delete runLoader;
- fileESD->Close();
-
- return events;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnReader::SetEfficiencyHistogram(AliPID::EParticleType type, TH3D *h, Bool_t pos)
-//
-// Sets efficiency histogram for a given particle species.
-// If third argument is "true", histo is assigned to positive particles,
-// otherwise it is assigned to negative particles.
-//
-{
- if (type >= AliPID::kElectron && type < AliPID::kPhoton) {
- if (pos) fEffPos[type] = h; else fEffNeg[type] = h;
- }
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnReader::SetPriorProbabilities(Double_t *prior)
-//
-// Set prior probabilities to be used in case of ESD PID.
-//
-{
- if (!prior) return;
-
- Int_t i = 0;
- for (i = 0; i < AliPID::kSPECIES; i++) fPrior[i] = prior[i];
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnReader::SetPriorProbability(AliPID::EParticleType type, Double_t value)
-//
-// Sets prior probability referred to a single particle species.
-//
-{
- if (type >= AliPID::kElectron && type < AliPID::kPhoton) fPrior[type] = value;
-}
-//--------------------------------------------------------------------------------------------------------
-void AliRsnReader::SmearMomentum(Double_t &px, Double_t &py, Double_t &pz)
-//
-// Use resolution histograms to do smearing of momentum
-// (to introduce reconstruction effects)
-//
-{
- Double_t pt = TMath::Sqrt(px*px + py*py);
- Double_t lambda = TMath::ATan2(pz, pt);
- Double_t phi = TMath::ATan2(py, px);
-
- Int_t ibin;
- Double_t ref;
- if (fResPt) {
- ibin = fResPt->FindBin(pt);
- ref = (Double_t)fResPt->GetBinContent(ibin);
- pt *= 1.0 + ref;
- }
- if (fResLambda) {
- ibin = fResLambda->FindBin(pt);
- ref = (Double_t)fResLambda->GetBinContent(ibin);
- lambda *= 1.0 + ref;
- }
- if (fResPhi) {
- ibin = fResPhi->FindBin(pt);
- ref = (Double_t)fResPhi->GetBinContent(ibin);
- phi *= 1.0 + ref;
- }
-
- px = pt * TMath::Cos(phi);
- py = pt * TMath::Sin(phi);
- pz = pt * TMath::Tan(lambda);
-}
-//--------------------------------------------------------------------------------------------------------
-AliPID::EParticleType AliRsnReader::FindType(Int_t pdg)
-//
-// Finds the enum value corresponding to a PDG code
-//
-{
- pdg = TMath::Abs(pdg);
- switch (pdg) {
- case 11: return AliPID::kElectron; break;
- case 13: return AliPID::kMuon; break;
- case 211: return AliPID::kPion; break;
- case 321: return AliPID::kKaon; break;
- case 2212: return AliPID::kProton; break;
- default : return AliPID::kPhoton;
- }
-}
-//--------------------------------------------------------------------------------------------------------
-AliRunLoader* AliRsnReader::OpenRunLoader(const char *path)
-//
-// Open the Run loader with events in a given path
-//
-{
- // clear gALICE
- if (gAlice) {
- delete gAlice;
- gAlice = 0;
- }
-
- // initialize run loader
- TString name(path);
- name += "/galice.root";
- AliRunLoader *runLoader = AliRunLoader::Open(name.Data());
- if (runLoader) {
- runLoader->LoadgAlice();
- gAlice = runLoader->GetAliRun();
- runLoader->LoadKinematics();
- runLoader->LoadHeader();
- }
-
- return runLoader;
-}
+++ /dev/null
-/**************************************************************************
- * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
- * See cxx source for full Copyright notice *
- **************************************************************************/
-
-//-------------------------------------------------------------------------
-// Class AliRsnEventReader
-//
-// Reader for conversion of ESD or Kinematics output into AliRsnEvent
-//
-// author: A. Pulvirenti (email: alberto.pulvirenti@ct.infn.it)
-//-------------------------------------------------------------------------
-
-#ifndef AliRSNREADER_H
-#define AliRSNREADER_H
-
-#include "AliPID.h"
-
-class TH3D;
-class TH1D;
-class TOrdCollection;
-class TTree;
-class AliRunLoader;
-
-class AliRsnReader : public TObject
-{
-public:
-
- enum EPIDMethod {
- kNoPID = 0, // no PID performed
- kPerfectPID = 1, // use Kinematics to simulate a 100% PID efficiency
- kESDPID = 2, // use experimental ESD weights
- kPerfectPIDwithRecSign = 3 // get particle type from Kine and charge sign from reconstruction
- };
-
- AliRsnReader();
- AliRsnReader(const AliRsnReader& copy);
- virtual ~AliRsnReader() {Clear("DELTREE");}
-
- void Clear(Option_t *option = "");
- Bool_t EffSim(Int_t pdg, Double_t pt, Double_t eta, Double_t z);
- TTree* GetEvents() const {return fEvents;}
- Double_t* GetPIDprobabilities(AliRsnDaughter track);
- void Identify(AliRsnDaughter &track);
- TTree* ReadParticles(const char *path, Option_t *option="");
- TTree* ReadTracks(const char *path, Option_t *option="R");
- TTree* ReadTracksAndParticles(const char *path, Option_t *option="R");
- void SetEfficiencyHistogram(AliPID::EParticleType type, TH3D *h, Bool_t pos=kTRUE);
- void SetMaxRadius(Double_t value) {fMaxRadius=value;}
- void SetPIDMethod(AliRsnReader::EPIDMethod pm) {fPIDMethod=pm;}
- void SetPriorProbabilities(Double_t *prior);
- void SetPriorProbability(AliPID::EParticleType type, Double_t p);
- void SetProbabilityThreshold(Double_t p) {fProbThreshold=p;}
- void SetPtLimit4PID(Double_t value) {fPtLimit4PID=value;}
- void SetResolutionPt(TH1D *histo) {fResPt=histo;}
- void SetResolutionLambda(TH1D *histo) {fResLambda=histo;}
- void SetResolutionPhi(TH1D *histo) {fResPhi=histo;}
- void SmearMomentum(Double_t &px, Double_t &py, Double_t &pz);
-
-protected:
-
- AliPID::EParticleType FindType(Int_t pdg);
- AliRunLoader* OpenRunLoader(const char *path);
-
- Bool_t fDoEffSim; // set to true to introduce efficiency effects
- Bool_t fDoSmearing; // when selecting only particles, set this to true
- // to smear particle momentum
-
- EPIDMethod fPIDMethod; // PID method
- Double_t fPrior[AliPID::kSPECIES]; // prior probabilities (in case of REAL pid)
- Double_t fPtLimit4PID; // maximum transverse momentum to accept realistic PID
- Double_t fProbThreshold; // minimum required probability to accept realistic PID
- Double_t fMaxRadius; // maximum allowed distance from primary vertex
-
- TTree *fEvents; //! tree of read events
- TH3D *fEffPos[AliPID::kSPECIES]; // tracking efficiency (in pt, eta and Zv) for positive particles
- TH3D *fEffNeg[AliPID::kSPECIES]; // tracking efficiency (in pt, eta and Zv) for negative particles
- TH1D *fResPt; // relative Pt resolution
- TH1D *fResLambda; // relative Lambda resolution
- TH1D *fResPhi; // relative Phi resolution
-
- // Rsn event reader implementation
- ClassDef(AliRsnReader,1)
-};
-
-#endif
-