// Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
//
+#include <Riostream.h>
+
#include <TSystem.h>
#include <TFile.h>
#include <TFolder.h>
ClassImp(AliRsnAnalysisSE)
//________________________________________________________________________
-AliRsnAnalysisSE::AliRsnAnalysisSE(const char * name)
- : AliRsnAnalysisTaskSEBase(name),
- fPairMgrs(0),fOutList(0x0),fRsnEventBuffer(0x0),
- fNumOfEventsInBuffer(1000)
+AliRsnAnalysisSE::AliRsnAnalysisSE(const char * name, Int_t bufferSize) :
+ AliRsnAnalysisTaskSEBase(name),
+ fDoesMixing(kFALSE),
+ fMixingNum(0),
+ fMixingCut(0x0),
+ fPairMgrs(0),
+ fOutList(0x0),
+ fBuffer(0x0),
+ fBufferSize(bufferSize)
{
-//=========================================================
+//
// Default constructor
-//=========================================================
+//
InitIOVars();
DefineOutput(1, TList::Class());
//________________________________________________________________________
AliRsnAnalysisSE::~AliRsnAnalysisSE()
{
-//=========================================================
+//
// Destructor
-//=========================================================
+//
}
//________________________________________________________________________
void AliRsnAnalysisSE::InitIOVars()
{
-//=========================================================
+//
// Init input output values
-//=========================================================
+//
- AliDebug(AliLog::kDebug, "<-");
AliRsnAnalysisTaskSEBase::InitIOVars();
- fRsnEventBuffer = 0;
+ fBuffer = 0;
fOutList = 0;
-
- AliDebug(AliLog::kDebug, "->");
}
//________________________________________________________________________
void AliRsnAnalysisSE::UserCreateOutputObjects()
{
-//=========================================================
+//
// UserCreateOutputObjects() of AliAnalysisTaskSE
-//=========================================================
+//
- AliDebug(AliLog::kDebug, "<-");
-// fPID.DumpPriors();
OpenFile(0);
fOutList = new TList();
fOutList->SetOwner();
AliRsnPair *def=0;
AliRsnPairMgr *mgr=0;
TList *listTmp;
+ fDoesMixing = kFALSE;
for (Int_t iMgr=0 ;iMgr< fPairMgrs.GetEntries();iMgr++)
{
mgr = (AliRsnPairMgr *) fPairMgrs.At(iMgr);
def = (AliRsnPair *) mgr->GetPairs()->At(i);
if (def)
{
- def->Init();
listTmp->Add(def->GenerateHistograms(mgr->GetName()));
- //def->GenerateHistograms(mgr->GetName(), listTmp);
- //def->Print();
+ if (def->IsMixed()) fDoesMixing = kTRUE;
}
}
fOutList->Add(listTmp);
}
- fRsnEventBuffer = new AliRsnEventBuffer(fNumOfEventsInBuffer);
-// fRsnEventBuffer = new AliRsnEventBuffer ( 10000 ,kFALSE );
- AliDebug(AliLog::kDebug, "->");
-
+ fBuffer = new AliRsnEventBuffer(fBufferSize);
}
//________________________________________________________________________
void AliRsnAnalysisSE::UserExec(Option_t *)
{
-//=========================================================
+//
// UserExec() of AliAnalysisTaskSE
-//=========================================================
+//
+
+ static UInt_t eventID = 0;
- if (fEntry++%1000==0)
- AliInfo(Form("Event %d",fEntry-1));
+ if (fEntry++ % 100 == 0) cout << "[" << GetName() << "] : processing entry " << fEntry-1 << endl;
AliRsnEvent *curEvent = GetRsnEventFromInputType();
if (!curEvent) return;
+ curEvent->SetUniqueID(eventID++);
ProcessEventAnalysis(curEvent);
PostEventProcess();
//________________________________________________________________________
void AliRsnAnalysisSE::Terminate(Option_t *)
{
-//=========================================================
+//
// Terminate() of AliAnalysisTask
-//=========================================================
+//
fOutList = dynamic_cast<TList*>(GetOutputData(1));
- if (!fOutList) { AliError(" fOutList not available"); return; }
- fOutList->Print();
+ if (!fOutList) { AliError("At end of analysis, output list is NULL"); return; }
+ //fOutList->Print();
}
//________________________________________________________________________
void AliRsnAnalysisSE::ProcessEventAnalysis(AliRsnEvent *curEvent)
{
-//=========================================================
+//
// Process of one event
-//=========================================================
-
+//
// Adds event to Event Buffer
- fRsnEventBuffer->AddEvent(curEvent);
+ fBuffer->AddEvent(curEvent);
+ Int_t index = fBuffer->GetEventsBufferIndex();
+
+ Int_t nmatches;
+ TArrayI matched(0);
+ if (fDoesMixing) matched = FindGoodMatches(index, nmatches);
// loop over all Pair managers
AliRsnPairMgr *mgr=0;
for (Int_t i=0;i< mgr->GetPairs()->GetEntriesFast();i++)
{
pair = (AliRsnPair *) mgr->GetPairs()->At(i);
- pair->ProcessPair(fRsnEventBuffer);
+ if (!pair->IsMixed()) {
+ pair->ProcessPair(curEvent, 0);
+ }
+ else {
+ Int_t i, iev;
+ for (i = 0; i < matched.GetSize(); i++) {
+ iev = matched[i];
+ if (iev < 0) continue;
+ AliRsnEvent *evmatch = fBuffer->GetEvent(iev);
+ pair->ProcessPair(curEvent, evmatch);
+ if (!pair->IsPairEqual()) pair->ProcessPair(evmatch, curEvent);
+ }
+ }
}
}
}
//________________________________________________________________________
void AliRsnAnalysisSE::PostEventProcess(const Short_t & index)
{
-//=========================================================
+//
// Post process of one event
-//=========================================================
+//
- switch (fInputType[index])
- {
- case kAOD:
- break;
- case kESD:
- break;
- case kESDMC:
- break;
- case kMC:
- break;
- case kRSN:
- {
- if (fRsnEventBuffer->GetDeleteBufferWhenReset() == kFALSE)
- {
- fRSN[index] = (AliRsnEvent*) fRsnEventBuffer->GetNextEvent();
- SetBranchAddress(0 , "rsnEvents", &fRSN[index]);
- }
- break;
- }
- default:
- break;
- }
+ if (fInputType[index] != kRSN) return;
+ if (fBuffer->GetDeleteBufferWhenReset() == kFALSE) {
+ fRSN[index] = (AliRsnEvent*) fBuffer->GetNextEvent();
+ SetBranchAddress(0 , "rsnEvents", &fRSN[index]);
+ }
}
+//________________________________________________________________________
void AliRsnAnalysisSE::AddPairMgr(AliRsnPairMgr * pairmgr)
{
fPairMgrs.Add(pairmgr);
}
-
+//________________________________________________________________________
void AliRsnAnalysisSE::AddPairMgrFromConfig(TString configfile)
{
gROOT->LoadMacro(configfile.Data());
configfile.ReplaceAll(".C","");
- AliRsnPairMgr *mgrRsn = (AliRsnPairMgr *) gROOT->ProcessLine(Form("%s();",configfile.Data()));
+ AliRsnPairMgr *mgrRsn = (AliRsnPairMgr *) gROOT->ProcessLine(Form("%s();", configfile.Data()));
if (!mgrRsn) return;
fPairMgrs.Add(mgrRsn);
}
+
+//________________________________________________________________________
+TArrayI AliRsnAnalysisSE::FindGoodMatches(Int_t iRef, Int_t &foundMatches)
+{
+ // initialize the output array to the size of required mixed events
+ // and initialize all members to -1
+ Int_t i;
+ TArrayI matched(fMixingNum);
+ for (i = 0; i < fMixingNum; i++) matched[i] = -1;
+ foundMatches = 0;
+
+ // starts from the position behind the reference index
+ // and goes backward; if it reaches the value 0, stops
+ AliRsnEvent *refEvent = fBuffer->GetEvent(iRef);
+ if (!refEvent) return matched;
+ AliRsnEvent *matchEvent = 0x0;
+ Int_t checkIndex;
+ for (checkIndex = iRef - 1; ; checkIndex--) {
+ if (checkIndex < 0) checkIndex = fBuffer->GetEventsBufferSize() - 1;
+ if (checkIndex == iRef) break;
+ matchEvent = fBuffer->GetEvent(checkIndex);
+ if (!matchEvent) continue;
+ if (fMixingCut) {
+ if (!fMixingCut->IsSelected(AliRsnCut::kMixEvent, refEvent, matchEvent)) continue;
+ }
+ // assign to current array slot the matched event
+ // and increment current slot and stops if it exceeds array size
+ matched[foundMatches++] = checkIndex;
+ if (foundMatches >= fMixingNum) break;
+ }
+
+ // returns the current index value,
+ // which is also the number of matched events found
+ return matched;
+}
#define ALIRSNANALYSISAT_H
#include <TH1.h>
+#include <TArrayI.h>
#include "AliRsnEventBuffer.h"
#include "AliRsnPair.h"
class AliRsnAnalysisSE : public AliRsnAnalysisTaskSEBase
{
public:
- AliRsnAnalysisSE(const char *name = "AliRsnAnalysisSE");
- AliRsnAnalysisSE(const AliRsnAnalysisSE& copy): AliRsnAnalysisTaskSEBase(copy),
- fPairMgrs(0),fOutList(0x0),fRsnEventBuffer(0x0),fNumOfEventsInBuffer(100) {}
- AliRsnAnalysisSE& operator= (const AliRsnAnalysisSE&) {return *this;}
+ AliRsnAnalysisSE(const char *name = "AliRsnAnalysisSE", Int_t bufferSize = 1000);
~AliRsnAnalysisSE();
virtual void InitIOVars();
-// virtual void LocalInit();
-// virtual Bool_t Notify();
virtual void UserCreateOutputObjects();
virtual void UserExec(Option_t *option);
virtual void Terminate(Option_t *);
void AddPairMgr(AliRsnPairMgr *pairmgr);
void AddPairMgrFromConfig(TString configfile);
- void SetNumOfEventsInBuffer(const Int_t& theValue) { fNumOfEventsInBuffer = theValue; }
- Int_t GetNumOfEventsInBuffer() const { return fNumOfEventsInBuffer; }
+ void SetBufferSize(const Int_t size) {fBufferSize = size;}
+ Int_t GetBufferSize() const {return fBufferSize;}
+
+ TArrayI FindGoodMatches(Int_t iRef, Int_t &found);
+ void SetMixingNum(Int_t i) {fMixingNum = i;}
+ void SetMixingCut(AliRsnCutSet *cut) {fMixingCut = cut;}
private:
- TObjArray fPairMgrs;
+ AliRsnAnalysisSE(const AliRsnAnalysisSE& copy) :
+ AliRsnAnalysisTaskSEBase(copy),
+ fDoesMixing(kFALSE),fMixingNum(0),fMixingCut(0x0),fPairMgrs(0),
+ fOutList(0x0),fBuffer(0x0),fBufferSize(0) {}
+ AliRsnAnalysisSE& operator= (const AliRsnAnalysisSE&) {return *this;}
- TList *fOutList; // List of output
- AliRsnEventBuffer *fRsnEventBuffer; // event buffer
- Int_t fNumOfEventsInBuffer; // number of events in buffer
+ void ProcessEventAnalysis(AliRsnEvent *curEvent);
+ void PostEventProcess(const Short_t &index = 0);
- void ProcessEventAnalysis(AliRsnEvent *curEvent);
- void PostEventProcess(const Short_t &index=0);
+ Bool_t fDoesMixing; // flag set to kTRUE if the task contains pairs for mixing
+ Int_t fMixingNum; // number of events to mix
+ AliRsnCutSet *fMixingCut; // mixing cut
+ TObjArray fPairMgrs; // collections of pairs
+ TList *fOutList; // List of output
+ AliRsnEventBuffer *fBuffer; // event buffer
+ Int_t fBufferSize; // number of events in buffer
ClassDef(AliRsnAnalysisSE, 1)
};
#include "AliRsnEvent.h"
#include "AliMCEvent.h"
+#include "AliRsnMCInfo.h"
+#include "AliRsnDaughter.h"
#include "AliRsnAnalysisTaskSEBase.h"
return IsBetween(daughter->Eta());
case kRadialImpactParam:
return IsBetween(daughter->Dr());
+ case kRadialImpactParamMC:
+ if (mcinfo) return IsBetween(mcinfo->Dr());
+ else return kTRUE;
case kMomentumMC:
if (mcinfo) return IsBetween(mcinfo->P());
else return kTRUE;
case kTransMomentumMC:
- if (mcinfo) return IsBetween(mcinfo->P());
+ if (mcinfo) return IsBetween(mcinfo->Pt());
else return kTRUE;
case kStatus:
return daughter->CheckFlag(fUIMin);
pidType = daughter->PIDType(prob);
if (fType == kPIDType) return MatchesValue((Int_t) pidType);
if (fType == kPIDProb) return IsBetween(prob);
+ case kPIDProbForSpecies:
+ return IsBetween(daughter->PIDProb()[(AliRsnPID::EType)fIMin]);
case kTruePID:
pdg = TMath::Abs(mcinfo->PDG());
cut = MatchesValue(pdg);
case kIsPrimary:
if (mcinfo) return (mcinfo->Mother() < 0);
else return kTRUE;
+ case kIsKinkMother:
+ return (!daughter->IsKink() || daughter->IsKinkMother());
+ case kIsKinkDaughter:
+ return daughter->IsKinkDaughter();
case kNSigma:
return IsBetween(daughter->NSigmaToVertex());
/*
{
case kMultiplicity:
return IsBetween((Int_t) event->GetMultiplicity());
+ case kVz:
+ return IsBetween((Double_t)event->GetPrimaryVertexZ());
default:
AliWarning("Requested a cut which cannot be applied to an event");
return kTRUE;
kTransMomentum,
kEta,
kRadialImpactParam,
+ kRadialImpactParamMC,
kMomentumMC,
kTransMomentumMC,
kEtaMC,
kIsLabelEqual,
kIsTruePair,
kIsPrimary,
+ kIsKinkDaughter,
+ kIsKinkMother,
+ kMCTracked,
kChargePos,
kChargeNeg,
kPIDType,
kPIDProb,
+ kPIDProbForSpecies,
kTruePID,
+ kVz,
kMultiplicity,
kMultiplicityDifference,
kMultiplicityRatio,
if (copy.fMCInfo) fMCInfo = new AliRsnMCInfo(*(copy.fMCInfo));
}
-//_____________________________________________________________________________
-AliRsnDaughter::AliRsnDaughter(AliESDtrack *track) :
- AliVParticle(),
- fIndex(-1),
- fLabel(-1),
- fCharge(0),
- fFlags(0),
- fKink(0),
- fMass(0.0),
- fChi2(0.0),
- fNSigmaToVertex(-1.0),
- fITSnum(0),
- fTPCnum(0),
- fRealisticPID(AliRsnPID::kUnknown),
- fMCInfo(0x0)
-{
-//
-// Constructor to get data from an ESD track.
-//
-
- Int_t i;
- for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDProb[i] = 0.0;
- AliRsnPIDDefESD pidDef;
- pidDef.UseESDWeights();
- Adopt(track, pidDef);
-}
-
-//_____________________________________________________________________________
-AliRsnDaughter::AliRsnDaughter(AliAODTrack *track) :
- AliVParticle(),
- fIndex(-1),
- fLabel(-1),
- fCharge(0),
- fFlags(0),
- fKink(0),
- fMass(0.0),
- fChi2(0.0),
- fNSigmaToVertex(-1.0),
- fITSnum(0),
- fTPCnum(0),
- fRealisticPID(AliRsnPID::kUnknown),
- fMCInfo(0x0)
-{
-//
-// Constructor to get data from an AOD track.
-//
-
- Int_t i;
- for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDProb[i] = 0.0;
- Adopt(track);
-}
-
-//_____________________________________________________________________________
-AliRsnDaughter::AliRsnDaughter(AliMCParticle *track) :
- AliVParticle(),
- fIndex(-1),
- fLabel(-1),
- fCharge(0),
- fFlags(0),
- fKink(0),
- fMass(0.0),
- fChi2(0.0),
- fNSigmaToVertex(-1.0),
- fITSnum(0),
- fTPCnum(0),
- fRealisticPID(AliRsnPID::kUnknown),
- fMCInfo(0x0)
-{
-//
-// Constructor to get data from an MC track.
-//
-
- Int_t i;
- for (i = 0; i < AliRsnPID::kSpecies; i++) fPIDProb[i] = 0.0;
- Adopt(track);
-}
-
//_____________________________________________________________________________
AliRsnDaughter& AliRsnDaughter::operator=(const AliRsnDaughter ©)
{
}
}
-//_____________________________________________________________________________
-Bool_t AliRsnDaughter::Adopt(AliESDtrack* esdTrack, AliRsnPIDDefESD pidDef)
-{
-//
-// Copies data from an AliESDtrack into "this":
-// - charge sign
-// - momentum
-// - point of closest approach to primary vertex
-// - ESD pid weights
-// In case of errors returns kFALSE, otherwise kTRUE.
-//
-
- if (!esdTrack) {
- AliError("Passed NULL object: nothing done");
- return kFALSE;
- }
-
- // copy values which don't need treatment:
- // momentum, vertex, chi2, flags, charge, number of TPC and ITS clusters
- esdTrack->GetPxPyPz(fP);
- esdTrack->GetXYZ(fV);
- fChi2 = esdTrack->GetConstrainedChi2();
- fFlags = esdTrack->GetStatus();
- fCharge = (Short_t)esdTrack->Charge();
- fITSnum = esdTrack->GetITSclusters(0x0);
- fTPCnum = esdTrack->GetTPCclusters(0x0);
-
- // define the kink index:
- // 0 = no kink
- // 1 = kink daughter
- // -1 = kink mother
- Int_t i, ik[3];
- for (i = 0; i < 3; i++) ik[i] = esdTrack->GetKinkIndex(i);
- if (ik[0] < 0 || ik[1] < 0 || ik[2] < 0) fKink = -1;
- else if (ik[0] > 0 || ik[1] > 0 || ik[2] > 0) fKink = 1;
- else fKink = 0;
-
- // store PID weights according to definition
- pidDef.ComputeWeights(esdTrack, fPIDWeight);
-
- // check:
- // if the sum of all weights is null, the adoption fails
- Double_t sum = 0.0;
- for (i = 0; i < AliRsnPID::kSpecies; i++) sum += fPIDWeight[i];
- if (sum <= 0.0) return kFALSE;
-
- // calculate N sigma to vertex
- AliESDtrackCuts trkCut;
- SetNSigmaToVertex(trkCut.GetSigmaToVertex(esdTrack));
-
- return kTRUE;
-}
-
-
-//_____________________________________________________________________________
-Bool_t AliRsnDaughter::Adopt(AliAODTrack* aodTrack)
-{
-//
-// Copies data from an AliAODtrack into "this":
-// - charge sign
-// - momentum
-// - point of closest approach to primary vertex
-// - ESD pid weights
-// In case of errors returns kFALSE, otherwise kTRUE.
-//
-
- if (!aodTrack)
- {
- AliError("Passed NULL object: nothing done");
- return kFALSE;
- }
-
- // copy momentum and vertex
- fP[0] = aodTrack->Px();
- fP[1] = aodTrack->Py();
- fP[2] = aodTrack->Pz();
- fV[0] = aodTrack->Xv();
- fV[1] = aodTrack->Yv();
- fV[2] = aodTrack->Zv();
-
- // chi2
- fChi2 = aodTrack->Chi2perNDF();
-
- // copy PID weights
- Int_t i;
- for (i = 0; i < 5; i++) fPIDWeight[i] = aodTrack->PID()[i];
-
- // copy flags
- fFlags = aodTrack->GetStatus();
-
- // copy sign
- fCharge = aodTrack->Charge();
-
- return kTRUE;
-}
-
-
-//_____________________________________________________________________________
-Bool_t AliRsnDaughter::Adopt(AliMCParticle *mcParticle)
-{
-//
-// Copies data from a MCParticle into "this":
-// - charge sign
-// - momentum
-// - point of closest approach to primary vertex
-// - ESD pid weights
-// - true PDG code
-// - mother
-// In case of errors returns kFALSE, otherwise kTRUE.
-//
-
- if (!mcParticle)
- {
- AliError("Passed NULL object: nothing done");
- return kFALSE;
- }
-
- // retrieve the TParticle object from the argument
- TParticle *particle = mcParticle->Particle();
- if (!particle)
- {
- AliError("AliMCParticle::Particle() returned NULL");
- return kFALSE;
- }
-
- // copy momentum and vertex
- fP[0] = particle->Px();
- fP[1] = particle->Py();
- fP[2] = particle->Pz();
- fV[0] = particle->Vx();
- fV[1] = particle->Vy();
- fV[2] = particle->Vz();
-
- // recognize charge sign from PDG code sign
- Int_t pdg = particle->GetPdgCode();
- Int_t absPDG = TMath::Abs(pdg);
- if (absPDG == 11 || absPDG == 13)
- {
- if (pdg > 0) fCharge = -1; else fCharge = 1;
- }
- else if (absPDG == 211 || absPDG == 321 || absPDG == 2212)
- {
- if (pdg > 0) fCharge = 1; else fCharge = -1;
- }
- else
- {
- // when trying to "adopt" a neutral track (photon, neutron, etc.)
- // for the moment a "failed" message is returned
- fCharge = 0;
- return kFALSE;
- }
-
- // flags and PID weights make no sense with MC tracks
- fFlags = 0;
- for (pdg = 0; pdg < AliRsnPID::kSpecies; pdg++) fPIDWeight[pdg] = 0.0;
- fPIDWeight[AliRsnPID::InternalType(absPDG)] = 1.0;
-
- // copy other MC info (mother PDG code cannot be retrieved here)
- InitMCInfo(particle);
-
- return kTRUE;
-}
-
//_____________________________________________________________________________
void AliRsnDaughter::Print(Option_t *option) const
{
return kTRUE;
}
-//_____________________________________________________________________________
-Bool_t AliRsnDaughter::InitMCInfo(AliMCParticle *mcParticle)
-{
-//
-// Copies data from an MC particle into the object which
-// contains all MC details taken from kinematics info.
-// If requested by second argument, momentum and vertex
-// of the Particle are copied into the 'fP' and 'fV'
-// data members, to simulate a perfect reconstruction.
-// If something goes wrong, returns kFALSE,
-// otherwise returns kTRUE.
-//
-
- // retrieve the TParticle object pointed by this MC track
- TParticle *particle = mcParticle->Particle();
- return InitMCInfo(particle);
-}
-
//_____________________________________________________________________________
Int_t AliRsnDaughter::Compare(const TObject* obj) const
{
#include "AliVParticle.h"
#include "AliRsnPID.h"
-#include "AliRsnPIDDefESD.h"
-#include "AliESDtrackCuts.h"
class TParticle;
AliRsnDaughter();
AliRsnDaughter(const AliRsnDaughter ©);
- AliRsnDaughter(AliESDtrack *track);
- AliRsnDaughter(AliAODTrack *track);
- AliRsnDaughter(AliMCParticle *track);
virtual ~AliRsnDaughter();
AliRsnDaughter& operator= (const AliRsnDaughter& copy);
void ShiftZero(Double_t x, Double_t y, Double_t z){fV[0]-=x;fV[1]-=y;fV[2]-=z;}
// Orientation
- virtual Double_t Phi() const {return TMath::ATan2(Py(), Px()) * TMath::RadToDeg();} // degrees
- virtual Double_t Theta() const {return TMath::ATan2(Pt(), Pz()) * TMath::RadToDeg();} // degrees
+ virtual Double_t Phi() const {return TMath::ATan2(Py(), Px());}
+ virtual Double_t PhiDeg() const {return Phi() * TMath::RadToDeg();}
+ virtual Double_t Theta() const {return TMath::ATan2(Pt(), Pz()) ;}
+ virtual Double_t ThetaDeg() const {return Theta() * TMath::RadToDeg();}
virtual Double_t Eta() const {return -TMath::Log(TMath::Tan(0.5*Theta()));}
virtual Double_t Y() const {return 0.5*TMath::Log((E() + Pz()) / (E() - Pz()));}
- // Charge
- virtual Short_t Charge() const {return fCharge;}
+ // Kink
virtual Char_t Kink() const {return fKink;}
virtual Bool_t IsKinkMother() const {return (fKink < 0);}
virtual Bool_t IsKinkDaughter() const {return (fKink > 0);}
virtual Bool_t IsKink() const {return (IsKinkMother() || IsKinkDaughter());}
- void SetCharge(Short_t value) {fCharge = value;}
void SetKink(Char_t kink) {fKink = kink;}
+ void SetKinkMother() {fKink = -1;}
+ void SetKinkDaughter() {fKink = 1;}
+ void SetNoKink() {fKink = 0;}
+
+ // Charge
+ virtual Short_t Charge() const {return fCharge;}
+ void SetCharge(Short_t value) {fCharge = value;}
// PID
virtual const Double_t* PID() const {return fPIDWeight;}
const Double_t* PIDProb() const {return fPIDProb;}
- void SetPIDProb(Int_t i, Double_t value) {if(i>=0&&i<AliRsnPID::kSpecies)fPIDProb[i]=value;}
- void SetPIDWeight(Int_t i, Double_t value) {if(i>=0&&i<AliRsnPID::kSpecies)fPIDWeight[i]=value;}
+ void SetPIDProb(Int_t i, Double_t value) {if
(i>=0&&i<AliRsnPID::kSpecies)fPIDProb[i]=value;}
+ void SetPIDWeight(Int_t i, Double_t value) {if
(i>=0&&i<AliRsnPID::kSpecies)fPIDWeight[i]=value;}
static void SetPIDMethod(EPIDMethod method) {fgPIDMethod = method;}
void AssignRealisticPID();
AliRsnPID::EType PIDType(Double_t &prob) const;
// check that contains a given ESD flag
+ void SetFlags(ULong_t flags) {fFlags = flags;}
UInt_t GetFlags() {return fFlags;}
Bool_t CheckFlag(ULong_t flag) {return ((fFlags & flag) == flag);}
- // information getters from objects
- Bool_t Adopt(AliESDtrack *track, AliRsnPIDDefESD pidDef);
- Bool_t Adopt(AliAODTrack *track);
- Bool_t Adopt(AliMCParticle *track);
-
// position in stack/array
Int_t Index() const {return fIndex;}
Int_t Label() const {return fLabel;}
Float_t NSigmaToVertex() const { return fNSigmaToVertex; }
void SetNSigmaToVertex(const Float_t& theValue) { fNSigmaToVertex = theValue; }
-
// ITS/TPC clusters
Int_t NumberOfITSClusters() const {return fITSnum;}
Int_t NumberOfTPCClusters() const {return fTPCnum;}
void Print(Option_t *option = "ALL") const;
void InitMCInfo();
Bool_t InitMCInfo(TParticle *particle);
- Bool_t InitMCInfo(AliMCParticle *mcParticle);
// MC info
AliRsnMCInfo* GetMCInfo() const { return fMCInfo; }
fPVz(0.0),
fPhiMean(0.0),
fMult(0),
+ fPVxMC(0.0),
+ fPVyMC(0.0),
+ fPVzMC(0.0),
fTracks(0x0),
fNoPID(0x0),
fPerfectPID(0x0),
fPVz(event.fPVz),
fPhiMean(event.fPhiMean),
fMult(event.fMult),
+ fPVxMC(event.fPVxMC),
+ fPVyMC(event.fPVyMC),
+ fPVzMC(event.fPVzMC),
fTracks(0x0),
fNoPID(0x0),
fPerfectPID(0x0),
fPVx = event.fPVx;
fPVy = event.fPVy;
fPVz = event.fPVz;
-
+ fPVxMC = event.fPVxMC;
+ fPVyMC = event.fPVyMC;
+ fPVzMC = event.fPVzMC;
+
// other data
fPhiMean = event.fPhiMean;
fMult = event.fMult;
// - mean phi of tracks
//
- if (!fTracks)
+ if (!fTracks)
{
- fMult = 0;
+ fMult = 0;
fPhiMean = 1000.0;
}
- else
+ else
{
fMult = fTracks->GetEntries();
if (fMult < 1) {
}
}
+//_____________________________________________________________________________
+void AliRsnEvent::CorrectTracks()
+{
+//
+// Corrects in all tracks the DCA vertex position using the primary vertex.
+// If present, also MC infos are corrected with primaryvertexMC.
+//
+
+ AliRsnDaughter *track = 0;
+ TObjArrayIter next(fTracks);
+
+ while ( (track = (AliRsnDaughter*)next()) ) {
+ track->ShiftZero(fPVx, fPVy, fPVz);
+ if (track->GetMCInfo()) {
+ track->GetMCInfo()->ShiftZero(fPVxMC, fPVyMC, fPVzMC);
+ }
+ }
+}
+
//_____________________________________________________________________________
Int_t AliRsnEvent::GetNCharged(Char_t sign)
{
pmean = 0.0;
}
}
-
+
return pmean;
}
AliRsnDaughter *leading = GetLeadingParticle(ptMin);
if (!leading) return kFALSE;
-
+
Int_t count = 0;
Double_t angle, angle2Mean;
AliRsnDaughter *trk = 0x0;
TObjArrayIter next(fTracks);
-
+
angleMean = angle2Mean = 0.0;
-
+
while ( (trk = (AliRsnDaughter*)next()) )
{
if (trk == leading) continue;
-
+
angle = leading->AngleTo(trk);
-
+
angleMean += angle;
angle2Mean += angle * angle;
count++;
}
-
+
if (!count) return kFALSE;
-
+
angleMean /= (Double_t)count;
angle2Mean /= (Double_t)count;
angleRMS = TMath::Sqrt(angle2Mean - angleMean*angleMean);
-
+
return kTRUE;
}
void SortTracks() {fTracks->Sort();}
void Print(Option_t *option = "") const;
void MakeComputations();
+ void CorrectTracks();
// Primary vertex
Double_t GetPrimaryVertexX() const {return fPVx;}
Double_t GetPrimaryVertexY() const {return fPVy;}
Double_t GetPrimaryVertexZ() const {return fPVz;}
+ Double_t GetPrimaryVertexXMC() const {return fPVxMC;}
+ Double_t GetPrimaryVertexYMC() const {return fPVyMC;}
+ Double_t GetPrimaryVertexZMC() const {return fPVzMC;}
void GetPrimaryVertex(Double_t &x, Double_t &y, Double_t &z) const {x=fPVx;y=fPVy;z=fPVz;}
Double_t GetVz() const {return GetPrimaryVertexZ();}
void SetPrimaryVertexX(Double_t value) {fPVx = value;}
void SetPrimaryVertexY(Double_t value) {fPVy = value;}
void SetPrimaryVertexZ(Double_t value) {fPVz = value;}
+ void SetPrimaryVertexXMC(Double_t value) {fPVxMC = value;}
+ void SetPrimaryVertexYMC(Double_t value) {fPVyMC = value;}
+ void SetPrimaryVertexZMC(Double_t value) {fPVzMC = value;}
void SetPrimaryVertex(Double_t x, Double_t y, Double_t z) {fPVx=x;fPVy=y;fPVz=z;}
+ void SetPrimaryVertexMC(Double_t x, Double_t y, Double_t z) {fPVxMC=x;fPVyMC=y;fPVzMC=z;}
// Multiplicity
Int_t GetMultiplicity() const {return fMult;}
Int_t GetNCharged(Char_t sign);
-
+
// Mean phi
Double_t GetPhiMean() const {return fPhiMean;}
-
+
// functions for event selection and computation (require on-fly setting of sel parameters)
void SetSelectionPIDType(AliRsnPID::EType type) {fSelPIDType = type;}
void SetSelectionCharge(Char_t charge) {fSelCharge = charge;}
void SetSelectionPIDMethod(AliRsnDaughter::EPIDMethod method) {fSelPIDMethod = method;}
void SetSelectionTrackCuts(AliRsnCutSet *cuts) {fSelCuts = cuts;}
void SetSelection(AliRsnPID::EType pid, Char_t charge, AliRsnDaughter::EPIDMethod meth, AliRsnCutSet *cuts = 0x0);
-
- Bool_t CutPass(AliRsnDaughter *d)
+
+ Bool_t CutPass(AliRsnDaughter *d)
{if (!fSelCuts) return kTRUE; else return fSelCuts->IsSelected(AliRsnCut::kParticle, d);}
-
+
AliRsnDaughter* GetLeadingParticle(Double_t ptMin = 0.0);
Double_t GetAverageMomentum(Int_t &count);
Bool_t GetAngleDistrWRLeading(Double_t &angleMean, Double_t &angleRMS, Double_t ptMin = 0.0);
Double_t fPVz; // vertex
Double_t fPhiMean; // mean "phi" coord of all tracks
Int_t fMult; // track multiplicity
+ Double_t fPVxMC; // position of
+ Double_t fPVyMC; // primary
+ Double_t fPVzMC; // vertex in MC
TClonesArray *fTracks; // collection of particles
AliRsnPIDIndex *fNoPID; // array index only for charged tracks
AliRsnPIDIndex *fPerfectPID; // array index for perfect PID
AliRsnPIDIndex *fRealisticPID; // array index for realistic PID (largest prob)
-
+
AliRsnPID::EType fSelPIDType; //! (for selection/cut functions) particle type
Char_t fSelCharge; //! (for selection/cut functions) particle charge ('0' = both)
AliRsnDaughter::EPIDMethod fSelPIDMethod; //! (for selection/cut functions) PID method used
AliRsnCutSet *fSelCuts; //! (for selection/cut functions) track cuts used
- ClassDef(AliRsnEvent, 2);
+ ClassDef(AliRsnEvent, 3);
};
#endif
//________________________________________________________________________________________
AliRsnEventFunction::AliRsnEventFunction
-(EType type, AliRsnHistoDef *hd, AliRsnDaughter::EPIDMethod pidMethod,
+(EType type, AliRsnHistoDef *hd, AliRsnDaughter::EPIDMethod pidMethod,
AliRsnPID::EType pidType, Char_t sign) :
fType(type),
fPIDMethod(pidMethod),
case kAngleLeadingRMS:
text = "ANGLEADRMS";
break;
+ case kVtResolution:
+ text = "VTRES";
+ break;
+ case kVzResolution:
+ text = "VZRES";
+ break;
default:
AliError("Type not defined");
}
-
+
switch (fPIDMethod)
{
case AliRsnDaughter::kNoPID:
default:
AliError("PID method not defined");
}
-
+
text += AliRsnPID::ParticleName(fPIDType);
text += fCharge;
-
+
if (fEventCuts) {
text += '_';
text += fEventCuts->GetName();
// Fillse the histogram with data contained in a defined pair.
// This method must be overidden by an appropriate definition in each inheriting class.
//
-
+
if (fEventCuts) if (!fEventCuts->IsSelected(AliRsnCut::kEvent, event)) return kFALSE;
-
+
// first of all, set all selection definitions, using the ones in this object
event->SetSelectionPIDType(fPIDType);
event->SetSelectionCharge(fCharge);
//
Int_t count;
- Double_t output, mean = 0.0, rms = 0.0;
+ Double_t output, v[3], vMC[3], vt, vtMC, mean = 0.0, rms = 0.0;
AliRsnDaughter *trk = 0x0;
switch (fType)
fAccept = event->GetAngleDistrWRLeading(mean, rms, fLeadPtMin);
if (fType == kAngleLeadingMean) output = mean; else output = rms;
break;
+ case kVzResolution:
+ case kVtResolution:
+ fAccept = kTRUE;
+ v[0] = event->GetPrimaryVertexX();
+ v[1] = event->GetPrimaryVertexY();
+ v[2] = event->GetPrimaryVertexZ();
+ vMC[0] = event->GetPrimaryVertexXMC();
+ vMC[1] = event->GetPrimaryVertexYMC();
+ vMC[2] = event->GetPrimaryVertexZMC();
+ vt = TMath::Sqrt(v[0]*v[0] + v[1]*v[1]);
+ vtMC = TMath::Sqrt(vMC[0]*vMC[0] + vMC[1]*vMC[1]);
+ if (fType == kVtResolution) return (vt - vtMC); else return (v[2] - vMC[2]);
+ break;
default:
AliError(Form("Type '%d' not supported for EVENT functions", fType));
fAccept = kFALSE;
output = 0.0;
}
-
+
return output;
}
kAverageMomentum,
kAngleLeadingMean,
kAngleLeadingRMS,
-
+ kPtResolution,
+ kVtResolution,
+ kVzResolution,
+
kTypes
};
-
+
AliRsnEventFunction();
AliRsnEventFunction(EType type, AliRsnHistoDef *hd,
AliRsnDaughter::EPIDMethod pidMethod = AliRsnDaughter::kNoPID,
const AliRsnEventFunction& operator=(const AliRsnEventFunction& /*copy*/) { return *this; }
EType fType; // function type
-
+
AliRsnDaughter::EPIDMethod fPIDMethod; // PID method to be used
AliRsnPID::EType fPIDType; // PID species to be used
Char_t fCharge; // charge sign
Double_t fLeadPtMin; // smallest acceptable momentum of leading particle
-
+
Bool_t fAccept; // internal flag to check a computed value
Bool_t fUseBins; // flag to choose if binning is used
TArrayD fBins; // low edge of each bin (upper is the low edge of next bin)
AliRsnCut fBinningCut; // binning cut
AliRsnCut::EType fBinningCutType; // binning cut type
-
+
AliRsnCutSet *fEventCuts; // selection cuts for events
AliRsnCutSet *fTrackCuts; // selection cuts for tracks in each event
// Alberto Pulvirenti (alberto.pulvirenti@ct.infn.it)
//
+#include <Riostream.h>
#include <TList.h>
#include "AliLog.h"
fOutList = new TList();
fOutList->SetOwner();
fOutList->SetName("EventFunctions");
-
+
AliRsnEventFunction *fcn = 0;
TObjArrayIter next(&fEventFunctions);
-
+
while ( (fcn = (AliRsnEventFunction*)next()) )
{
fcn->Init(fOutList);
}
-
- AliInfo("List");
+
+ AliInfo("*************************** List");
fOutList->Print();
- AliInfo("End List");
+ AliInfo("*********************** End List");
}
//________________________________________________________________________
// UserExec() of AliAnalysisTaskSE
//
- if (fEntry++ % 1000 == 0) AliInfo(Form("Event %d", fEntry-1));
+ if (fEntry++ % 100 == 0) cout << "[" << GetName() << "]: event " << fEntry-1 << endl;
- AliRsnEvent *curEvent = GetRsnEventFromInputType();return;
+ AliRsnEvent *curEvent = GetRsnEventFromInputType();
if (!curEvent) return;
AliRsnEventFunction *fcn = 0;
TObjArrayIter next(&fEventFunctions);
-
+
while ( (fcn = (AliRsnEventFunction*)next()) )
{
fcn->Fill(curEvent);
}
-
+
PostData(1, fOutList);
}
fRotAngle(0.0),
fUseBins(kFALSE),
fSkipOutsideInterval(kFALSE),
- fBins(0),
- fBinningCut(),
- fBinningCutType(AliRsnCut::kLastCutType),
+ fNumberOfBinTypes(0),
+// fBins(0),
+// fBinningCut(),
+// fBinningCutType(AliRsnCut::kLastCutType),
fHistoDef(0x0)
{
//
// which must be overridden in any derivate implementation.
//
- Int_t i;
- for (i = 0; i < 100; i++)
- {
- fHisto[i] = 0x0;
- }
+ Int_t i, j;
+ for (j = 0 ; j < kFcnBinTypes; j++)
+ for (i = 0; i < 100; i++)
+ fHisto[j][i] = 0x0;
+
}
//________________________________________________________________________________________
fRotAngle(0.0),
fUseBins(kFALSE),
fSkipOutsideInterval(skipOut),
- fBins(0),
- fBinningCut(),
- fBinningCutType(AliRsnCut::kLastCutType),
+ fNumberOfBinTypes(0),
+// fBins(0),
+// fBinningCut(),
+// fBinningCutType(AliRsnCut::kLastCutType),
fHistoDef(hd)
{
//
// which must be overridden in any derivate implementation.
//
- Int_t i;
- for (i = 0; i < 100; i++)
- {
- fHisto[i] = 0x0;
- }
+ Int_t i, j;
+ for (j = 0 ; j < kFcnBinTypes; j++)
+ for (i = 0; i < 100; i++)
+ fHisto[j][i] = 0x0;
}
//________________________________________________________________________________________
fRotAngle(copy.fRotAngle),
fUseBins(copy.fUseBins),
fSkipOutsideInterval(copy.fSkipOutsideInterval),
- fBins(0),
- fBinningCut(),
- fBinningCutType(AliRsnCut::kLastCutType),
+ fNumberOfBinTypes(copy.fNumberOfBinTypes),
+// fBins(0),
+// fBinningCut(),
+// fBinningCutType(AliRsnCut::kLastCutType),
fHistoDef(copy.fHistoDef)
{
//
// Calls the function to define binning.
//
- Int_t i, n = 100;
- for (i = 0; i < n; i++)
- {
- fHisto[i] = 0x0;
- }
+ Int_t i, j, n;
+ for (j = 0 ; j < kFcnBinTypes; j++)
+ for (i = 0; i < 100; i++)
+ fHisto[j][i] = 0x0;
if (fUseBins)
{
- n = copy.fBins.GetSize();
- Double_t *array = new Double_t[n];
- for (i = 0; i < n; i++) array[i] = copy.fBins[i];
- SetBinningCut(copy.fBinningCutType, copy.fBins.GetSize(), array);
- delete [] array;
+ for (i = 0 ; i < kFcnBinTypes; i++){
+ if (fNumberOfBinTypes<=i) continue;
+ n = copy.fBins[i].GetSize();
+ Double_t *array = new Double_t[n];
+ for (j = 0; j < n; j++) array[j] = copy.fBins[i][j];
+ SetBinningCut(copy.fBinningCutType[i], copy.fBins[i].GetSize(), array,i,kTRUE);
+ delete [] array;
+ }
}
}
//________________________________________________________________________________________
// For the sake of security, all pointers are also set explicitly to NULL.
//
- Int_t i;
- for (i = 0; i < 100; i++)
- {
- delete fHisto[i];
- fHisto[i] = 0x0;
- }
+ Int_t i, j;
+ for (j = 0 ; j < kFcnBinTypes; j++)
+ for (i = 0; i < 100; i++)
+ {
+ delete fHisto[j][i];
+ fHisto[j][i] = 0x0;
+ }
}
//________________________________________________________________________________________
// a general histogram is always added,
// which overrides the binning and collects everything
- fHisto[0] = new TH1D(histoName, histoTitle, nbins, min, max);
- fHisto[0]->Sumw2();
- histos->AddLast(fHisto[0]);
+
+ fHisto[0][0] = new TH1D(histoName, histoTitle, nbins, min, max);
+ fHisto[0][0]->Sumw2();
+ histos->AddLast(fHisto[0][0]);
// if requested a binning w.r. to some cut variable, histograms are added
// for that in this part of the method (one per each bin)
Char_t hName[255];
Char_t hTitle[255];
+ Int_t j;
if (fUseBins)
{
- for (ibin = 0, i = 1; ibin < fBins.GetSize() - 1; ibin++, i++)
- {
- sprintf(hName, "%s[%.2f-%.2f]", histoName, fBins[ibin], fBins[ibin+1]);
- sprintf(hTitle, "%s [%.2f-%.2f]", histoTitle, fBins[ibin], fBins[ibin+1]);
- fHisto[i] = new TH1D(hName, hTitle, nbins, min, max);
- fHisto[i]->Sumw2();
- histos->AddLast(fHisto[i]);
+ for (j = 0 ; j < kFcnBinTypes; j++){
+ if (fNumberOfBinTypes<=j) continue;
+ for (ibin = 0, i = 1; ibin < fBins[j].GetSize() - 1; ibin++, i++)
+ {
+ sprintf(hName, "%s_%d%02d_[%.2f-%.2f]", histoName, j,i,fBins[j][ibin], fBins[j][ibin+1]);
+ sprintf(hTitle, "%s [%.2f-%.2f]", histoTitle, fBins[j][ibin], fBins[j][ibin+1]);
+// AliInfo(Form("Adding %s",hName));
+ fHisto[j][i] = new TH1D(hName, hTitle, nbins, min, max);
+ fHisto[j][i]->Sumw2();
+ histos->AddLast(fHisto[j][i]);
+ }
}
}
// a general histogram is always added,
// which overrides the binning and collects everything
- fHisto[0] = new TH1D(histoName, histoTitle, nbins, min, max);
- histos->AddLast(fHisto[0]);
+ fHisto[0][0] = new TH1D(histoName, histoTitle, nbins, min, max);
+ histos->AddLast(fHisto[0][0]);
// if requested a binning w.r. to some cut variable, histograms are added
// for that in this part of the method (one per each bin)
Char_t hName[255];
Char_t hTitle[255];
+ Int_t j;
if (fUseBins)
{
- for (ibin = 0, i = 1; ibin < fBins.GetSize() - 1; ibin++, i++)
- {
- sprintf(hName, "%s[%.2f-%.2f]", histoName, fBins[ibin], fBins[ibin+1]);
- sprintf(hTitle, "%s [%.2f-%.2f]", histoTitle, fBins[ibin], fBins[ibin+1]);
- fHisto[i] = new TH1D(hName, hTitle, nbins, min, max);
- histos->AddLast(fHisto[i]);
+ for (j = 0 ; j < kFcnBinTypes; j++){
+ if (fNumberOfBinTypes<=j) continue;
+ for (ibin = 0, i = 1; ibin < fBins[j].GetSize() - 1; ibin++, i++)
+ {
+
+ sprintf(hName, "%s_%d_%02d[%.2f-%.2f]", histoName,j,i, fBins[j][ibin], fBins[j][ibin+1]);
+ sprintf(hTitle, "%s [%.2f-%.2f]", histoTitle, fBins[j][ibin], fBins[j][ibin+1]);
+// AliInfo(Form("Adding %s",hName));
+ fHisto[j][i] = new TH1D(hName, hTitle, nbins, min, max);
+ histos->AddLast(fHisto[j][i]);
+ }
}
}
}
//________________________________________________________________________________________
void AliRsnFunction::SetBinningCut
-(AliRsnCut::EType type, Double_t min, Double_t max, Double_t step)
+(AliRsnCut::EType type, Double_t min, Double_t max, Double_t step,Int_t index,Bool_t IsCopyConstructor)
{
//
// Set fixed bins
//
- fUseBins = kTRUE;
+ if (index >= kFcnBinTypes) {
+ AliError(Form("We support only %d Binning cuts(0-%d). Skipping...",kFcnBinTypes,kFcnBinTypes-1));
+ return;
+ }
+
+ if (!IsCopyConstructor){
+ // TODO if some one sets indexes 0,2,3 it is a bug here(i'll solve it)
+ if (index == fNumberOfBinTypes)
+ fNumberOfBinTypes++;
+ else {
+ AliError(Form("Wrong index %d. fUseBins is set to kFALSE",index));
+// fUseBins = kFALSE;
+ return;
+ }
+ }
+ fUseBins = kTRUE;
Int_t i, nBins = (Int_t)((max - min) / step) + 1;
- fBinningCutType = type;
- fBins.Set(nBins);
+ fBinningCutType[index] = type;
+ fBins[index].Set(nBins);
for (i = 0; i < nBins; i++)
{
- fBins[i] = min + (Double_t)i * step;
+ fBins[index][i] = min + (Double_t)i * step;
}
}
//________________________________________________________________________________________
void AliRsnFunction::SetBinningCut
-(AliRsnCut::EType type, Int_t nbins, Double_t *bins)
+(AliRsnCut::EType type, Int_t nbins, Double_t *bins,Int_t index,Bool_t IsCopyConstructor)
{
//
// Set variable bins
//
- fUseBins = kTRUE;
+ if (index >= kFcnBinTypes) {
+ AliError(Form("We support only %d Binning cuts(0-%d). Skipping...",kFcnBinTypes,kFcnBinTypes-1));
+ return;
+ }
+ if (!IsCopyConstructor){
+ // TODO if some one sets indexes 0,2,3 it is a bug here(i'll solve it)
+ if (index >= fNumberOfBinTypes)
+ fNumberOfBinTypes++;
+ else {
+ AliError(Form("Wrong index %d (%d). fUseBins is set to kFALSE",index,fNumberOfBinTypes));
+ // fUseBins = kFALSE;
+ return;
+ }
+ }
+ fUseBins = kTRUE;
Int_t i;
- fBinningCutType = type;
- fBins.Set(nbins);
+ fBinningCutType[index] = type;
+ fBins[index].Set(nbins);
for (i = 0; i < nbins; i++)
{
- fBins[i] = bins[i];
+ fBins[index][i] = bins[i];
}
}
}
//________________________________________________________________________________________
-Bool_t AliRsnFunction::Fill(AliRsnPairParticle *pair, AliRsnPairDef *ref, Double_t weight)
+Bool_t AliRsnFunction::Fill(AliRsnPairParticle *pair, AliRsnPairDef *ref)
{
//
// Fillse the histogram with data contained in a defined pair.
}
// fill global histogram
- if (weight == 0.0) fHisto[0]->Fill(value);
- else fHisto[0]->Fill(value, weight);
+ fHisto[0][0]->Fill(value);
// if bins are allocated, find right one and fill it
if (fUseBins)
{
- Int_t i, ibin;
- for (ibin = 0, i = 1; ibin < fBins.GetSize() - 1; ibin++, i++)
- {
- if (!fHisto[i]) continue;
- fBinningCut.SetCutValues(fBinningCutType, fBins[ibin], fBins[ibin+1]);
- if (fBinningCut.IsSelected(AliRsnCut::kPair, pair))
+ Int_t i, j, ibin;
+ for (j = 0 ; j < kFcnBinTypes; j++){
+ if (fNumberOfBinTypes<=j) continue;
+ for (ibin = 0, i = 1; ibin < fBins[j].GetSize() - 1; ibin++, i++)
{
- if (weight == 0.0) fHisto[i]->Fill(value);
- else fHisto[i]->Fill(value, weight);
- break;
+ if (!fHisto[j][i]) continue;
+ fBinningCut[j].SetCutValues(fBinningCutType[j], fBins[j][ibin], fBins[j][ibin+1]);
+ if (fBinningCut[j].IsSelected(AliRsnCut::kPair, pair))
+ {
+ fHisto[j][i]->Fill(value);
+ break;
+ }
}
}
}
kEtaSpectrum,
kFcnTypes
};
+ enum EFcnBinTypes
+ {
+ kFcnBinTypes=5
+ };
AliRsnFunction();
AliRsnFunction(EFcnType type, AliRsnHistoDef *hd, Bool_t skipOut = kTRUE);
TString GetFcnName();
TString GetFcnTitle();
- void SetBinningCut(AliRsnCut::EType type, Double_t min, Double_t max, Double_t step);
- void SetBinningCut(AliRsnCut::EType type, Int_t nbins, Double_t *bins);
+ void SetBinningCut(AliRsnCut::EType type, Double_t min, Double_t max, Double_t step,Int_t index=0,Bool_t IsCopyConstructor=kFALSE);
+ void SetBinningCut(AliRsnCut::EType type, Int_t nbins, Double_t *bins,Int_t index=0,Bool_t IsCopyConstructor=kFALSE);
void SetHistoDef(AliRsnHistoDef *def) {fHistoDef = def;}
void SetRotationAngle(Double_t rotAngle) {fRotAngle = rotAngle;}
// working routines
TList* Init(const char *histoName, const char *histoTitle);
void Init(const char *histoName, const char *histoTitle, TList *tgt);
- Bool_t Fill(AliRsnPairParticle *pair, AliRsnPairDef *ref, Double_t weight = 0.0);
+ Bool_t Fill(AliRsnPairParticle *pair, AliRsnPairDef *ref);
Double_t FcnValue(AliRsnPairParticle *pair, AliRsnPairDef *ref);
private:
Double_t FcnResolution(AliRsnPairParticle *pair, AliRsnPairDef *pd);
- EFcnType fFcnType; // function type
+ EFcnType fFcnType; // function type
- Double_t fRotAngle; // rotation angle (for "rotated" invMass)
+ Double_t fRotAngle; // rotation angle (for "rotated" invMass)
- Bool_t fUseBins; // flag to choose if binning is used
- Bool_t fSkipOutsideInterval; // skip pairs which fall outside histogram interval
+ Bool_t fUseBins; // flag to choose if binning is used
+ Bool_t fSkipOutsideInterval; // skip pairs which fall outside histogram interval
- TArrayD fBins; // low edge of each bin (upper is the low edge of next bin)
- AliRsnCut fBinningCut; // binning cut
- AliRsnCut::EType fBinningCutType;// binning cut type
+ Int_t fNumberOfBinTypes; // number of binning types
+ TArrayD fBins[kFcnBinTypes]; // low edge of each bin (upper is the low edge of next bin)
+ AliRsnCut fBinningCut[kFcnBinTypes]; // binning cut
+ AliRsnCut::EType fBinningCutType[kFcnBinTypes]; // binning cut type
- AliRsnHistoDef *fHistoDef; // definitions for histogram
- TH1D *fHisto[100]; // binned histograms
+ AliRsnHistoDef *fHistoDef; // definitions for histogram
+ TH1D *fHisto[kFcnBinTypes][100]; // binned histograms
// ROOT dictionary
ClassDef(AliRsnFunction, 1)
//_____________________________________________________________________________
AliRsnMCInfo::AliRsnMCInfo() :
- TObject(), fEnergy(0), fPDG(0), fMother(-1), fMotherPDG(0)
+ TObject(), fEnergy(0), fPDG(0), fMother(-1), fMotherPDG(0)
{
//
// Default constructor.
//
Int_t i;
- for (i = 0; i < 3; i++) fP[i] = 0.0;
+ for (i = 0; i < 3; i++) fP[i] = fV[i] = 0.0;
}
//_____________________________________________________________________________
AliRsnMCInfo::AliRsnMCInfo(const AliRsnMCInfo & copy) :
- TObject(copy),
- fEnergy(copy.fEnergy),
- fPDG(copy.fPDG),
- fMother(copy.fMother),
- fMotherPDG(copy.fMotherPDG)
+ TObject(copy),
+ fEnergy(copy.fEnergy),
+ fPDG(copy.fPDG),
+ fMother(copy.fMother),
+ fMotherPDG(copy.fMotherPDG)
{
//
// Copy constructor.
//
Int_t i;
- for (i = 0; i < 3; i++) fP[i] = copy.fP[i];
+ for (i = 0; i < 3; i++) {
+ fP[i] = copy.fP[i];
+ fV[i] = copy.fV[i];
+ }
}
//_____________________________________________________________________________
// If the argument is NULL, nothing is done and an error message is returned.
//
- if (!particle)
- {
- AliError("NULL argument passed. Nothing done.");
- return;
- }
+ if (!particle) return;
fP[0] = particle->Px();
fP[1] = particle->Py();
fP[2] = particle->Pz();
+ fV[0] = particle->Vx();
+ fV[1] = particle->Vy();
+ fV[2] = particle->Vz();
+
fEnergy = particle->Energy();
fPDG = particle->GetPdgCode();
void SetP(Double_t px, Double_t py, Double_t pz) {SetPx(px); SetPy(py); SetPz(pz);}
void SetE(Double_t e) {fEnergy = e;}
+ Double_t Vx() const {return fV[0];}
+ Double_t Vy() const {return fV[1];}
+ Double_t Vz() const {return fV[2];}
+ Double_t Dr() const {return TMath::Sqrt(Vx()*Vx() + Vy()*Vy());}
+ void ShiftZero(Double_t x, Double_t y, Double_t z){fV[0]-=x;fV[1]-=y;fV[2]-=z;}
+
+ void SetVx(Double_t value) {fV[0] = value;}
+ void SetVy(Double_t value) {fV[1] = value;}
+ void SetVz(Double_t value) {fV[2] = value;}
+ void SetV(Double_t x, Double_t y, Double_t z) {SetVx(x); SetVy(y); SetVz(z);}
+
Int_t PDG() const {return fPDG;}
Int_t Mother() const {return fMother;}
Short_t MotherPDG() const {return fMotherPDG;}
private:
Double_t fP[3]; // MC momentum
+ Double_t fV[3]; // MC position
Double_t fEnergy; // MC energy
Int_t fPDG; // PDG code
Int_t fMother; // GEANT label of mother particle
//
#include <Riostream.h>
-#include "TObjArray.h"
+#include <TObjArray.h>
#include "AliLog.h"
ClassImp(AliRsnPair)
//_____________________________________________________________________________
-AliRsnPair::AliRsnPair
-(EPairType type, AliRsnPairDef *def, Int_t mixNum) :
+AliRsnPair::AliRsnPair(EPairType type, AliRsnPairDef *def) :
TObject(),
fIsMixed(kFALSE),
- fUseMC(kFALSE),
- fIsLikeSign(kFALSE),
- fMixNum(mixNum),
- fMixingCut(0x0),
- fMatched(mixNum),
- fPairDef(def),
fPairType(type),
- fTypePID(AliRsnDaughter::kRealistic),
+ fPIDMethod(AliRsnDaughter::kRealistic),
+ fPairDef(def),
fCutMgr(0),
fFunctions("AliRsnFunction", 0)
{
//
SetUp(type);
- if (!fIsMixed) fMatched.Set(1);
}
//_____________________________________________________________________________
AliRsnPair::~AliRsnPair()
switch (type)
{
case kNoPID:
- SetAllFlags(AliRsnDaughter::kNoPID, kFALSE, kFALSE);
+ SetAllFlags(AliRsnDaughter::kNoPID, kFALSE);
break;
case kNoPIDMix:
- SetAllFlags(AliRsnDaughter::kNoPID, kTRUE, kFALSE);
+ SetAllFlags(AliRsnDaughter::kNoPID, kTRUE);
break;
case kRealisticPID:
- SetAllFlags(AliRsnDaughter::kRealistic, kFALSE, kFALSE);
+ SetAllFlags(AliRsnDaughter::kRealistic, kFALSE);
break;
case kRealisticPIDMix:
- SetAllFlags(AliRsnDaughter::kRealistic, kTRUE, kFALSE);
+ SetAllFlags(AliRsnDaughter::kRealistic, kTRUE);
break;
case kPerfectPID:
- // SetAllFlags (AliRsnDaughter::kPerfect, kFALSE, kFALSE);
- SetAllFlags(AliRsnDaughter::kPerfect, kFALSE, kTRUE);
+ SetAllFlags (AliRsnDaughter::kPerfect, kFALSE);
break;
case kPerfectPIDMix:
- // SetAllFlags (AliRsnDaughter::kPerfect, kTRUE, kFALSE);
- SetAllFlags(AliRsnDaughter::kPerfect, kTRUE, kTRUE);
+ SetAllFlags (AliRsnDaughter::kPerfect, kTRUE);
break;
default :
- AliWarning("Wrong type selected: setting up for kRealisticPID.");
- SetAllFlags(AliRsnDaughter::kRealistic, kFALSE, kFALSE);
+ AliWarning("Wrong type selected: setting up for realistic PID - no mixing.");
+ SetAllFlags(AliRsnDaughter::kRealistic, kFALSE);
break;
}
}
-//_____________________________________________________________________________
-void AliRsnPair::SetAllFlags(AliRsnDaughter::EPIDMethod pidType, Bool_t isMix, Bool_t useMC)
-{
-//
-// Sets up all flags values
-//
-
- fTypePID = pidType;
- fIsMixed = isMix;
- fUseMC = useMC;
-}
-
-//_____________________________________________________________________________
-void AliRsnPair::Init()
-{
-//
-// Init pair
-//
-
- fIsLikeSign = fPairDef->IsLikeSign();
- Print();
-}
-
//_____________________________________________________________________________
void AliRsnPair::Print(Option_t* /*option*/) const
{
//
// Prints info about pair
//
+
AliInfo(Form("%s", GetPairHistTitle(0x0).Data()));
AliInfo(Form("PDG %d %d", AliRsnPID::PDGCode(fPairDef->GetType(0)),
AliRsnPID::PDGCode(fPairDef->GetType(1))));
AliInfo(Form("Masses %f %f", fPairDef->GetMass(0), fPairDef->GetMass(1)));
AliInfo(Form("Number of functions %d", fFunctions.GetEntries()));
- switch(fTypePID) {
+
+ switch(fPIDMethod)
+ {
case AliRsnDaughter::kNoPID:
AliInfo("PID method: none");
break;
break;
default:
AliInfo("PID method: undefined");
- }
-}
-
-//_____________________________________________________________________________
-void AliRsnPair::ProcessPair(AliRsnEventBuffer *buf)
-{
-//
-// Process current event in the passed buffer.
-// If this pair is a single-event pair, only that event is processed,
-// otherwise, if it is a mixing pair, all good matches are found and mixed.
-//
-
- if (fIsMixed && buf->GetEventsBufferIndex() >= fMixNum) ProcessPairMix(buf);
- else ProcessPairSingle(buf);
-
- /*
- // track type/charge #0 in pairDef is taken from this event,
- // track tipe/charge #1 is taken from the matched event
- AliRsnEvent *e1 = buf->GetCurrentEvent();
- if (!e1) return;
- TArrayI* array1 = e1->GetTracksArray(fTypePID, fPairDef->GetCharge(0), fPairDef->GetType(0));
-
- // reset array of matched events
- Int_t i;
- for (i = 0; i < fMixNum; i++) fMatched[i] = -1;
-
- // here searches for all matched events:
- // in case of single-event pairs, there will be only one, equal to current,
- // otherwise there will be many
- Int_t lastOkEvent = buf->IndexOf(e1) - 1;
- Double_t nMatched = 0.0;
- if (fIsMixed) {
- for (i = 0; i < fMixNum; i++) {
- // find other event by event cut
- AliRsnEvent *e2 = FindEventByEventCut(buf, lastOkEvent);
- if (!e2) break;
- fMatched[i] = lastOkEvent;
- lastOkEvent--;
- nMatched++;
- }
- if (fMatched.GetSize() < 0) {
- AliWarning(Form("Event #%d: found no events to match", buf->IndexOf(e1)));
- return;
- }
}
- else {
- fMatched[0] = buf->IndexOf(e1);
- nMatched = 1;
- }
-
- // in order to balance the problem that different events could be matched
- // a different number of times (between 1 and the maximum allowd quantity),
- // for each event, the histograms are filled with a variable weight depending
- // on the number of matches found:
- Double_t weight = (Double_t)fMixNum / (Double_t)nMatched;
-
- // now that matched events are found, they are used to fill the histograms
- TArrayI* array2 = 0;
- for (i = 0; i < fMixNum; i++) {
- if (fMatched[i] < 0) break;
- AliRsnEvent *e2 = buf->GetEvent(fMatched[i]);
- array2 = e2->GetTracksArray(fTypePID, fPairDef->GetCharge(1), fPairDef->GetType(1));
- if (fIsMixed) LoopPair(e1, array1, e2, array2, weight);
- else LoopPair(e1, array1, e2, array2);
- }
- */
}
//_____________________________________________________________________________
-void AliRsnPair::ProcessPairSingle(AliRsnEventBuffer *buf)
+void AliRsnPair::ProcessPair(AliRsnEvent *ev1, AliRsnEvent *ev2)
{
//
-// SINGLE EVENT PROCESSING
-// This function fills the functions' histograms using tracks
-// from the same event for both components of the defined pair.
-// This function is used for signal, like-sign and rotated background.
-//
-
- AliRsnEvent *event = buf->GetCurrentEvent();
- if (!event) return;
-
- TArrayI* array1 = event->GetTracksArray(fTypePID, fPairDef->GetCharge(0), fPairDef->GetType(0));
- TArrayI* array2 = event->GetTracksArray(fTypePID, fPairDef->GetCharge(1), fPairDef->GetType(1));
-
- LoopPair(event, array1, event, array2);
-}
-
-//_____________________________________________________________________________
-void AliRsnPair::ProcessPairMix(AliRsnEventBuffer *buf)
-{
+// Fills the functions' histograms using tracks from passed events.
+// What tracks are taken in each event depend from the order of
+// track types defined in the AliRsnPairDef for this object:
+// - tracks of type 0 are the ones stored as pairDef data members with index [0]
+// ---> taken from first argument (ev1)
+// - tracks of type 1 are the ones stored as pairDef data members with index [1]
+// ---> taken from second argument (ev2)
//
-// MIXED EVENT PROCESSING
-// This function fills the functions' histograms using
-// tracks from different events.
-// Used for event mixing.
+// When doing single-event analysis (e.g. signal, like-sign), second argument
+// can be NULL, and it will be forced to point to the same object of first one.
//
- // track type/charge #0 in pairDef is taken from this event,
- // track tipe/charge #1 is taken from the matched events
- AliRsnEvent *e1 = buf->GetCurrentEvent();
- if (!e1) return;
- TArrayI* array10 = e1->GetTracksArray(fTypePID, fPairDef->GetCharge(0), fPairDef->GetType(0));
- TArrayI* array11 = e1->GetTracksArray(fTypePID, fPairDef->GetCharge(1), fPairDef->GetType(1));
+ if (!ev2) ev2 = ev1;
- // find matched events
- Int_t i, iev = buf->GetEventsBufferIndex();
- Int_t nMatched = FindMatchedEvents(iev, buf);
- if (!nMatched) {
- AliWarning(Form("Event #%d: found no events to match", iev));
- return;
- }
- else if (nMatched < fMixNum) {
- AliWarning(Form("Event #%d: found only %d events to match", iev, nMatched));
- return;
- }
+ TArrayI *array1 = ev1->GetTracksArray(fPIDMethod, fPairDef->GetCharge(0), fPairDef->GetType(0));
+ TArrayI *array2 = ev2->GetTracksArray(fPIDMethod, fPairDef->GetCharge(1), fPairDef->GetType(1));
- /*
- else {
- TString str("Matched events: ");
- for (i = 0; i < nMatched; i++) str.Append(Form("%d ", fMatched[i]));
- AliInfo(Form("Event #%d: %s", iev, str.Data()));
- }
- */
-
- // in order to balance the problem that different events could be matched
- // a different number of times (between 1 and the maximum allowd quantity),
- // for each event, the histograms are filled with a variable weight depending
- // on the number of matches found:
- // Double_t weight = 0.5 * (Double_t)fMixNum / (Double_t)nMatched;
-
- // now that matched events are found, they are used to fill the histograms
- TArrayI *array20 = 0, *array21 = 0;
- for (i = 0; i < fMixNum; i++) {
- if (fMatched[i] < 0) break;
- AliRsnEvent *e2 = buf->GetEvent(fMatched[i]);
- array20 = e2->GetTracksArray(fTypePID, fPairDef->GetCharge(0), fPairDef->GetType(0));
- array21 = e2->GetTracksArray(fTypePID, fPairDef->GetCharge(1), fPairDef->GetType(1));
- // track type/charge #0 from event1, track type/charge #1 from event2
- LoopPair(e1, array10, e2, array21);
- // track type/charge #1 from event1, track type/charge #0 from event2
- LoopPair(e1, array11, e2, array20);
- }
-}
-
-//_____________________________________________________________________________
-Int_t AliRsnPair::FindMatchedEvents(Int_t evIndex, AliRsnEventBuffer *buf)
-{
-//
-// Resets the 'fMatched' data member and stores into it
-// a number of well-matched events found in the buffer,
-// (maximum amount = fMixNum).
-// The argument is the index of the event to be matched,
-// in the event buffer.
-//
-
- // reset array of matched events
- Int_t i;
- for (i = 0; i < fMixNum; i++) fMatched[i] = -1;
-
- // starts from the position behind the one
- // of the event to be matched and goes backward
- // if it reaches the value 0, stops
- AliRsnEvent *eventToBeMatched = buf->GetCurrentEvent();
- AliRsnEvent *matchEvent = 0x0;
- Int_t checkIndex = evIndex - 1, curIndex = 0;
- for (;;checkIndex--) {
- if (checkIndex < 0) checkIndex = buf->GetEventsBufferSize() - 1;
- if (checkIndex == evIndex) break;
- matchEvent = buf->GetEvent(checkIndex);
- if (!matchEvent) continue;
- if (fMixingCut) {
- if (fMixingCut->IsSelected(AliRsnCut::kMixEvent, eventToBeMatched, matchEvent)) continue;
- }
- // assign to current array slot the matched event
- // and increment current slot and stops if it exceeds array size
- fMatched[curIndex++] = checkIndex;
- if (curIndex >= fMixNum) break;
- }
-
- // returns the current index value,
- // which is also the number of matched events found
- return curIndex;
-}
-
-//_____________________________________________________________________________
-AliRsnEvent * AliRsnPair::FindEventByEventCut(AliRsnEventBuffer *buf, Int_t& num)
-{
-//
-// For now it just returns num events before current event
-// in buffer (buf)
-// TODO event cut selection
-//
-
- AliRsnEvent *returnEvent = 0x0;
-
- if (fIsMixed)
- {
- //returnEvent = buf->GetEvent(buf->GetEventsBufferIndex() - num);
- returnEvent = buf->GetNextGoodEvent(num, fMixingCut);
- }
- else
- {
- returnEvent = buf->GetCurrentEvent();
- }
-
- return returnEvent;
+ LoopPair(ev1, array1, ev2, array2);
}
//_____________________________________________________________________________
void AliRsnPair::LoopPair
-(AliRsnEvent * ev1, TArrayI * a1, AliRsnEvent * ev2, TArrayI * a2, Double_t weight)
+(AliRsnEvent * ev1, TArrayI * a1, AliRsnEvent * ev2, TArrayI * a2)
{
//
// Loop on all pairs of tracks of the defined types/charges,
// using the arrays of indexes and the events containing them.
+// This method is private, for safety reasons.
//
if (!a1) {AliDebug(4, "No TArrayI 1 from currentEvent->GetTracksArray(...)"); return;}
if (!a2) {AliDebug(4, "No TArrayI 2 from currentEvent->GetTracksArray(...)"); return;}
- AliRsnDaughter::SetPIDMethod(fTypePID);
+ // cuts on events
+ if (!CutPass(ev1) || !CutPass(ev2)) return;
+
+ AliRsnDaughter::SetPIDMethod(fPIDMethod);
AliRsnDaughter *daughter1 = 0;
AliRsnDaughter *daughter2 = 0;
AliRsnFunction *fcn = 0;
+
+ Bool_t isLikeSign = fPairDef->IsLikeSign();
Int_t j, startj = 0;
+
for (Int_t i = 0; i < a1->GetSize(); i++)
{
// get track #1
daughter2 = 0;
// check starting index for searching the event:
// for like-sign pairs we avoid duplicating the pairs
- if (fIsLikeSign) startj = i+1; else startj = 0;
+ if (isLikeSign) startj = i+1; else startj = 0;
// AliInfo(Form("%d",startj));
// loop on event for all track #2 to be combined with the found track #1
for (j = startj; j < a2->GetSize(); j++)
// fill all histograms
TObjArrayIter nextFcn(&fFunctions);
while ( (fcn = (AliRsnFunction*)nextFcn()) ) {
- fcn->Fill(&pairParticle, fPairDef, weight);
+ fcn->Fill(&pairParticle, fPairDef);
}
}
}
TList * AliRsnPair::GenerateHistograms(TString prefix)
{
//
-// Generates needed histograms
+// Generates needed histograms, giving them a name based on
+// the flags defined here, on the pair definition, and attaches
+// a prefix to it, according to the argument.
+//
+// All generated histograms are stored into the output TList.
//
TList *list = new TList();
sprintf(hName, "%s_%s", prefix.Data(), GetPairHistName(fcn).Data());
sprintf(hTitle, "%s", GetPairHistTitle(fcn).Data());
TList *histos = fcn->Init(hName, hTitle);
- histos->Print();
+ //histos->Print();
list->Add(histos);
}
void AliRsnPair::GenerateHistograms(TString prefix, TList *tgt)
{
//
-// Generates needed histograms
+// Generates needed histograms, giving them a name based on
+// the flags defined here, on the pair definition, and attaches
+// a prefix to it, according to the argument.
+//
+// All generated histograms are stored into the TList passed as argument
//
if (!tgt) {
//
// Returns type name, made with particle names ant chosen PID
//
+
switch (type)
{
case kNoPID : return ("NOPID_");break;
case kRealisticPIDMix : return ("REALISTICMIX_");break;
case kPerfectPID : return ("PERFECT_");break;
case kPerfectPIDMix : return ("PERFECTMIX_");break;
- case kTruePairs : return ("TRUEPAIRS_"); break;
default:
AliWarning("Unrecognized value of EPairTypeName argument");
break;
//
// Retruns pair name
//
+
TString sName;
sName += GetPairTypeName(fPairType);
sName += fPairDef->GetPairName();
TString AliRsnPair::GetPairHistName(AliRsnFunction *fcn, TString text) const
{
//
-// Returns eff. mass histogram name
+// Returns definitive histogram name
//
TString sName;
TString AliRsnPair::GetPairHistTitle(AliRsnFunction *fcn, TString text) const
{
//
-// Returns eff. mass histogram title
+// Returns definitive histogram title
//
TString sTitle;
//
// Adds a new computing function
//
+
Int_t size = fFunctions.GetEntries();
- new(fFunctions[size]) AliRsnFunction(*fcn);
+ new (fFunctions[size]) AliRsnFunction(*fcn);
}
//________________________________________________________________________________________
// Check if the AliRsnDaughter argument pass its cuts.
// If the cut data member is not initialized for it, returns kTRUE.
//
+
if (!fCutMgr) return kTRUE;
else return fCutMgr->IsSelected(AliRsnCut::kParticle, d);
}
//
// Check if the AliRsnPairParticle argument pass its cuts.
// If the cut data member is not initialized for it, returns kTRUE.
-// In this case, another separate check which could be necessary
-// concerns the possibility that the two tracks are a "true pair" of
-// daughters of the same resonance. If the corresponding flag is set,
-// this further check is done, and the method returns kTRUE only
-// when also this check is passed.
//
if (!fCutMgr) return kTRUE;
// Check if the AliRsnEvent argument pass its cuts.
// If the cut data member is not initialized for it, returns kTRUE.
//
+
if (!fCutMgr) return kTRUE;
else return fCutMgr->IsSelected(AliRsnCut::kEvent, e);
}
kNoPID = 0, kNoPIDMix,
kRealisticPID, kRealisticPIDMix,
kPerfectPID, kPerfectPIDMix,
- kTruePairs,
kPairTypes
};
- AliRsnPair(EPairType type = kRealisticPID, AliRsnPairDef *def = 0, Int_t mixNum = 1);
+ AliRsnPair(EPairType type = kRealisticPID, AliRsnPairDef *def = 0);
~AliRsnPair();
- void Init();
void Print(Option_t *option = "") const;
- void ProcessPair(AliRsnEventBuffer *buf);
+ void ProcessPair(AliRsnEvent *ev1, AliRsnEvent *ev2 = 0);
void SetCutMgr(AliRsnCutMgr* theValue) { fCutMgr = theValue; }
- void SetMixingCut(AliRsnCutSet* theValue) { fMixingCut = theValue; }
void AddFunction(AliRsnFunction *fcn);
TList* GenerateHistograms(TString prefix = "");
void GenerateHistograms(TString prefix, TList *tgt);
+ Bool_t IsMixed() {return fIsMixed;}
+ Bool_t IsPairEqual() {if (fPIDMethod == AliRsnDaughter::kNoPID) return (fPairDef->IsLikeSign());
+ else return (fPairDef->IsLikeSign() && fPairDef->HasEqualTypes());}
+
TString GetPairTypeName(EPairType type) const;
TString GetPairName() const;
TString GetPairHistName(AliRsnFunction *fcn, TString text = "") const;
private:
AliRsnPair (const AliRsnPair ©) : TObject(copy),
- fIsMixed(kFALSE),fUseMC(kFALSE),fIsLikeSign(kFALSE),fMixNum(1),fMixingCut(0x0),fMatched(0),
- fPairDef(0x0),fPairType(kPairTypes),fTypePID(AliRsnDaughter::kRealistic),
- fCutMgr(0x0),fFunctions("AliRsnFunction",0) {}
+ fIsMixed(kFALSE),fPairType(kPairTypes),fPIDMethod(AliRsnDaughter::kRealistic),
+ fPairDef(0x0),fCutMgr(0x0),fFunctions("AliRsnFunction",0) {}
AliRsnPair& operator=(const AliRsnPair&) {return *this;}
- void SetUp(EPairType type); // sets up all flags
- void SetAllFlags(AliRsnDaughter::EPIDMethod pidType,Bool_t isMix, Bool_t useMC);
- AliRsnEvent* FindEventByEventCut(AliRsnEventBuffer *buf, Int_t & num);
- Int_t FindMatchedEvents(Int_t evIndex, AliRsnEventBuffer *buf);
- void ProcessPairSingle(AliRsnEventBuffer *buf);
- void ProcessPairMix(AliRsnEventBuffer *buf);
- void LoopPair(AliRsnEvent *ev1, TArrayI *a1, AliRsnEvent *ev2, TArrayI *a2, Double_t weight = 0.0);
-
- Bool_t CutPass(AliRsnDaughter *d);
- Bool_t CutPass(AliRsnPairParticle *p);
- Bool_t CutPass(AliRsnEvent *e);
-
- // flags & integer data
- Bool_t fIsMixed; // doing event-mixing ?
- Bool_t fUseMC; // using MC inv. mass ?
- Bool_t fIsLikeSign; // is a like-sign pair ?
- Int_t fMixNum; // number of mixed events
- AliRsnCutSet *fMixingCut; // cut for event mixing
-
- // work management
- TArrayI fMatched; // array with indexes of matched events (used in mixing)
- AliRsnPairDef *fPairDef; // pair definition (particles, charges)
- EPairType fPairType; // pair type (PID + mixing or not)
- AliRsnDaughter::EPIDMethod fTypePID; // pid type variable for single track
- AliRsnCutMgr *fCutMgr; // cut manager
- TClonesArray fFunctions; // functions
-
- ClassDef (AliRsnPair, 1)
+ void SetUp(EPairType type);
+ void SetAllFlags(AliRsnDaughter::EPIDMethod pid, Bool_t mix) {fPIDMethod = pid; fIsMixed = mix;}
+ void LoopPair(AliRsnEvent *ev1, TArrayI *a1, AliRsnEvent *ev2, TArrayI *a2);
+
+ Bool_t CutPass(AliRsnDaughter *d);
+ Bool_t CutPass(AliRsnPairParticle *p);
+ Bool_t CutPass(AliRsnEvent *e);
+
+ Bool_t fIsMixed; // doing event-mixing ?
+ EPairType fPairType; // pair type (PID + mixing or not)
+ AliRsnDaughter::EPIDMethod fPIDMethod; // pid type variable for single track
+
+ AliRsnPairDef *fPairDef; // pair definition (particles, charges)
+ AliRsnCutMgr *fCutMgr; // cut manager
+ TClonesArray fFunctions; // functions
+
+ ClassDef (AliRsnPair, 2)
};
#endif
// author: A. Pulvirenti (alberto.pulvirenti@ct.infn.it)
//
+#include <Riostream.h>
#include <TString.h>
#include "AliLog.h"
#include "AliESDEvent.h"
#include "AliESDtrack.h"
#include "AliESDVertex.h"
+#include "AliESDtrackCuts.h"
#include "AliAODEvent.h"
#include "AliAODTrack.h"
#include "AliRsnDaughter.h"
#include "AliRsnEvent.h"
#include "AliRsnPIDDefESD.h"
+#include "AliRsnCut.h"
#include "AliRsnReader.h"
fCheckSplit(kFALSE),
fRejectFakes(kFALSE),
fTPCOnly(kFALSE),
+ fUseESDTrackCuts(kFALSE),
+ fUseRsnTrackCuts(kFALSE),
fPIDDef(),
fITSClusters(0),
fTPCClusters(0),
fTRDClusters(0),
fTrackRefs(0),
fTrackRefsITS(0),
- fTrackRefsTPC(0)
+ fTrackRefsTPC(0),
+ fESDTrackCuts(),
+ fRsnTrackCuts("primaryCuts")
{
//
// Constructor.
}
}
+//_____________________________________________________________________________
+Bool_t AliRsnReader::ConvertTrack(AliRsnDaughter *daughter, AliESDtrack *esdTrack)
+{
+//
+// Translates an ESD track into a RSN daughter
+//
+
+ if (!esdTrack) {
+ AliError("Passed NULL object: nothing done");
+ return kFALSE;
+ }
+
+ Double_t p[3], v[3], pid[AliRsnPID::kSpecies];
+
+ // copy values which don't need treatment:
+ // momentum, vertex, chi2, flags, charge, number of TPC and ITS clusters
+ esdTrack->GetPxPyPz(p);
+ daughter->SetP(p[0], p[1], p[2]);
+ esdTrack->GetXYZ(v);
+ daughter->SetV(v[0], v[1], v[2]);
+ daughter->SetChi2(esdTrack->GetConstrainedChi2());
+ daughter->SetFlags(esdTrack->GetStatus());
+ daughter->SetCharge((Short_t)esdTrack->Charge());
+ daughter->SetNumberOfITSClusters(esdTrack->GetITSclusters(0x0));
+ daughter->SetNumberOfTPCClusters(esdTrack->GetTPCclusters(0x0));
+
+ // define the kink index:
+ // 0 = no kink
+ // 1 = kink daughter
+ // -1 = kink mother
+ Int_t i, ik[3];
+ for (i = 0; i < 3; i++) ik[i] = esdTrack->GetKinkIndex(i);
+ if (ik[0] < 0 || ik[1] < 0 || ik[2] < 0) daughter->SetKinkMother();
+ else if (ik[0] > 0 || ik[1] > 0 || ik[2] > 0) daughter->SetKinkDaughter();
+ else daughter->SetNoKink();
+
+ // store PID weights according to definition
+ // check: if the sum of all weights is null, the adoption fails
+ Double_t sum = 0.0;
+ fPIDDef.ComputeWeights(esdTrack, pid);
+ for (i = 0; i < AliRsnPID::kSpecies; i++)
+ {
+ daughter->SetPIDWeight(i, pid[i]);
+ sum += pid[i];
+ }
+ if (sum <= 0.0) return kFALSE;
+
+ // calculate N sigma to vertex
+ AliESDtrackCuts trkCut;
+ daughter->SetNSigmaToVertex(trkCut.GetSigmaToVertex(esdTrack));
+
+ return kTRUE;
+}
+
+//_____________________________________________________________________________
+Bool_t AliRsnReader::ConvertTrack(AliRsnDaughter *daughter, AliAODTrack *aodTrack)
+{
+//
+// Translates an AOD track into a RSN daughter
+//
+
+ if (!aodTrack)
+ {
+ AliError("Passed NULL object: nothing done");
+ return kFALSE;
+ }
+
+ // copy momentum and vertex
+ daughter->SetP(aodTrack->Px(), aodTrack->Py(), aodTrack->Pz());
+ daughter->SetV(aodTrack->Xv(), aodTrack->Yv(), aodTrack->Zv());
+
+ // chi2
+ daughter->SetChi2(aodTrack->Chi2perNDF());
+
+ // copy PID weights
+ Int_t i;
+ for (i = 0; i < 5; i++) daughter->SetPIDWeight(i, aodTrack->PID()[i]);
+
+ // copy flags
+ daughter->SetFlags(aodTrack->GetStatus());
+
+ // copy sign
+ daughter->SetCharge(aodTrack->Charge());
+
+ return kTRUE;
+}
+
+//_____________________________________________________________________________
+Bool_t AliRsnReader::ConvertTrack(AliRsnDaughter *daughter, TParticle *particle)
+{
+//
+// Translates an MC particle into a RSN daughter
+//
+
+ if (!particle) return kFALSE;
+
+ // copy other MC info (mother PDG code cannot be retrieved here)
+ daughter->InitMCInfo(particle);
+
+ // copy momentum and vertex
+ daughter->SetP(particle->Px(), particle->Py(), particle->Pz());
+ daughter->SetV(particle->Vx(), particle->Vy(), particle->Vz());
+
+ // recognize charge sign from PDG code sign
+ Int_t pdg = particle->GetPdgCode();
+ Int_t absPDG = TMath::Abs(pdg);
+ if (absPDG == 11 || absPDG == 13)
+ {
+ if (pdg > 0) daughter->SetCharge(-1); else daughter->SetCharge(1);
+ }
+ else if (absPDG == 211 || absPDG == 321 || absPDG == 2212)
+ {
+ if (pdg > 0) daughter->SetCharge(1); else daughter->SetCharge(-1);
+ }
+ else
+ {
+ // when trying to "adopt" a neutral track (photon, neutron, etc.)
+ // for the moment a "failed" message is returned
+ return kFALSE;
+ }
+
+ // flags and PID weights make no sense with MC tracks
+ daughter->SetFlags(0);
+ for (pdg = 0; pdg < AliRsnPID::kSpecies; pdg++) daughter->SetPIDWeight(pdg, 0.0);
+ daughter->SetPIDWeight(AliRsnPID::InternalType(absPDG), 1.0);
+
+ return kTRUE;
+}
+
//_____________________________________________________________________________
Bool_t AliRsnReader::Fill
(AliRsnEvent *rsn, AliVEvent *event, AliMCEvent *mc)
AliStack *stack = 0x0;
if (mc) stack = mc->Stack();
+ // set MC primary vertex if present
+ TArrayF fvertex(3);
+ Double_t mcvx = 0., mcvy = 0., mcvz = 0.;
+ if (mc) {
+ mc->GenEventHeader()->PrimaryVertex(fvertex);
+ mcvx = (Double_t)fvertex[0];
+ mcvy = (Double_t)fvertex[1];
+ mcvz = (Double_t)fvertex[2];
+ rsn->SetPrimaryVertexMC(mcvx, mcvy, mcvz);
+ }
+
// get number of tracks
Int_t ntracks = esd->GetNumberOfTracks();
if (!ntracks) return kFALSE;
}
// get ESD track
esdTrack = esd->GetTrack(index);
+ // check for ESD track cuts (if required)
+ if (fUseESDTrackCuts && (!fESDTrackCuts.AcceptTrack(esdTrack))) continue;
// check for fake tracks
label = esdTrack->GetLabel();
if (fRejectFakes && (label < 0)) continue;
}
// try to get information from this track
// output value tells if this was successful or not
- check = temp.Adopt(esdTrack, fPIDDef);
+ check = ConvertTrack(&temp, esdTrack);
if (!check) {
AliDebug(10, Form("Failed adopting track #%d", index));
continue;
temp.SetIndex(index);
temp.SetLabel(label);
- // shifts the track DCA to the found vertex
- temp.ShiftZero(vertex[0], vertex[1], vertex[2]);
+ // check this object against the Rsn cuts (if required)
+ if (fUseRsnTrackCuts) {
+ if (!fRsnTrackCuts.IsSelected(AliRsnCut::kParticle, &temp)) continue;
+ }
// try to add track to collection and returns an error in case of problems
AliRsnDaughter *ptr = rsn->AddTrack(temp);
}
// sort tracks w.r. to Pt (from largest to smallest)
- rsn->SortTracks();
+ //rsn->SortTracks();
+
+ // correct tracks for primary vertex
+ //rsn->CorrectTracks();
+
return kTRUE;
}
vertex[2] = aod->GetPrimaryVertex()->GetZ();
rsn->SetPrimaryVertex(vertex[0], vertex[1], vertex[2]);
+ // set MC primary vertex if present
+ TArrayF fvertex(3);
+ Double_t mcvx = 0., mcvy = 0., mcvz = 0.;
+ if (mc) {
+ mc->GenEventHeader()->PrimaryVertex(fvertex);
+ mcvx = (Double_t)fvertex[0];
+ mcvy = (Double_t)fvertex[1];
+ mcvz = (Double_t)fvertex[2];
+ rsn->SetPrimaryVertexMC(mcvx, mcvy, mcvz);
+ }
+
// store tracks from ESD
Int_t index, label, labmum;
Bool_t check;
if (fRejectFakes && (label < 0)) continue;
// copy ESD track data into RsnDaughter
// if unsuccessful, this track is skipped
- check = temp.Adopt(aodTrack);
+ check = ConvertTrack(&temp, aodTrack);
if (!check) continue;
// if stack is present, copy MC info
if (stack)
// set index and label and add this object to the output container
temp.SetIndex(index);
temp.SetLabel(label);
+
+ // check this object against the Rsn cuts (if required)
+ if (fUseRsnTrackCuts) {
+ if (!fRsnTrackCuts.IsSelected(AliRsnCut::kParticle, &temp)) continue;
+ }
+
AliRsnDaughter *ptr = rsn->AddTrack(temp);
// if problems occurred while storin, that pointer is NULL
if (!ptr) AliWarning(Form("Failed storing track#%d"));
return kFALSE;
}
+ // correct tracks for primary vertex
+ rsn->CorrectTracks();
+
return kTRUE;
}
// get primary vertex
TArrayF fvertex(3);
- Double_t vertex[3];
+ Double_t vx, vy, vz;
mc->GenEventHeader()->PrimaryVertex(fvertex);
- vertex[0] = (Double_t)fvertex[0];
- vertex[1] = (Double_t)fvertex[1];
- vertex[2] = (Double_t)fvertex[2];
- rsn->SetPrimaryVertex(vertex[0], vertex[1], vertex[2]);
+ vx = (Double_t)fvertex[0];
+ vy = (Double_t)fvertex[1];
+ vz = (Double_t)fvertex[2];
+ rsn->SetPrimaryVertex(vx, vy, vz);
+ rsn->SetPrimaryVertexMC(vx, vy, vz);
// store tracks from MC
- Int_t i, index, labmum, nHitsITS, nHitsTPC, nRef;
- Bool_t check;
+ Int_t i, index, labmum, nRef, nRefITS, nRefTPC, detectorID;
AliRsnDaughter temp;
- for (index = 0; index < ntracks; index++) {
- // get and check MC track
+ for (index = 0; index < ntracks; index++)
+ {
+ // get MC track & take index and label
AliMCParticle *mcTrack = mc->GetTrack(index);
+ temp.SetIndex(index);
+ temp.SetLabel(mcTrack->Label());
+
// check particle track references
nRef = mcTrack->GetNumberOfTrackReferences();
- if (fTrackRefs > 0 && nRef < fTrackRefs) continue;
- else if (fTrackRefsITS > 0 || fTrackRefsTPC > 0) {
- nHitsITS = nHitsTPC = 0;
- for (i = 0; i < nRef; i++) {
- AliTrackReference *trackRef = mcTrack->GetTrackReference(i);
- if (trackRef) {
- Int_t detectorId = trackRef->DetectorId();
- switch(detectorId) {
- case AliTrackReference::kITS : nHitsITS++ ; break ;
- case AliTrackReference::kTPC : nHitsTPC++ ; break ;
- default : break ;
- }
- }
+ for (i = 0, nRefITS = 0, nRefTPC = 0; i < nRef; i++) {
+ AliTrackReference *trackRef = mcTrack->GetTrackReference(i);
+ if (!trackRef) continue;
+ detectorID = trackRef->DetectorId();
+ switch (detectorID) {
+ case AliTrackReference::kITS : nRefITS++; break;
+ case AliTrackReference::kTPC : nRefTPC++; break;
+ default: break;
}
- if (fTrackRefsITS > 0 && nHitsITS < fTrackRefsITS) continue;
- if (fTrackRefsTPC > 0 && nHitsTPC < fTrackRefsTPC) continue;
}
+ if (fTrackRefs > 0 && nRef < fTrackRefs) continue;
+ if (fTrackRefsITS > 0 && nRefITS < fTrackRefsITS) continue;
+ if (fTrackRefsTPC > 0 && nRefTPC < fTrackRefsTPC) continue;
+
// try to insert in the RsnDaughter its data
- check = temp.Adopt(mcTrack);
- if (!check) continue;
+ TParticle *mcPart = mcTrack->Particle();
+ if (!ConvertTrack(&temp, mcPart)) continue;
+
+ // assign mother label and PDG code (if any)
labmum = temp.GetMCInfo()->Mother();
- if (labmum >= 0)
- {
+ if (labmum >= 0) {
TParticle *mum = stack->Particle(labmum);
temp.GetMCInfo()->SetMotherPDG(mum->GetPdgCode());
}
- // if successful, set other data and stores it
- temp.SetIndex(index);
- temp.SetLabel(mcTrack->Label());
+
+ // check this object against the Rsn cuts (if required)
+ if (fUseRsnTrackCuts) {
+ if (!fRsnTrackCuts.IsSelected(AliRsnCut::kParticle, &temp)) continue;
+ }
+
AliRsnDaughter *ptr = rsn->AddTrack(temp);
- // if problems occurred while storin, that pointer is NULL
- if (!ptr) AliWarning(Form("Failed storing track#%d", index));
+ if (!ptr) AliWarning(Form("Failed storing track #%d", index));
}
// compute total multiplicity
// sort tracks w.r. to Pt (from largest to smallest)
rsn->SortTracks();
+
+ // correct tracks for primary vertex
+ rsn->CorrectTracks();
+
return kTRUE;
}
#include <TNamed.h>
+#include "AliESDtrackCuts.h"
+
#include "AliRsnDaughter.h"
#include "AliRsnPIDDefESD.h"
+#include "AliRsnCutSet.h"
class AliVEvent;
class AliESDEvent;
public:
AliRsnReader();
- virtual ~AliRsnReader() {}
+ virtual ~AliRsnReader() { }
void SetCheckSplit(Bool_t doit = kTRUE) {fCheckSplit = doit;}
Bool_t AreSplitted(AliESDtrack *track1, AliESDtrack *track2);
void SetTPCOnly(Bool_t doit = kTRUE);
Bool_t DoesTPCOnly() {return fTPCOnly;}
- AliRsnPIDDefESD& GetPIDDef() {return fPIDDef;}
+ void SetUseESDTrackCuts(Bool_t doit = kTRUE) {fUseESDTrackCuts = doit;}
+ Bool_t DoesESDTrackCuts() {return fUseESDTrackCuts;}
+ AliESDtrackCuts* GetESDTrackCuts() {return &fESDTrackCuts;}
+
+ void SetUseRsnTrackCuts(Bool_t doit = kTRUE) {fUseRsnTrackCuts = doit;}
+ Bool_t DoesRsnTrackCuts() {return fUseRsnTrackCuts;}
+ AliRsnCutSet* GetRsnTrackCuts() {return &fRsnTrackCuts;}
+
+ AliRsnPIDDefESD* GetPIDDef() {return &fPIDDef;}
void SetTrackRefs(Int_t value) {fTrackRefs = value;}
void SetTrackRefsITS(Int_t value) {fTrackRefsITS = value;}
void SetMinTRDClusters(Int_t value) {fTRDClusters = value;}
void SetITSTPCTRDSectors(const Int_t& its = -1, const Int_t& tpc = -1, const Int_t& trd = -1);
+ Bool_t ConvertTrack(AliRsnDaughter *daughter, AliESDtrack *track);
+ Bool_t ConvertTrack(AliRsnDaughter *daughter, AliAODTrack *track);
+ Bool_t ConvertTrack(AliRsnDaughter *daughter, TParticle *particle);
+
Bool_t Fill(AliRsnEvent *rsn, AliVEvent *event, AliMCEvent *refMC = 0);
Bool_t FillFromESD(AliRsnEvent *rsn, AliESDEvent *event, AliMCEvent *refMC = 0);
Bool_t FillFromAOD(AliRsnEvent *rsn, AliAODEvent *event, AliMCEvent *refMC = 0);
// dummy copy methods
AliRsnReader(const AliRsnReader ©) :
- TObject(copy),fCheckSplit(0),fRejectFakes(0),fTPCOnly(0),fPIDDef(copy.fPIDDef),
- fITSClusters(0),fTPCClusters(0),fTRDClusters(0),
- fTrackRefs(0),fTrackRefsITS(0),fTrackRefsTPC(0) { /*nothing*/ }
+ TObject(copy),fCheckSplit(0),fRejectFakes(0),fTPCOnly(0),fUseESDTrackCuts(0),fUseRsnTrackCuts(0),
+ fPIDDef(copy.fPIDDef),fITSClusters(0),fTPCClusters(0),fTRDClusters(0),
+ fTrackRefs(0),fTrackRefsITS(0),fTrackRefsTPC(0),fESDTrackCuts(),fRsnTrackCuts("") { /*nothing*/ }
AliRsnReader& operator=(const AliRsnReader&) {return (*this);}
- Bool_t fCheckSplit; // flag to check and remove split tracks
- Bool_t fRejectFakes; // flag to reject fake tracks (negative label)
- Bool_t fTPCOnly; // flag to use only the TPC for reading data
+ Bool_t fCheckSplit; // flag to check and remove split tracks
+ Bool_t fRejectFakes; // flag to reject fake tracks (negative label)
+ Bool_t fTPCOnly; // flag to use only the TPC for reading data
+ Bool_t fUseESDTrackCuts; // flag to use ESD track cuts
+ Bool_t fUseRsnTrackCuts; // flag to use ESD track cuts
+
+ AliRsnPIDDefESD fPIDDef; // manager for alternative PID weights (ESD only)
- AliRsnPIDDefESD fPIDDef; // manager for alternative PID weights (ESD only)
+ Int_t fITSClusters; // minimum number of ITS clusters to accept a track
+ Int_t fTPCClusters; // minimum number of TPC clusters to accept a track
+ Int_t fTRDClusters; // minimum number of TRD clusters to accept a track
- Int_t fITSClusters; // minimum number of ITS clusters to accept a track
- Int_t fTPCClusters; // minimum number of TPC clusters to accept a track
- Int_t fTRDClusters; // minimum number of TRD clusters to accept a track
+ Int_t fTrackRefs; // minimum required track references for MC reading
+ Int_t fTrackRefsITS; // minimum required track references for MC reading (ITS)
+ Int_t fTrackRefsTPC; // minimum required track references for MC reading (TPC)
- Int_t fTrackRefs; // minimum required track references for MC reading
- Int_t fTrackRefsITS; // minimum required track references for MC reading (ITS)
- Int_t fTrackRefsTPC; // minimum required track references for MC reading (TPC)
+ AliESDtrackCuts fESDTrackCuts; // object for ESD track cuts
+ AliRsnCutSet fRsnTrackCuts; // other local cuts used in preliminary track selection
ClassDef(AliRsnReader, 1);
};
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff