* 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. *
- **************************************************************************/
+ **********************************************************i****************/
/* $Id$ */
// Data selection for flow framework
//
// origin: Mikolaj Krzewicki (mikolaj.krzewicki@cern.ch)
+// mods: Redmer A. Bertens (rbertens@cern.ch)
//
// This class gurantees consistency of cut methods, trackparameter
// selection (global tracks, TPC only, etc..) and parameter mixing
#include "AliESDEvent.h"
#include "AliAODEvent.h"
#include "AliVParticle.h"
+#include "AliVVZERO.h"
#include "AliMCParticle.h"
#include "AliESDtrack.h"
+#include "AliESDMuonTrack.h" // XZhang 20120604
#include "AliMultiplicity.h"
#include "AliAODTrack.h"
+#include "AliAODTracklets.h" // XZhang 20120615
#include "AliFlowTrackSimple.h"
#include "AliFlowTrack.h"
#include "AliFlowTrackCuts.h"
AliFlowTrackCuts::AliFlowTrackCuts():
AliFlowTrackSimpleCuts(),
fAliESDtrackCuts(NULL),
+ fMuonTrackCuts(NULL), // XZhang 20120604
fQA(NULL),
fCutMC(kFALSE),
fCutMChasTrackReferences(kFALSE),
fRequireStrictTOFTPCagreement(kFALSE),
fCutRejectElectronsWithTPCpid(kFALSE),
fProbBayes(0.0),
- fCurrCentr(0.0)
+ fCurrCentr(0.0),
+ fV0gainEqualization(NULL),
+ fApplyRecentering(kFALSE),
+ fV0gainEqualizationPerRing(kFALSE)
{
//io constructor
SetPriors(); //init arrays
+
+ // New PID procedure (Bayesian Combined PID)
+ // allocating here is necessary because we don't
+ // stream this member
+ // TODO: fix streaming problems AliFlowBayesianPID
+ fBayesianResponse = new AliFlowBayesianPID();
+ fBayesianResponse->SetNewTrackParam();
+ for(Int_t i(0); i < 4; i++) {
+ fV0Apol[i] = 0;
+ fV0Cpol[i] = 0;
+ }
+ for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = kTRUE;
+
}
//-----------------------------------------------------------------------
AliFlowTrackCuts::AliFlowTrackCuts(const char* name):
AliFlowTrackSimpleCuts(name),
fAliESDtrackCuts(NULL),
+ fMuonTrackCuts(NULL), // XZhang 20120604
fQA(NULL),
fCutMC(kFALSE),
fCutMChasTrackReferences(kFALSE),
fRequireStrictTOFTPCagreement(kFALSE),
fCutRejectElectronsWithTPCpid(kFALSE),
fProbBayes(0.0),
- fCurrCentr(0.0)
+ fCurrCentr(0.0),
+ fV0gainEqualization(NULL),
+ fApplyRecentering(kFALSE),
+ fV0gainEqualizationPerRing(kFALSE)
{
//constructor
SetTitle("AliFlowTrackCuts");
// New PID procedure (Bayesian Combined PID)
fBayesianResponse = new AliFlowBayesianPID();
fBayesianResponse->SetNewTrackParam();
+ for(Int_t i(0); i < 4; i++) {
+ fV0Apol[i] = 0;
+ fV0Cpol[i] = 0;
+ }
+ for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = kTRUE;
}
//-----------------------------------------------------------------------
AliFlowTrackCuts::AliFlowTrackCuts(const AliFlowTrackCuts& that):
AliFlowTrackSimpleCuts(that),
fAliESDtrackCuts(NULL),
+ fMuonTrackCuts(NULL), // XZhang 20120604
fQA(NULL),
fCutMC(that.fCutMC),
fCutMChasTrackReferences(that.fCutMChasTrackReferences),
fRequireStrictTOFTPCagreement(that.fRequireStrictTOFTPCagreement),
fCutRejectElectronsWithTPCpid(that.fCutRejectElectronsWithTPCpid),
fProbBayes(0.0),
- fCurrCentr(0.0)
+ fCurrCentr(0.0),
+ fV0gainEqualization(NULL),
+ fApplyRecentering(that.fApplyRecentering),
+ fV0gainEqualizationPerRing(that.fV0gainEqualizationPerRing)
{
//copy constructor
+ printf(" \n\n claling copy ctor \n\n" );
if (that.fTPCpidCuts) fTPCpidCuts = new TMatrixF(*(that.fTPCpidCuts));
if (that.fTOFpidCuts) fTOFpidCuts = new TMatrixF(*(that.fTOFpidCuts));
if (that.fAliESDtrackCuts) fAliESDtrackCuts = new AliESDtrackCuts(*(that.fAliESDtrackCuts));
+ if (that.fMuonTrackCuts) fMuonTrackCuts = new AliMuonTrackCuts(*(that.fMuonTrackCuts)); // XZhang 20120604
SetPriors(); //init arrays
if (that.fQA) DefineHistograms();
// New PID procedure (Bayesian Combined PID)
fBayesianResponse = new AliFlowBayesianPID();
fBayesianResponse->SetNewTrackParam();
+
+ // V0 gain calibration
+ // no reason to init fV0gainEqualizationPerRing, will be initialized on node if necessary
+ // pointer is set to NULL in initialization list of this constructor
+// if (that.fV0gainEqualization) fV0gainEqualization = new TH1(*(that.fV0gainEqualization));
+ for(Int_t i(0); i < 4; i++) { // no use to copy these guys since they're only initialized on worked node
+ fV0Apol[i] = 0.;
+ fV0Cpol[i] = 0.;
+ }
+ for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = that.fUseVZERORing[i];
+
}
//-----------------------------------------------------------------------
if (that.fAliESDtrackCuts && fAliESDtrackCuts) *fAliESDtrackCuts=*(that.fAliESDtrackCuts);
if (that.fAliESDtrackCuts && !fAliESDtrackCuts) fAliESDtrackCuts=new AliESDtrackCuts(*(that.fAliESDtrackCuts));
if (!that.fAliESDtrackCuts) delete fAliESDtrackCuts; fAliESDtrackCuts=NULL;
+
+ if ( that.fMuonTrackCuts && fMuonTrackCuts) *fMuonTrackCuts = *(that.fMuonTrackCuts); // XZhang 20120604
+ if ( that.fMuonTrackCuts && !fMuonTrackCuts) fMuonTrackCuts = new AliMuonTrackCuts(*(that.fMuonTrackCuts)); // XZhang 20120604
+ if (!that.fMuonTrackCuts) delete fMuonTrackCuts; fMuonTrackCuts = NULL; // XZhang 20120604
+ if (!that.fV0gainEqualization) delete fV0gainEqualization; fV0gainEqualization = NULL;
//these guys we don't need to copy, just reinit
if (that.fQA) {fQA->Delete(); delete fQA; fQA=NULL; DefineHistograms();}
fCutMC=that.fCutMC;
fCutRejectElectronsWithTPCpid=that.fCutRejectElectronsWithTPCpid;
fProbBayes = that.fProbBayes;
fCurrCentr = that.fCurrCentr;
-
- // New PID procedure (Bayesian Combined PID)
- fBayesianResponse = new AliFlowBayesianPID();
- fBayesianResponse->SetNewTrackParam();
+
+ fApplyRecentering = that.fApplyRecentering;
+ fV0gainEqualizationPerRing = that.fV0gainEqualizationPerRing;
+#if ROOT_VERSION_CODE < ROOT_VERSION(5,99,0)
+ if (that.fV0gainEqualization) fV0gainEqualization = new TH1(*(that.fV0gainEqualization));
+#else
+ //PH Lets try Clone, however the result might be wrong
+ if (that.fV0gainEqualization) fV0gainEqualization = (TH1*)that.fV0gainEqualization->Clone();
+#endif
+ for(Int_t i(0); i < 4; i++) { // no use to copy these guys since they're only initialized on worked node
+ fV0Apol[i] = that.fV0Apol[i];
+ fV0Cpol[i] = that.fV0Cpol[i];
+ }
+ for(Int_t i(0); i < 8; i++) fUseVZERORing[i] = that.fUseVZERORing[i];
return *this;
}
delete fAliESDtrackCuts;
delete fTPCpidCuts;
delete fTOFpidCuts;
+ if (fMuonTrackCuts) delete fMuonTrackCuts; // XZhang 20120604
if (fQA) { fQA->SetOwner(); fQA->Delete(); delete fQA; }
+ if (fV0gainEqualization) delete fV0gainEqualization;
}
//-----------------------------------------------------------------------
//do the magic for ESD
AliESDEvent* myESD = dynamic_cast<AliESDEvent*>(event);
+ AliAODEvent* myAOD = dynamic_cast<AliAODEvent*>(event);
if (fCutPID && myESD)
{
//TODO: maybe call it only for the TOF options?
fESDpid.SetTOFResponse(myESD,AliESDpid::kTOF_T0);
// End F. Noferini added part
}
-
- //TODO: AOD
+ if (fCutPID && myAOD){
+ fBayesianResponse->SetDetResponse(myAOD, fCurrCentr,AliESDpid::kTOF_T0); // centrality = PbPb centrality class (0-100%) or -1 for pp collisions
+ if(myAOD->GetTOFHeader()){
+ fESDpid.SetTOFResponse(myAOD,AliESDpid::kTOF_T0);
+ }
+ else{ // corrected on the fly track by track if tof header is not present (old AODs)
+ }
+ // End F. Noferini added part
+ }
+
+ if(fPIDsource==kTOFbayesian) fBayesianResponse->SetDetAND(1);
+ else if(fPIDsource==kTPCbayesian) fBayesianResponse->ResetDetOR(1);
}
//-----------------------------------------------------------------------
{
//check cuts
AliVParticle* vparticle = dynamic_cast<AliVParticle*>(obj);
- if (vparticle) return PassesCuts(vparticle);
+//if (vparticle) return PassesCuts(vparticle); // XZhang 20120604
+ if (vparticle) { // XZhang 20120604
+ if (fParamType==kMUON) return PassesMuonCuts(vparticle); // XZhang 20120604
+ return PassesCuts(vparticle); // XZhang 20120604
+ } // XZhang 20120604
+
AliFlowTrackSimple* flowtrack = dynamic_cast<AliFlowTrackSimple*>(obj);
if (flowtrack) return PassesCuts(flowtrack);
AliMultiplicity* tracklets = dynamic_cast<AliMultiplicity*>(obj);
if (tracklets) return PassesCuts(tracklets,id);
+ AliAODTracklets* trkletAOD = dynamic_cast<AliAODTracklets*>(obj); // XZhang 20120615
+ if (trkletAOD) return PassesCuts(trkletAOD,id); // XZhang 20120615
AliESDPmdTrack* pmdtrack = dynamic_cast<AliESDPmdTrack*>(obj);
if (pmdtrack) return PassesPMDcuts(pmdtrack);
AliVEvent* vvzero = dynamic_cast<AliVEvent*>(obj); // should be removed; left for protection only
return kTRUE;
}
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesCuts(const AliAODTracklets* tracklet, Int_t id)
+{
+ // XZhang 20120615
+ //check cuts on a tracklets
+
+ if (id<0) return kFALSE;
+
+ //clean up from last iteration, and init label
+ fTrack = NULL;
+ fMCparticle=NULL;
+ fTrackLabel=-999;
+
+ fTrackPhi = tracklet->GetPhi(id);
+//fTrackEta = tracklet->GetEta(id);
+ fTrackEta = -1.*TMath::Log(TMath::Tan(tracklet->GetTheta(id)/2.));
+ fTrackWeight = 1.0;
+ if (fCutEta) {if ( fTrackEta < fEtaMin || fTrackEta >= fEtaMax ) return kFALSE;}
+ if (fCutPhi) {if ( fTrackPhi < fPhiMin || fTrackPhi >= fPhiMax ) return kFALSE;}
+
+ //check MC info if available
+ //if the 2 clusters have different label track cannot be good
+ //and should therefore not pass the mc cuts
+ Int_t label0 = tracklet->GetLabel(id,0);
+ Int_t label1 = tracklet->GetLabel(id,1);
+ //if possible get label and mcparticle
+ fTrackLabel = (label0==label1)?tracklet->GetLabel(id,1):-1;
+ if (!fFakesAreOK && fTrackLabel<0) return kFALSE;
+ if (fTrackLabel>=0 && fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
+ //check MC cuts
+ if (fCutMC && !PassesMCcuts()) return kFALSE;
+ return kTRUE;
+}
+
//-----------------------------------------------------------------------
Bool_t AliFlowTrackCuts::PassesMCcuts(AliMCEvent* mcEvent, Int_t label)
{
//the case of ESD or AOD
if (esdTrack) { if (!PassesESDcuts(esdTrack)) { pass=kFALSE; } }
- if (aodTrack) { if (!PassesAODcuts(aodTrack)) { pass=kFALSE; } }
+ if (aodTrack) { if (!PassesAODcuts(aodTrack,pass)) { pass=kFALSE; } }
if (fQA)
{
}
//_______________________________________________________________________
-Bool_t AliFlowTrackCuts::PassesAODcuts(const AliAODTrack* track)
+Bool_t AliFlowTrackCuts::PassesAODcuts(const AliAODTrack* track, Bool_t passedFid)
{
//check cuts for AOD
- Bool_t pass = kTRUE;
+ Bool_t pass = passedFid;
if (fCutNClustersTPC)
{
if (fCutDCAToVertexZ && track->ZAtDCA()>GetMaxDCAToVertexZ()) pass=kFALSE;
- if (track->GetTPCsignal() < fMinimalTPCdedx) pass=kFALSE;
+ Double_t dedx = track->GetTPCsignal();
+ if (dedx < fMinimalTPCdedx) pass=kFALSE;
+ Double_t time[5];
+ track->GetIntegratedTimes(time);
+ if (fQA) {
+ Double_t momTPC = track->GetTPCmomentum();
+ QAbefore( 1)->Fill(momTPC,dedx);
+ QAbefore( 5)->Fill(track->Pt(),track->DCA());
+ QAbefore( 6)->Fill(track->Pt(),track->ZAtDCA());
+ if (pass) QAafter( 1)->Fill(momTPC,dedx);
+ if (pass) QAafter( 5)->Fill(track->Pt(),track->DCA());
+ if (pass) QAafter( 6)->Fill(track->Pt(),track->ZAtDCA());
+ QAbefore( 8)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kElectron]));
+ if (pass) QAafter( 8)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kElectron]));
+ QAbefore( 9)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kMuon]));
+ if (pass) QAafter( 9)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kMuon]));
+ QAbefore(10)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kPion]));
+ if (pass) QAafter( 10)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kPion]));
+ QAbefore(11)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kKaon]));
+ if (pass) QAafter( 11)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kKaon]));
+ QAbefore(12)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kProton]));
+ if (pass) QAafter( 12)->Fill(track->P(),(track->GetTOFsignal()-time[AliPID::kProton]));
+ }
+ if (fCutPID && (fParticleID!=AliPID::kUnknown)) //if kUnknown don't cut on PID
+ {
+ if (!PassesAODpidCut(track)) pass=kFALSE;
+ }
+
return pass;
}
//-----------------------------------------------------------------------
AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts()
+{
+ //returns the lhc10h vzero track cuts, this function
+ //is left here for backward compatibility
+ //if a run is recognized as 11h, the calibration method will
+ //switch to 11h calbiration, which means that the cut
+ //object is updated but not replaced.
+ //calibratin is only available for PbPb runs
+ return GetStandardVZEROOnlyTrackCuts2010();
+}
+//-----------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts2010()
{
//get standard V0 cuts
- AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard vzero flow cuts");
+ AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard vzero flow cuts 2010");
cuts->SetParamType(kV0);
cuts->SetEtaRange( -10, +10 );
cuts->SetPhiMin( 0 );
cuts->SetPhiMax( TMath::TwoPi() );
+ // options for the reweighting
+ cuts->SetV0gainEqualizationPerRing(kFALSE);
+ cuts->SetApplyRecentering(kTRUE);
+ // to exclude a ring , do e.g.
+ // cuts->SetUseVZERORing(7, kFALSE);
+ return cuts;
+}
+//-----------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardVZEROOnlyTrackCuts2011()
+{
+ //get standard V0 cuts for 2011 data
+ //in this case, the vzero segments will be weighted by
+ //VZEROEqMultiplicity,
+ //if recentering is enableded, the sub-q vectors
+ //will be taken from the event header, so make sure to run
+ //the VZERO event plane selection task before this task !
+ //recentering replaces the already evaluated q-vectors, so
+ //when chosen, additional settings (e.g. excluding rings)
+ //have no effect. recentering is true by default
+ //
+ //NOTE user is responsible for running the vzero event plane
+ //selection task in advance, e.g. add to your launcher macro
+ //
+ // gROOT->LoadMacro("$ALICE_ROOT/ANALYSIS/macros/AddTaskVZEROEPSelection.C");
+ // AddTaskVZEROEPSelection();
+ //
+ AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard vzero flow cuts 2011");
+ cuts->SetParamType(kV0);
+ cuts->SetEtaRange( -10, +10 );
+ cuts->SetPhiMin( 0 );
+ cuts->SetPhiMax( TMath::TwoPi() );
+ cuts->SetApplyRecentering(kTRUE);
return cuts;
}
-
//-----------------------------------------------------------------------
AliFlowTrackCuts* AliFlowTrackCuts::GetStandardGlobalTrackCuts2010()
{
return cuts;
}
+//-----------------------------------------------------------------------------
+AliFlowTrackCuts* AliFlowTrackCuts::GetStandardMuonTrackCuts(Bool_t isMC, Int_t passN)
+{
+// XZhang 20120604
+ AliFlowTrackCuts* cuts = new AliFlowTrackCuts("standard muon track cuts");
+ cuts->SetParamType(kMUON);
+ cuts->SetStandardMuonTrackCuts();
+ cuts->SetIsMuonMC(isMC);
+ cuts->SetMuonPassNumber(passN);
+ return cuts;
+}
+
//-----------------------------------------------------------------------
Bool_t AliFlowTrackCuts::FillFlowTrackGeneric(AliFlowTrack* flowtrack) const
{
flowtrack->SetSource(AliFlowTrack::kFromESD);
flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
}
- else if (dynamic_cast<AliAODTrack*>(fTrack))
+ else if (dynamic_cast<AliESDMuonTrack*>(fTrack)) // XZhang 20120604
+ { // XZhang 20120604
+ flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
+ flowtrack->SetID((Int_t)static_cast<AliESDMuonTrack*>(fTrack)->GetUniqueID()); // XZhang 20120604
+ } // XZhang 20120604
+ else if (dynamic_cast<AliAODTrack*>(fTrack))
{
- flowtrack->SetSource(AliFlowTrack::kFromAOD);
+ if (fParamType==kMUON) // XZhang 20120604
+ flowtrack->SetSource(AliFlowTrack::kFromMUON); // XZhang 20120604
+ else // XZhang 20120604
+ flowtrack->SetSource(AliFlowTrack::kFromAOD); // XZhang 20120604
flowtrack->SetID(static_cast<AliVTrack*>(fTrack)->GetID());
}
else if (dynamic_cast<AliMCParticle*>(fTrack))
before->Add(hb);//7
after->Add(ha);//7
+ before->Add(new TH2F("TOFkElectron",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//8
+ after->Add( new TH2F("TOFkElectron",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//8
+
+ before->Add(new TH2F("TOFkMuon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//9
+ after->Add( new TH2F("TOFkMuon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//9
+
+ before->Add(new TH2F("TOFkPion",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//10
+ after->Add( new TH2F("TOFkPion",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//10
+
+ before->Add(new TH2F("TOFkKaon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//11
+ after->Add( new TH2F("TOFkKaon",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//11
+
+ before->Add(new TH2F("TOFkProton",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//12
+ after->Add( new TH2F("TOFkProton",";p_{t}[GeV/c];TOF signal - IT", kNbinsP,binsP,1000,-2e4, 2e4));//12
+
TH1::AddDirectory(adddirstatus);
}
//get the number of tracks in the input event according source
//selection (ESD tracks, tracklets, MC particles etc.)
AliESDEvent* esd=NULL;
+ AliAODEvent* aod=NULL; // XZhang 20120615
switch (fParamType)
{
case kSPDtracklet:
+ if (!fEvent) return 0; // XZhang 20120615
esd = dynamic_cast<AliESDEvent*>(fEvent);
- if (!esd) return 0;
- return esd->GetMultiplicity()->GetNumberOfTracklets();
+ aod = dynamic_cast<AliAODEvent*>(fEvent); // XZhang 20120615
+// if (!esd) return 0; // XZhang 20120615
+// return esd->GetMultiplicity()->GetNumberOfTracklets(); // XZhang 20120615
+ if (esd) return esd->GetMultiplicity()->GetNumberOfTracklets(); // XZhang 20120615
+ if (aod) return aod->GetTracklets()->GetNumberOfTracklets(); // XZhang 20120615
case kMC:
if (!fMCevent) return 0;
return fMCevent->GetNumberOfTracks();
return esd->GetNumberOfPmdTracks();
case kV0:
return fgkNumberOfV0tracks;
+ case kMUON: // XZhang 20120604
+ if (!fEvent) return 0; // XZhang 20120604
+ esd = dynamic_cast<AliESDEvent*>(fEvent); // XZhang 20120604
+ if (esd) return esd->GetNumberOfMuonTracks(); // XZhang 20120604
+ return fEvent->GetNumberOfTracks(); // if AOD // XZhang 20120604
default:
if (!fEvent) return 0;
return fEvent->GetNumberOfTracks();
//get the input object according the data source selection:
//(esd tracks, traclets, mc particles,etc...)
AliESDEvent* esd=NULL;
+ AliAODEvent* aod=NULL; // XZhang 20120615
switch (fParamType)
{
case kSPDtracklet:
+ if (!fEvent) return NULL; // XZhang 20120615
esd = dynamic_cast<AliESDEvent*>(fEvent);
- if (!esd) return NULL;
- return const_cast<AliMultiplicity*>(esd->GetMultiplicity());
+ aod = dynamic_cast<AliAODEvent*>(fEvent); // XZhang 20120615
+// if (!esd) return NULL; // XZhang 20120615
+// return const_cast<AliMultiplicity*>(esd->GetMultiplicity()); // XZhang 20120615
+ if (esd) return const_cast<AliMultiplicity*>(esd->GetMultiplicity()); // XZhang 20120615
+ if (aod) return const_cast<AliAODTracklets*>(aod->GetTracklets()); // XZhang 20120615
case kMC:
if (!fMCevent) return NULL;
return fMCevent->GetTrack(i);
//}
//return esd->GetVZEROData();
return fEvent; // left only for compatibility
+ case kMUON: // XZhang 20120604
+ if (!fEvent) return NULL; // XZhang 20120604
+ esd = dynamic_cast<AliESDEvent*>(fEvent); // XZhang 20120604
+ if (esd) return esd->GetMuonTrack(i); // XZhang 20120604
+ return fEvent->GetTrack(i); // if AOD // XZhang 20120604
default:
if (!fEvent) return NULL;
return fEvent->GetTrack(i);
fTrackEta=0.0;
fTrackPhi=0.0;
}
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesAODpidCut(const AliAODTrack* track )
+{
+ if(!track->GetAODEvent()->GetTOFHeader()){
+ AliAODPid *pidObj = track->GetDetPid();
+ if (!pidObj) fESDpid.GetTOFResponse().SetTimeResolution(84.);
+ else{
+ Double_t sigmaTOFPidInAOD[10];
+ pidObj->GetTOFpidResolution(sigmaTOFPidInAOD);
+ if(sigmaTOFPidInAOD[0] > 84.){
+ fESDpid.GetTOFResponse().SetTimeResolution(sigmaTOFPidInAOD[0]); // use the electron TOF PID sigma as time resolution (including the T0 used)
+ }
+ }
+ }
+
+ //check if passes the selected pid cut for ESDs
+ Bool_t pass = kTRUE;
+ switch (fPIDsource)
+ {
+ case kTOFbeta:
+ if (!PassesTOFbetaCut(track)) pass=kFALSE;
+ break;
+ case kTOFbayesian:
+ if (!PassesTOFbayesianCut(track)) pass=kFALSE;
+ break;
+ case kTPCbayesian:
+ if (!PassesTPCbayesianCut(track)) pass=kFALSE;
+ break;
+ default:
+ return kTRUE;
+ break;
+ }
+ return pass;
+}
//-----------------------------------------------------------------------
Bool_t AliFlowTrackCuts::PassesESDpidCut(const AliESDtrack* track )
{
}
//-----------------------------------------------------------------------
-Float_t AliFlowTrackCuts::GetBeta(const AliESDtrack* track)
+Float_t AliFlowTrackCuts::GetBeta(const AliVTrack* track)
{
//get beta
+ Double_t integratedTimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
+ track->GetIntegratedTimes(integratedTimes);
+
const Float_t c = 2.99792457999999984e-02;
- Float_t p = track->GetP();
- Float_t l = track->GetIntegratedLength();
+ Float_t p = track->P();
+ Float_t l = integratedTimes[0]*c;
Float_t trackT0 = fESDpid.GetTOFResponse().GetStartTime(p);
Float_t timeTOF = track->GetTOFsignal()- trackT0;
return l/timeTOF/c;
}
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesTOFbetaCut(const AliAODTrack* track )
+{
+ //check if track passes pid selection with an asymmetric TOF beta cut
+ if (!fTOFpidCuts)
+ {
+ //printf("no TOFpidCuts\n");
+ return kFALSE;
+ }
+
+ //check if passes PID cut using timing in TOF
+ Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFpid) &&
+ (track->GetTOFsignal() > 12000) &&
+ (track->GetTOFsignal() < 100000);
+
+ if (!goodtrack) return kFALSE;
+
+ const Float_t c = 2.99792457999999984e-02;
+ Double_t integratedTimes[5] = {-1.0,-1.0,-1.0,-1.0,-1.0};
+ track->GetIntegratedTimes(integratedTimes);
+ Float_t l = integratedTimes[0]*c;
+
+ goodtrack = goodtrack && (l > 365);
+
+ if (!goodtrack) return kFALSE;
+
+ if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
+
+ Bool_t statusMatchingHard = TPCTOFagree(track);
+ if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
+ return kFALSE;
+
+
+ Float_t beta = GetBeta(track);
+
+ //construct the pid index because it's not AliPID::EParticleType
+ Int_t pid = 0;
+ switch (fParticleID)
+ {
+ case AliPID::kPion:
+ pid=2;
+ break;
+ case AliPID::kKaon:
+ pid=3;
+ break;
+ case AliPID::kProton:
+ pid=4;
+ break;
+ default:
+ return kFALSE;
+ }
+
+ //signal to cut on
+ Float_t p = track->P();
+ Float_t betahypothesis = l/integratedTimes[pid]/c;
+ Float_t betadiff = beta-betahypothesis;
+
+ Float_t* arr = fTOFpidCuts->GetMatrixArray();
+ Int_t col = TMath::BinarySearch(fTOFpidCuts->GetNcols(),arr,static_cast<Float_t>(p));
+ if (col<0) return kFALSE;
+ Float_t min = (*fTOFpidCuts)(1,col);
+ Float_t max = (*fTOFpidCuts)(2,col);
+ Bool_t pass = (betadiff>min && betadiff<max);
+
+ return pass;
+}
//-----------------------------------------------------------------------
Bool_t AliFlowTrackCuts::PassesTOFbetaCut(const AliESDtrack* track )
{
if (!fAllowTOFmismatchFlag) {if ((track->GetStatus() & AliESDtrack::kTOFmismatch)) return kFALSE;}
+ if (!goodtrack) return kFALSE;
+
Bool_t statusMatchingHard = TPCTOFagree(track);
if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
return kFALSE;
-
- if (!goodtrack) return kFALSE;
Float_t beta = GetBeta(track);
}
}
+//-----------------------------------------------------------------------
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesTPCbayesianCut(const AliAODTrack* track)
+{
+ fBayesianResponse->ComputeProb(track,track->GetAODEvent()); // fCurrCentr is needed for mismatch fraction
+ Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mism$
+
+ Int_t kTPC = fBayesianResponse->GetCurrentMask(0); // is TPC on
+
+ if(! kTPC) return kFALSE;
+
+ fProbBayes = 0.0;
+
+switch (fParticleID)
+ {
+ case AliPID::kPion:
+ fProbBayes = probabilities[2];
+ break;
+ case AliPID::kKaon:
+ fProbBayes = probabilities[3];
+ break;
+ case AliPID::kProton:
+ fProbBayes = probabilities[4];
+ break;
+ case AliPID::kElectron:
+ fProbBayes = probabilities[0];
+ break;
+ case AliPID::kMuon:
+ fProbBayes = probabilities[1];
+ break;
+ case AliPID::kDeuteron:
+ fProbBayes = probabilities[5];
+ break;
+ case AliPID::kTriton:
+ fProbBayes = probabilities[6];
+ break;
+ case AliPID::kHe3:
+ fProbBayes = probabilities[7];
+ break;
+ default:
+ return kFALSE;
+ }
+
+ if(fProbBayes > fParticleProbability){
+ if(!fCutCharge)
+ return kTRUE;
+ else if (fCutCharge && fCharge * track->Charge() > 0)
+ return kTRUE;
+ }
+ return kFALSE;
+}
//-----------------------------------------------------------------------
Bool_t AliFlowTrackCuts::PassesTPCbayesianCut(const AliESDtrack* track)
{
// return kFALSE;
fBayesianResponse->ComputeProb(track,fCurrCentr); // fCurrCentr is needed for mismatch fraction
Int_t kTPC = fBayesianResponse->GetCurrentMask(0); // is TPC on
- Int_t kTOF = fBayesianResponse->GetCurrentMask(1); // is TOF on
+ //Int_t kTOF = fBayesianResponse->GetCurrentMask(1); // is TOF on
if(! kTPC) return kFALSE;
- Bool_t statusMatchingHard = 1;
- Float_t mismProb = 0;
- if(kTOF){
- statusMatchingHard = TPCTOFagree(track);
- mismProb = fBayesianResponse->GetTOFMismProb();
- }
- if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
- return kFALSE;
+ // Bool_t statusMatchingHard = 1;
+ // Float_t mismProb = 0;
+ // if(kTOF){
+ // statusMatchingHard = TPCTOFagree(track);
+ // mismProb = fBayesianResponse->GetTOFMismProb();
+ // }
+ // if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
+ // return kFALSE;
Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mismatch parameterization
return kFALSE;
}
- if(fProbBayes > fParticleProbability && mismProb < 0.5)
+ if(fProbBayes > fParticleProbability)
{
if(!fCutCharge)
return kTRUE;
}
return kFALSE;
}
+//-----------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(const AliAODTrack* track)
+{
+//check is track passes bayesian combined TOF+TPC pid cut
+ Bool_t goodtrack = (track->GetStatus() & AliESDtrack::kTOFout) &&
+ (track->GetStatus() & AliESDtrack::kTIME) &&
+ (track->GetTOFsignal() > 12000) &&
+ (track->GetTOFsignal() < 100000);
+
+ if (! goodtrack)
+ return kFALSE;
+
+ Bool_t statusMatchingHard = TPCTOFagree(track);
+//ciao
+ if (fRequireStrictTOFTPCagreement && (!statusMatchingHard))
+ return kFALSE;
+
+ fBayesianResponse->ComputeProb(track,track->GetAODEvent()); // fCurrCentr is needed for mismatch fraction
+ Float_t *probabilities = fBayesianResponse->GetProb(); // Bayesian Probability (from 0 to 4) (Combined TPC || TOF) including a tuning of priors and TOF mism$
+
+ Float_t mismProb = fBayesianResponse->GetTOFMismProb(); // mismatch Bayesian probabilities
+
+ fProbBayes = 0.0;
+
+ switch (fParticleID)
+ {
+ case AliPID::kPion:
+ fProbBayes = probabilities[2];
+ break;
+ case AliPID::kKaon:
+ fProbBayes = probabilities[3];
+ break;
+ case AliPID::kProton:
+ fProbBayes = probabilities[4];
+ break;
+ case AliPID::kElectron:
+ fProbBayes = probabilities[0];
+ break;
+ case AliPID::kMuon:
+ fProbBayes = probabilities[1];
+ break;
+ case AliPID::kDeuteron:
+ fProbBayes = probabilities[5];
+ break;
+ case AliPID::kTriton:
+ fProbBayes = probabilities[6];
+ break;
+ case AliPID::kHe3:
+ fProbBayes = probabilities[7];
+ break;
+ default:
+ return kFALSE;
+ }
+ if(fProbBayes > fParticleProbability && mismProb < 0.5){
+ if(!fCutCharge)
+ return kTRUE;
+ else if (fCutCharge && fCharge * track->Charge() > 0)
+ return kTRUE;
+ }
+ return kFALSE;
+
+}
//-----------------------------------------------------------------------
// part added by F. Noferini (some methods)
Bool_t AliFlowTrackCuts::PassesTOFbayesianCut(const AliESDtrack* track)
}
//---------------------------------------------------------------//
-Bool_t AliFlowTrackCuts::TPCTOFagree(const AliESDtrack *track)
+Bool_t AliFlowTrackCuts::TPCTOFagree(const AliVTrack *track)
{
//check pid agreement between TPC and TOF
Bool_t status = kFALSE;
-
+
+ const Float_t c = 2.99792457999999984e-02;
+
Float_t mass[5] = {5.10998909999999971e-04,1.05658000000000002e-01,1.39570000000000000e-01,4.93676999999999977e-01,9.38271999999999995e-01};
Float_t p = track->P();
Float_t time = track->GetTOFsignal()- fESDpid.GetTOFResponse().GetStartTime(p);
- Float_t tl = track->GetIntegratedLength();
+ Float_t tl = exptimes[0]*c; // to work both for ESD and AOD
Float_t betagammares = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[4], mass[4]);
else betagamma2 = 100;
- Double_t ptpc[3];
- track->GetInnerPxPyPz(ptpc);
- Float_t momtpc=TMath::Sqrt(ptpc[0]*ptpc[0] + ptpc[1]*ptpc[1] + ptpc[2]*ptpc[2]);
+ Float_t momtpc=track->GetTPCmomentum();
for(Int_t i=0;i < 5;i++){
Float_t resolutionTOF = fESDpid.GetTOFResponse().GetExpectedSigma(p, exptimes[i], mass[i]);
return "PMD";
case kV0:
return "V0";
+ case kMUON: // XZhang 20120604
+ return "MUON"; // XZhang 20120604
default:
return "unknown";
}
fTrackPhi = TMath::PiOver4()*(0.5+id%8);
- if(id<32)
- fTrackEta = -3.45+0.5*(id/8);
- else
- fTrackEta = +4.8-0.6*((id/8)-4);
+ // 10102013 weighting vzero tiles - rbertens@cern.ch
+ if(!fV0gainEqualization) {
+ // if for some reason the equalization is not initialized (e.g. 2011 data)
+ // the fV0xpol[] weights are used to enable or disable vzero rings
+ if(id<32) { // v0c side
+ fTrackEta = -3.45+0.5*(id/8);
+ if(id < 8) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Cpol[0];
+ else if (id < 16 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Cpol[1];
+ else if (id < 24 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Cpol[2];
+ else if (id < 32 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Cpol[3];
+ } else { // v0a side
+ fTrackEta = +4.8-0.6*((id/8)-4);
+ if( id < 40) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Apol[0];
+ else if ( id < 48 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Apol[1];
+ else if ( id < 56 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Apol[2];
+ else if ( id < 64 ) fTrackWeight = fEvent->GetVZEROEqMultiplicity(id)*fV0Apol[3];
+ }
+ } else { // the equalization is initialized
+ // note that disabled rings have already been excluded on calibration level in
+ // AliFlowEvent (so for a disabled ring, fV0xpol is zero
+ if(id<32) { // v0c side
+ fTrackEta = -3.45+0.5*(id/8);
+ if(id < 8) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Cpol[0]/fV0gainEqualization->GetBinContent(1+id);
+ else if (id < 16 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Cpol[1]/fV0gainEqualization->GetBinContent(1+id);
+ else if (id < 24 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Cpol[2]/fV0gainEqualization->GetBinContent(1+id);
+ else if (id < 32 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Cpol[3]/fV0gainEqualization->GetBinContent(1+id);
+ } else { // v0a side
+ fTrackEta = +4.8-0.6*((id/8)-4);
+ if( id < 40) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Apol[0]/fV0gainEqualization->GetBinContent(1+id);
+ else if ( id < 48 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Apol[1]/fV0gainEqualization->GetBinContent(1+id);
+ else if ( id < 56 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Apol[2]/fV0gainEqualization->GetBinContent(1+id);
+ else if ( id < 64 ) fTrackWeight = fEvent->GetVZEROData()->GetMultiplicity(id)*fV0Apol[3]/fV0gainEqualization->GetBinContent(1+id);
+ }
+ // printf ( " tile %i and weight %.2f \n", id, fTrackWeight);
+ }
- fTrackWeight = fEvent->GetVZEROEqMultiplicity(id);
-
- if (fLinearizeVZEROresponse)
+ if (fLinearizeVZEROresponse && id < 64)
{
//this is only needed in pass1 of LHC10h
Float_t multV0[fgkNumberOfV0tracks];
return pass;
}
+//-----------------------------------------------------------------------------
+Bool_t AliFlowTrackCuts::PassesMuonCuts(AliVParticle* vparticle)
+{
+// XZhang 20120604
+ fTrack=NULL;
+ fTrackLabel = (fFakesAreOK)?TMath::Abs(vparticle->GetLabel()):vparticle->GetLabel();
+ if (fMCevent) fMCparticle = static_cast<AliMCParticle*>(fMCevent->GetTrack(fTrackLabel));
+ else fMCparticle=NULL;
+
+ AliESDMuonTrack *esdTrack = dynamic_cast<AliESDMuonTrack*>(vparticle);
+ AliAODTrack *aodTrack = dynamic_cast<AliAODTrack*>(vparticle);
+ if ((!esdTrack) && (!aodTrack)) return kFALSE;
+ if (!fMuonTrackCuts->IsSelected(vparticle)) return kFALSE;
+ HandleVParticle(vparticle); if (!fTrack) return kFALSE;
+ return kTRUE;
+}
+
//----------------------------------------------------------------------------//
Double_t AliFlowTrackCuts::GetPmdEta(Float_t xPos, Float_t yPos, Float_t zPos)
{