//
// Event Loop
//
+
AliMCEventHandler* mcHandler = (dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()));
- AliESDInputHandler *inh = dynamic_cast<AliESDInputHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
- AliESDpid *workingPID = NULL;
- if(inh && (workingPID = inh->GetESDpid())) workingPID = inh->GetESDpid();
- else workingPID = AliHFEtools::GetDefaultPID(mcHandler ? kTRUE : kFALSE);
+
//AliESDInputHandler *inh = dynamic_cast<AliESDInputHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+ //AliESDpid *workingPID = NULL;
+ //if(inh && (workingPID = inh->GetESDpid())) workingPID = inh->GetESDpid();
+ //else workingPID = AliHFEtools::GetDefaultPID(mcHandler ? kTRUE : kFALSE);
// check the MC data
if(fMCEvent && !mcHandler ) return;
Int_t nV0s = fInputEvent->GetNumberOfV0s();
for(Int_t i=0; i<nV0s; ++i){
- AliESDv0 *esdV0 = (
dynamic_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
+ AliESDv0 *esdV0 = (
static_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
if(!esdV0) continue;
if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
Int_t pid = GetTenderPidV0(esdV0);
// V0 loop
for(Int_t i=0; i<nV0s; ++i){
Bool_t id = kFALSE;
- AliESDv0 *esdV0 = (
dynamic_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
+ AliESDv0 *esdV0 = (
static_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
if(!esdV0) continue;
if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
Int_t pid = GetTenderPidV0(esdV0);
class TH1F;
class TList;
class AliHFEcollection;
-class AliESDv0;
-class AliESDtrack;
class AliAnalysisTaskCheckV0tender : public AliAnalysisTaskSE{
public:
// Event Loop
//
AliMCEventHandler* mcHandler = (dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()));
- AliESDInputHandler *inh = dynamic_cast<AliESDInputHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
- AliESDpid *workingPID = NULL;
- if(inh && (workingPID = inh->GetESDpid())) workingPID = inh->GetESDpid();
- else workingPID = AliHFEtools::GetDefaultPID(mcHandler ? kTRUE : kFALSE);
+
+ //AliESDInputHandler *inh = dynamic_cast<AliESDInputHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+ //AliESDpid *workingPID = NULL;
+ //if(inh && (workingPID = inh->GetESDpid())) workingPID = inh->GetESDpid();
+ //else workingPID = AliHFEtools::GetDefaultPID(mcHandler ? kTRUE : kFALSE);
// check the MC data
if(fMCEvent && !mcHandler ) return;
Int_t nV0s = fInputEvent->GetNumberOfV0s();
for(Int_t i=0; i<nV0s; ++i){
- AliESDv0 *esdV0 = (
dynamic_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
+ AliESDv0 *esdV0 = (
static_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
if(!esdV0) continue;
if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
// check the identity of the V0tender selected V0 daughters
// !!! for positive check only the true identity plays a role here,
// not the true V0 mother identity (e.g. selected electron could come
- // from primary vertex or pion dalitz deca or true gamma conversion) !!!
+ // from primary vertex or pion dalitz deca or true gamma conversion) !!!!
//
const char * type[3] = {"Electron", "Pion", "Proton"};
class AliHFEcollection;
class AliESDv0KineCuts;
class AliKFVertex;
-class AliESDv0;
class AliAnalysisTaskCheckV0tenderII : public AliAnalysisTaskSE{
public:
Bool_t CheckSigns(AliESDv0 * const v0);
Int_t PDGtoPIDv0(Int_t pdgV0) const;
- Int_t PDGtoPID(Int_t pdg) const;
-
+ Int_t PDGtoPID(Int_t pdg) const;
TList *fOutput; //! Container for output histos
AliHFEcollection *fColl; //! collection of Data output
DefineOutput(1, TH1I::Class());
DefineOutput(2, TList::Class());
- printf(" ---> Dummy constructor used!\n");
+ //printf(" ---> Dummy constructor used!\n");
AliInfo("Dummy constructor used!");
fDefaultPID = new AliESDpid();
fHFEpid = new AliHFEpid("dummyPID");
hfecuts->SetCheckITSLayerStatus(kFALSE);
hfecuts->SetMaxImpactParam(maxDcaXY, maxDcaZ);
SetHFECuts(hfecuts);
-
- printf(" ---> Default constructor used!\n");
- cout<<" [!] DEFAULT CONSTRUCTOR! [!]"<<endl;
- AliInfo(" ---> Default constructor used!");
+
fDefaultPID = new AliESDpid();
fHFEpid = new AliHFEpid("PIDforDCAanalysis");
// create output objects
// fNEvents
// residual and pull
- printf("\n=====UserCreateOutputObjects=====\n");
+ //printf("\n=====UserCreateOutputObjects=====\n");
// Automatic determination of the analysis mode
AliVEventHandler *inputHandler = dynamic_cast<AliVEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
if(!fHFEpid) AliWarning("Hello, fHFEpid is not available");
cout<<" Hello this is a cout "<<endl<<endl;
- if(GetPlugin(kHFEpid)) {
- AliWarning(" ---> Adding TPC and TOF to the PID");
+ if(GetPlugin(kHFEpid)) {
fHFEpid->SetHasMCData(HasMCData());
fHFEpid->AddDetector("TOF", 0);
fHFEpid->AddDetector("TPC", 1);
//
// Run the analysis
//
- printf("\n=====UserExec=====\n");
+ //printf("\n=====UserExec=====\n");
if(HasMCData()) printf("WITH MC!\n");
AliDebug(3, "Processing ESD events");
//____________________________________________________________
void AliAnalysisTaskDCA::ProcessDcaAnalysis(){
- printf("\n=====ProcessDcaAnalysis=====\n");
+ //printf("\n=====ProcessDcaAnalysis=====\n");
//
// Loop ESD
//
-
- AliMCEvent *fMC = 0x0;
- AliESDVertex *vtxESDSkip = 0x0;
-
AliESDEvent *fESD = dynamic_cast<AliESDEvent *>(fInputEvent);
if(!fESD){
AliError("ESD Event required for ESD Analysis");
return;
}
+ AliMCEvent *fMC = 0x0;
if(HasMCData()){
fMC = dynamic_cast<AliMCEvent*>(fMCEvent);
if(!fMC){
}
}
-
fNEvents->Fill(1); // original event number before cut
fDCA->ApplyExtraCuts(fESD,fMinNprimVtxContrbutor); // cut on primVtx contributors
fNEvents->Fill(3); // events number after cut
-
- AliESDtrack *track = 0x0;
- AliMCParticle *mctrack = 0x0;
-
fCFM->SetRecEventInfo(fESD);
// event cut level
if(!fCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
+ AliESDtrack *track = 0x0;
+ AliMCParticle *mctrack = 0x0;
+ AliESDVertex *vtxESDSkip = 0x0;
+ AliHFEpidObject hfetrack;
+
for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
track = fESD->GetTrack(itrack);
- if(HasMCData())mctrack = dynamic_cast<AliMCParticle *>(fMC->GetTrack(TMath::Abs(track->GetLabel())));
+ if(HasMCData())
mctrack = dynamic_cast<AliMCParticle *>(fMC->GetTrack(TMath::Abs(track->GetLabel())));
// RecPrim: primary cuts
if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepRecPrim + AliHFEcuts::kNcutStepsMCTrack, track)) continue;
if(track->GetITSclusters(0)<=fNclustersITS) continue; // require number of ITS clusters
// track accepted, do PID
- AliHFEpidObject hfetrack;
hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
hfetrack.SetRecTrack(track);
if(HasMCData()) hfetrack.SetMCTrack(mctrack);
+ //printf("Track %d passed all the cuts!\n",itrack);
+
if(HasMCData()){
if(GetPlugin(kPrimVtx))
fDCA->FillHistogramsVtx(fESD, fMC);
if(GetPlugin(kCombinedPid))
fDCA->FillHistogramsPid(track, fMC);
if(GetPlugin(kHFEpid)) { // data-like
- if(fHFEpid->IsSelected(&hfetrack))
-
+ if(fHFEpid->IsSelected(&hfetrack)){
+
// printf("Found an electron in p+p collision! from HFE pid \n");
if(!vtxESDSkip){
// method from Andrea D 28.05.2010
delete vertexer; vertexer = NULL;
if(vtxESDSkip->GetNContributors()<fMinNprimVtxContrbutor) continue;
}
-
- fDCA->FillHistogramsHfeDataDca(fESD, track,vtxESDSkip);
-
+ //printf("\n[ABOUT TO FILL HFE DCA: MC!]\n");
+ fDCA->FillHistogramsHfeDataDca(fESD, track, vtxESDSkip);
+ }
} // plugin for hfepid
} // MC
delete vertexer; vertexer = NULL;
if(vtxESDSkip->GetNContributors()<fMinNprimVtxContrbutor) continue;
}
-
+ printf("\n[ABOUT TO FILL HFE DCA: DATA!]\n");
fDCA->FillHistogramsHfeDataDca(fESD, track,vtxESDSkip);
}
} // plugin for hfepid
//
// Terminate not implemented at the moment
//
- printf("\n=====Terminate=====\n");
+ //printf("\n=====Terminate=====\n");
if(GetPlugin(kPostProcess)){
fOutput = dynamic_cast<TList *>(GetOutputData(1));
//____________________________________________________________
void AliAnalysisTaskDCA::Load(TString filename){
- printf("\n=====Load=====\n");
+ //printf("\n=====Load=====\n");
// no need for postprocessing for the moment
TFile *input = TFile::Open(filename.Data());
// should do fitting here for dca resolution
// moved to an external macro to do the job
- printf("\n=====PostProcess=====\n");
+ //printf("\n=====PostProcess=====\n");
Load("HFEdca.root");
TCanvas *c1 = new TCanvas("c1", "number of analyzed events", 300, 400);
fNEvents->Draw();
//____________________________________________________________
void AliAnalysisTaskDCA::MakeParticleContainer(){
- printf("\n=====MakeParticleContainer=====\n");
+ //printf("\n=====MakeParticleContainer=====\n");
//
// Create the particle container (borrowed from AliAnalysisTaskHFE)
//
//
// Adding PID detector to the task
//
- printf("\n=====AddPIDdetector=====\n");
+ //printf("\n=====AddPIDdetector=====\n");
if(!fPIDdetectors.Length())
fPIDdetectors = detector;
class AliESDtrackCuts;
class AliMCEvent;
class AliVParticle;
-class AliESDpid;
class AliVEvent;
class AliVertexerTracks;
class AliHFEpid;
class AliHFEcuts;
class AliHFEextraCuts;
-class AliHFEdca;
class AliAnalysisTaskDCA : public AliAnalysisTaskSE{
public:
fDeCFM->SetMCEventInfo(fMCEvent);
- Double_t nContributor = 0;
+ Double_t nContributor[1] = {0};
const AliVVertex *mcPrimVtx = fMCEvent->GetPrimaryVertex();
- if(mcPrimVtx) nContributor = mcPrimVtx->GetNContributors();
+ if(mcPrimVtx) nContributor[0] = mcPrimVtx->GetNContributors();
//
// cut at MC event level
//
if(!fDeCFM->CheckEventCuts(AliHFEcuts::kEventStepGenerated, fMCEvent)) return;
- if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(&nContributor,AliHFEcuts::kEventStepGenerated);
+ if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(nContributor,AliHFEcuts::kEventStepGenerated);
AliStack *stack = 0x0;
}
fDeCFM->SetRecEventInfo(fESD);
- Double_t nContrib = fESD->GetPrimaryVertex()->GetNContributors();
+ Double_t nContrib[1] = {fESD->GetPrimaryVertex()->GetNContributors()};
Bool_t alreadyseen = kFALSE;
AliLabelContainer cont(fESD->GetNumberOfTracks());
// cut at ESD event level
//
if(!fDeCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
- if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(&nContrib, AliHFEcuts::kEventStepReconstructed);
+ if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(nContrib, AliHFEcuts::kEventStepReconstructed);
for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
track = fESD->GetTrack(itrack);
Int_t nHFEelectrons= 0;
fDeCFM->SetRecEventInfo(fESD);
- Double_t nContrib = fESD->GetPrimaryVertex()->GetNContributors();
+ Double_t nContrib[1] = {fESD->GetPrimaryVertex()->GetNContributors()};
if(!fDeCFM->CheckEventCuts(AliHFEcuts::kEventStepReconstructed, fESD)) return;
- if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(&nContrib, AliHFEcuts::kEventStepReconstructed);
+ if(GetPlugin(kCorrection)) fDeCFM->GetEventContainer()->Fill(nContrib, AliHFEcuts::kEventStepReconstructed);
for(Int_t itrack = 0; itrack < fESD->GetNumberOfTracks(); itrack++){
DefineOutput(2, TList::Class());
fPID = new AliHFEpid("hfePid");
+ fPIDqa = new AliHFEpidQAmanager;
fVarManager = new AliHFEvarManager("hfeVarManager");
}
fPID->InitializePID();
if(IsQAOn(kPIDqa)){
AliInfo("PID QA switched on");
- fPIDqa = new AliHFEpidQAmanager;
fPIDqa->Initialize(fPID);
fQA->Add(fPIDqa->MakeList("HFEpidQA"));
}
fOutput->Add(fTaggedTrackAnalysis->GetContainer());
fQA->Add(fTaggedTrackAnalysis->GetPIDQA());
fQA->Add(fTaggedTrackAnalysis->GetCutQA());
+ fQA->Add(fTaggedTrackAnalysis->GetQAcollection());
}
PrintStatus();
}
// need the centrality for everything (MC also)
fCentralityF = 99.0;
ReadCentrality();
+ //printf("fCentralityF %f\n",fCentralityF);
// See if pile up and z in the range
RejectionPileUpVertexRangeEventCut();
// Protect agains missing
if(HasMCData()){
+ //printf("Has MC data\n");
fSignalCuts->SetMCEvent(fMCEvent);
ProcessMC(); // Run the MC loop + MC QA in case MC Data are available
}
Double_t eventContainer [3];
eventContainer[0] = fMCEvent->GetPrimaryVertex()->GetZ();
eventContainer[2] = fCentralityF;
- if(fCFM->CheckEventCuts(AliHFEcuts::kEventStepGenerated, fMCEvent))
- fCFM->GetEventContainer()->Fill(eventContainer,AliHFEcuts::kEventStepGenerated);
+ //printf("z position is %f\n",eventContainer[0]);
+ //if(fCFM->CheckEventCuts(AliHFEcuts::kEventStepGenerated, fMCEvent))
+ fCFM->GetEventContainer()->Fill(eventContainer,AliHFEcuts::kEventStepGenerated);
Int_t nElectrons = 0;
if(IsESDanalysis()){
if (HasMCData() && IsQAOn(kMCqa)) {
// Do event Normalization
Double_t eventContainer[3];
- eventContainer[0] = fInputEvent->GetPrimaryVertex()->GetZ();
+ eventContainer[0] = 0.0;
+ if(fESD->GetPrimaryVertexTracks()) eventContainer[0] = fESD->GetPrimaryVertexTracks()->GetZ();
eventContainer[1] = 0.;
eventContainer[2] = fCentralityF;
if(fTriggerAnalysis->IsOfflineTriggerFired(fESD, AliTriggerAnalysis::kV0AND))
// Function is still in development
//
AliDebug(3, "Processing AOD Event");
- Double_t eventContainer[2];
+ Double_t eventContainer[3];
eventContainer[0] = fInputEvent->GetPrimaryVertex()->GetZ();
eventContainer[1] = 1.; // No Information available in AOD analysis, assume all events have V0AND
-
+ eventContainer[2] = fCentralityF;
+
AliAODEvent *fAOD = dynamic_cast<AliAODEvent *>(fInputEvent);
if(!fAOD){
AliError("AOD Event required for AOD Analysis")
Double_t vertex[3]; // Production vertex cut to mask gammas which are NOT supposed to have hits in the first ITS layer(s)
if(IsESDanalysis()){
AliMCParticle *mctrack = dynamic_cast<AliMCParticle *>(track);
- vertex[0] = mctrack->Particle()->Vx();
- vertex[1] = mctrack->Particle()->Vy();
+ if(mctrack){
+ vertex[0] = mctrack->Particle()->Vx();
+ vertex[1] = mctrack->Particle()->Vy();
+ }
} else {
AliAODMCParticle *aodmctrack = dynamic_cast<AliAODMCParticle *>(track);
- aodmctrack->XvYvZv(vertex);
+ if(aodmctrack) aodmctrack->XvYvZv(vertex);
}
if(!fCFM->CheckParticleCuts(AliHFEcuts::kStepMCGenerated, track)) return kFALSE;
Double_t *vertexBins = AliHFEtools::MakeLinearBinning(nBins[0], binMin[0], binMax[0]);
Double_t *v0andBins = AliHFEtools::MakeLinearBinning(nBins[1], binMin[1], binMax[1]);
+ Double_t *centralityBins = AliHFEtools::MakeLinearBinning(nBins[2], binMin[2], binMax[2]);
evCont->SetBinLimits(0, vertexBins);
evCont->SetBinLimits(1, v0andBins);
- delete[] vertexBins; delete[] v0andBins;
+ evCont->SetBinLimits(2, centralityBins);
+ delete[] vertexBins; delete[] v0andBins; delete[] centralityBins;
fCFM->SetEventContainer(evCont);
}
if(TString(mctrack->IsA()->GetName()).CompareTo("AliMCParticle") == 0){
// case MCParticle
- xv = (dynamic_cast<const AliMCParticle *>(mctrack)->Xv());
- yv = (dynamic_cast<const AliMCParticle *>(mctrack)->Yv());
-
+ const AliMCParticle *mcpart = dynamic_cast<const AliMCParticle *>(mctrack);
+ if(mcpart){
+ xv = mcpart->Xv();
+ yv = mcpart->Yv();
+ }
} else {
// case AODMCParticle
- xv = (dynamic_cast<const AliAODMCParticle *>(mctrack)->Xv());
- yv = (dynamic_cast<const AliAODMCParticle *>(mctrack)->Yv());
+ const AliAODMCParticle *mcpart = dynamic_cast<const AliAODMCParticle *>(mctrack);
+ if(mcpart){
+ xv = mcpart->Xv();
+ yv = mcpart->Yv();
+ }
}
//printf("xv %f, yv %f\n",xv,yv);
//AliAODCentrality *aodCentrality = fAOD->GetCentrality();
//Double_t fCentralityF = aodCentrality->GetCentralityPercentile("V0M");
fCentralityF = 99.0; // Fake for the moment
-
+
} else {
AliError("ESD Event required for ESD Analysis")
return;
}
+ const char* type = fESD->GetBeamType();
+
+ if (strstr(type,"Pb-Pb")) {
+
// Centrality
AliCentrality *esdCentrality = fESD->GetCentrality();
Float_t fCentralityF_temp = esdCentrality->GetCentralityPercentile("V0M");
else if ( fCentralityF_temp >= 80. && fCentralityF_temp < 90.) fCentralityF = 9;
else if ( fCentralityF_temp >= 90. && fCentralityF_temp <=100.) fCentralityF = 10;
+ }
+
+
+ if (strstr(type,"p-p")) {
+ fCentralityF = 0;
+ Int_t centralityF_temp = 0;
+ const AliESDVertex *vtxESD = fESD->GetPrimaryVertexTracks();
+ if(vtxESD && vtxESD->GetStatus()) centralityF_temp = vtxESD->GetNContributors();
+
+ //printf("centralityF_temp %d\n",centralityF_temp);
+
+ if( centralityF_temp <= 0) fCentralityF = 0;
+ else if ( centralityF_temp > 0 && centralityF_temp < 2) fCentralityF = 1;
+ else if ( centralityF_temp >= 2 && centralityF_temp < 3) fCentralityF = 2;
+ else if ( centralityF_temp >= 3 && centralityF_temp < 4) fCentralityF = 3;
+ else if ( centralityF_temp >= 4 && centralityF_temp < 5) fCentralityF = 4;
+ else if ( centralityF_temp >= 5 && centralityF_temp < 10) fCentralityF = 5;
+ else if ( centralityF_temp >= 10 && centralityF_temp < 20) fCentralityF = 6;
+ else if ( centralityF_temp >= 20 && centralityF_temp < 30) fCentralityF = 7;
+ else if ( centralityF_temp >= 30 && centralityF_temp < 40) fCentralityF = 8;
+ else if ( centralityF_temp >= 40 && centralityF_temp < 50) fCentralityF = 9;
+ else if ( centralityF_temp >= 50) fCentralityF = 10;
+
+
+ }
+
//printf("centrality %f\n",fCentralityF);
}
class AliMCEvent;
class AliVParticle;
class AliTriggerAnalysis;
-class AliESDtrack;
class TH1I;
class TList;
// Get Components for configuration
AliHFEvarManager *GetVarManager() const { return fVarManager; }
+ AliHFEpidQAmanager *GetPIDQAManager() const { return fPIDqa; }
AliHFEpid *GetPID() const { return fPID; }
void SetHFECuts(AliHFEcuts * const cuts) { fCuts = cuts; };
AliCFManager *fCFM; //! Correction Framework Manager
AliTriggerAnalysis *fTriggerAnalysis; //! Trigger Analysis for Normalisation
AliHFEpid *fPID; // PID
- AliHFEpidQAmanager *fPIDqa; //! PID QA
+ AliHFEpidQAmanager *fPIDqa; // PID QA
AliHFEpid *fPIDpreselect; // PID oject for pre-selected tracks (without QA)
AliHFEcuts *fCuts; // Cut Collection
AliHFEcuts *fTaggedTrackCuts; // Cut Collection for V0 tagged tracks
fV0(0x0)
, fEvent(0x0)
, fPrimaryVertex(0x0)
+ , fGcutChi2NDF(0)
+ , fGcutInvMass(0)
+ , fK0cutChi2NDF(0)
+ , fLcutChi2NDF(0)
{
//
// Default constructor
//
+ // default gamma cuts values
+ fGcutChi2NDF = 40; // Chi2NF cut value for the AliKFparticle gamma
+ fGcutCosPoint[0] = 0; // cos of the pointing angle [min, max]
+ fGcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
+ fGcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
+ fGcutDCA[1] = 0.25; // DCA between the daughter tracks [min, max]
+ fGcutVertexR[0] = 8.; // radius of the conversion point [min, max]
+ fGcutVertexR[1] = 90.; // radius of the conversion point [min, max]
+ fGcutPsiPair[0] = 0.; // value of the psi pair cut [min, max]
+ fGcutPsiPair[1] = 0.05; // value of the psi pair cut [min, max]
+ fGcutInvMass = 0.05; // upper value on the gamma invariant mass
+
+ fK0cutChi2NDF = 40; // Chi2NF cut value for the AliKFparticle K0
+ fK0cutCosPoint[0] = 0.; // cos of the pointing angle [min, max]
+ fK0cutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
+ fK0cutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
+ fK0cutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]
+ fK0cutVertexR[0] = 2.0; // radius of the decay point [min, max]
+ fK0cutVertexR[1] = 30.0; // radius of the decay point [min, max]
+ fK0cutInvMass[0] = 0.486; // invariant mass window
+ fK0cutInvMass[1] = 0.508; // invariant mass window
+ // Lambda & anti-Lambda cut values
+ fLcutChi2NDF = 40; // Chi2NF cut value for the AliKFparticle K0
+ fLcutCosPoint[0] = 0.; // cos of the pointing angle [min, max]
+ fLcutCosPoint[1] = 0.02; // cos of the pointing angle [min, max]
+ fLcutDCA[0] = 0.; // DCA between the daughter tracks [min, max]
+ fLcutDCA[1] = 0.2; // DCA between the daughter tracks [min, max]
+ fLcutVertexR[0] = 2.0; // radius of the decay point [min, max]
+ fLcutVertexR[1] = 40.0; // radius of the decay point [min, max]
+ fLcutInvMass[0] = 1.11; // invariant mass window
+ fLcutInvMass[1] = 1.12; // invariant mass window
+
}
//____________________________________________________________________
AliESDv0KineCuts::~AliESDv0KineCuts(){
// Gamma cuts
const Double_t cutAlphaG = 0.35;
- const Double_t cutAlphaG2[2] = {0.6, 0.8};
const Double_t cutQTG = 0.05;
+ const Double_t cutAlphaG2[2] = {0.6, 0.8};
const Double_t cutQTG2 = 0.04;
// K0 cuts
Float_t iMass = v0->GetEffMass(0, 0);
- // Cut values
- const Double_t cutMass = 0.05; // old [0.05]
- const Double_t cutChi2NDF = 40.; // old [7.]
- const Double_t cutCosPoint[2] = {0., 0.02}; // old [0., 0.03]
- const Double_t cutDCA[2] = {0., 0.25}; // old [0., 0.25]
- const Double_t cutProdVtxR[2] = {8., 90.}; // old [6., 9999]
- const Double_t cutPsiPair[2] = {0., 0.05}; // old [0. 0.05]
- const Double_t cutOAngle[2] = {0, 0.1}; // old [0., 0.1]
-
// cos pointing angle
Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
cosPoint = TMath::ACos(cosPoint);
r2 = TMath::Sqrt(xy[0]*xy[0] + xy[1]*xy[1]);
}
- // Opening angle
- Double_t oAngle = OpenAngle(v0);
-
// psi pair
Double_t psiPair = PsiPair(v0);
// apply the cuts
- if(iMass > cutMass) return kFALSE;
+ if(iMass > fGcutInvMass) return kFALSE;
- if(chi2ndf > cutChi2NDF) return kFALSE;
+ if(chi2ndf > fGcutChi2NDF) return kFALSE;
- if(cosPoint < cutCosPoint[0] || cosPoint > cutCosPoint[1]) return kFALSE;
+ if(cosPoint < fGcutCosPoint[0] || cosPoint > fGcutCosPoint[1]) return kFALSE;
- if(dca < cutDCA[0] || dca > cutDCA[1]) return kFALSE;
+ if(dca < fGcutDCA[0] || dca > fGcutDCA[1]) return kFALSE;
- if(r < cutProdVtxR[0] || r > cutProdVtxR[1]) return kFALSE;
+ if(r < fGcutVertexR[0] || r > fGcutVertexR[1]) return kFALSE;
- if(psiPair < cutPsiPair[0] || psiPair > cutPsiPair[1]) return kFALSE;
-
- if(oAngle < cutOAngle[0] || oAngle > cutOAngle[1]) return kFALSE;
+ if(psiPair < fGcutPsiPair[0] || psiPair > fGcutPsiPair[1]) return kFALSE;
// all cuts passed
Float_t iMass = v0->GetEffMass(2, 2);
- const Double_t cutMass[2] = {0.486, 0.508}; // ORG [0.485, 0.51]
- const Double_t cutChi2NDF = 40.; // ORG [7.]
- const Double_t cutCosPoint[2] = {0., 0.02}; // ORG [0., 0.03]
- const Double_t cutDCA[2] = {0., 0.2}; // ORG [0., 0.1]
- const Double_t cutProdVtxR[2] = {2.0, 30.}; // ORG [0., 8.1]
-
// cos pointing angle
Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
cosPoint = TMath::ACos(cosPoint);
//
// apply the cuts
//
- if(iMass < cutMass[0] || iMass > cutMass[1]) return kFALSE;
+ if(iMass < fK0cutInvMass[0] || iMass > fK0cutInvMass[1]) return kFALSE;
- if(chi2ndf > cutChi2NDF) return kFALSE;
+ if(chi2ndf > fK0cutChi2NDF) return kFALSE;
- if(cosPoint < cutCosPoint[0] || cosPoint > cutCosPoint[1]) return kFALSE;
+ if(cosPoint < fK0cutCosPoint[0] || cosPoint > fK0cutCosPoint[1]) return kFALSE;
- if(dca < cutDCA[0] || dca > cutDCA[1]) return kFALSE;
+ if(dca < fK0cutDCA[0] || dca > fK0cutDCA[1]) return kFALSE;
- if(r < cutProdVtxR[0] || r > cutProdVtxR[1]) return kFALSE;
+ if(r < fK0cutVertexR[0] || r > fK0cutVertexR[1]) return kFALSE;
// all cuts passed
pdgV0 = 310;
iMass = (type == 0) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
}
- // Cuts
- const Double_t cutMass[2] = {1.11, 1.12}; // ORG [1.11, 1.12]
- const Double_t cutChi2NDF = 40.; // ORG [5.]
- const Double_t cutCosPoint[2] = {0., 0.02}; // ORG [0., 0.03]
- const Double_t cutDCA[2] = {0., 0.2}; // ORG [0., 0.2]
- const Double_t cutProdVtxR[2] = {2., 40.}; // ORG [0., 24.]
-
// cos pointing angle
Double_t cosPoint = v0->GetV0CosineOfPointingAngle();
cosPoint = TMath::ACos(cosPoint);
// apply the cuts
//
- if(iMass < cutMass[0] || iMass > cutMass[1]) return kFALSE;
+ if(iMass < fLcutInvMass[0] || iMass > fLcutInvMass[1]) return kFALSE;
- if(chi2ndf > cutChi2NDF) return kFALSE;
+ if(chi2ndf > fLcutChi2NDF) return kFALSE;
- if(cosPoint < cutCosPoint[0] || cosPoint > cutCosPoint[1]) return kFALSE;
+ if(cosPoint < fLcutCosPoint[0] || cosPoint > fLcutCosPoint[1]) return kFALSE;
- if(dca < cutDCA[0] || dca > cutDCA[1]) return kFALSE;
+ if(dca < fLcutDCA[0] || dca > fLcutDCA[1]) return kFALSE;
- if(r < cutProdVtxR[0] || r > cutProdVtxR[1]) return kFALSE;
+ if(r < fLcutVertexR[0] || r > fLcutVertexR[1]) return kFALSE;
// all cuts passed
SetEvent(dynamic_cast<AliESDEvent*>(event));
}
//________________________________________________________________
-Double_t AliESDv0KineCuts::OpenAngle(AliESDv0 *v0) const {
- //
- // Opening angle between two daughter tracks
- //
- Double_t mn[3] = {0,0,0};
- Double_t mp[3] = {0,0,0};
-
-
- v0->GetNPxPyPz(mn[0],mn[1],mn[2]);//reconstructed cartesian momentum components of negative daughter;
- v0->GetPPxPyPz(mp[0],mp[1],mp[2]);//reconstructed cartesian momentum components of positive daughter;
-
-
- Double_t openAngle = TMath::ACos((mp[0]*mn[0] + mp[1]*mn[1] + mp[2]*mn[2])/(TMath::Sqrt(mp[0]*mp[0] + mp[1]*mp[1] + mp[2]*mp[2])*TMath::Sqrt(mn[0]*mn[0] + mn[1]*mn[1] + mn[2]*mn[2])));
-
- return TMath::Abs(openAngle);
-}
-//________________________________________________________________
Double_t AliESDv0KineCuts::PsiPair(AliESDv0* const v0) {
//
// Angle between daughter momentum plane and plane
void SetEvent(AliVEvent* const event);
void SetPrimaryVertex(AliKFVertex* const v) { fPrimaryVertex = v; };
- Double_t OpenAngle(AliESDv0 *v0) const;//opening angle between V0 daughters; close to zero for conversions
+ // setter functions for cut values
+ void SetGammaCutChi2NDF(Float_t val) { fGcutChi2NDF = val; };
+ void SetGammaCutCosPoint(Float_t *val) {
+ fGcutCosPoint[0] = val[0];
+ fGcutCosPoint[1] = val[1];
+ };
+ void SetGammaCutDCA(Float_t *val){
+ fGcutDCA[0] = val[0];
+ fGcutDCA[1] = val[1];
+ };
+ void SetGammaCutVertexR(Float_t *val){
+ fGcutVertexR[0] = val[0];
+ fGcutVertexR[1] = val[1];
+ };
+ void SetGammaCutPsiPair(Float_t *val){
+ fGcutPsiPair[0] = val[0];
+ fGcutPsiPair[1] = val[1];
+ };
+ void SetGammaCutInvMass(Float_t val){
+ fGcutInvMass = val;
+ };
+
Double_t PsiPair(AliESDv0* const v0);
Bool_t GetConvPosXY(AliESDtrack * const ptrack, AliESDtrack * const ntrack, Double_t convpos[2]);
AliESDEvent *fEvent; // current event
AliKFVertex *fPrimaryVertex; // primary vertex
+ // gamma cut values
+ Float_t fGcutChi2NDF; // Chi2NF cut value for the AliKFparticle gamma
+ Float_t fGcutCosPoint[2]; // cos of the pointing angle [min, max]
+ Float_t fGcutDCA[2]; // DCA between the daughter tracks [min, max]
+ Float_t fGcutVertexR[2]; // radius of the conversion point [min, max]
+ Float_t fGcutPsiPair[2]; // value of the psi pair cut [min, max]
+ Float_t fGcutInvMass; // upper value on the gamma invariant mass
+ // K0 cut values
+ Float_t fK0cutChi2NDF; // Chi2NF cut value for the AliKFparticle K0
+ Float_t fK0cutCosPoint[2]; // cos of the pointing angle [min, max]
+ Float_t fK0cutDCA[2]; // DCA between the daughter tracks [min, max]
+ Float_t fK0cutVertexR[2]; // radius of the decay point [min, max]
+ Float_t fK0cutInvMass[2]; // invariant mass window
+ // Lambda & anti-Lambda cut values
+ Float_t fLcutChi2NDF; // Chi2NF cut value for the AliKFparticle K0
+ Float_t fLcutCosPoint[2]; // cos of the pointing angle [min, max]
+ Float_t fLcutDCA[2]; // DCA between the daughter tracks [min, max]
+ Float_t fLcutVertexR[2]; // radius of the decay point [min, max]
+ Float_t fLcutInvMass[2]; // invariant mass window
+
+
ClassDef(AliESDv0KineCuts, 0);
};
if(!TString(fInputEvent->IsA()->GetName()).CompareTo("AliESDEvent")){
// case ESD
SetESDanalysis();
- AliESDv0 *esdV0 = (
dynamic_cast<AliESDEvent *>(fInputEvent))->GetV0(iv0);
+ AliESDv0 *esdV0 = (
static_cast<AliESDEvent *>(fInputEvent))->GetV0(iv0);
if(!esdV0) continue;
if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
v0status = ProcessV0(esdV0);
} else {
// case AOD
SetAODanalysis();
- AliAODv0 *aodV0 = (
dynamic_cast<AliAODEvent *>(fInputEvent))->GetV0(iv0);
+ AliAODv0 *aodV0 = (
static_cast<AliAODEvent *>(fInputEvent))->GetV0(iv0);
if(aodV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
v0status = ProcessV0(aodV0);
if(AliHFEV0cuts::kUndef != v0status){
// Process single V0
// Apply general cut and special cuts for gamma, K0s, Lambda
//
+ if(!v0) return AliHFEV0cuts::kUndef;
AliVTrack* daughter[2];
- daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack((
dynamic_cast<AliESDv0 *>(v0))->GetPindex()));
- daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack((
dynamic_cast<AliESDv0 *>(v0))->GetNindex()));
+ daughter[0] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack((
static_cast<AliESDv0 *>(v0))->GetPindex()));
+ daughter[1] = dynamic_cast<AliVTrack *>(fInputEvent->GetTrack((
static_cast<AliESDv0 *>(v0))->GetNindex()));
if(!daughter[0] || !daughter[1]) return AliHFEV0cuts::kUndef;
- if(fMCEvent) fMCon = kTRUE;
+ if(fMCEvent != NULL) fMCon = kTRUE;
//printf("-D: fMCEvent %x, fMCon: %i\n", fMCEvent, fMCon);
Int_t dMC[2] = {-1, -1};
}
// check common single track cuts
for(Int_t i=0; i<2; ++i){
- if(!fV0cuts->TrackCutsCommon(
dynamic_cast<AliESDtrack*>(daughter[i]))) return AliHFEV0cuts::kUndef;
+ if(!fV0cuts->TrackCutsCommon(
static_cast<AliESDtrack*>(daughter[i]))) return AliHFEV0cuts::kUndef;
}
// check commom V0 cuts
if(!fV0cuts->V0CutsCommon(dynamic_cast<AliESDv0 *>(v0))) return AliHFEV0cuts::kUndef;
//
// Identify Gamma
//
+
+ if(!v0) return kFALSE;
+
AliVTrack* daughter[2];
Int_t pIndex = 0, nIndex = 0;
Double_t invMass = 0.;
Int_t v0id = -1;
if(IsESDanalysis()){
// ESD - cut V0
- AliESDv0 *esdV0 =
dynamic_cast<AliESDv0 *>(v0);
+ AliESDv0 *esdV0 =
static_cast<AliESDv0 *>(v0);
v0id = esdV0->GetLabel();
// apply FULL gamma cuts
if(!fV0cuts->GammaCuts(esdV0)) return kFALSE;
mPt = esdV0->Pt();
} else {
// AOD Analysis - not possible to cut
- AliAODv0 *aodV0 =
dynamic_cast<AliAODv0 *>(v0);
+ AliAODv0 *aodV0 =
static_cast<AliAODv0 *>(v0);
v0id = aodV0->GetID();
pIndex = aodV0->GetPosID();
nIndex = aodV0->GetNegID();
//
// Identify K0s
//
+
+ if(!v0) return kFALSE;
+
AliVTrack* daughter[2];
Int_t pIndex = 0, nIndex = 0;
Int_t v0id = -1;
Double_t mPt = 0.;
if(IsESDanalysis()){
// ESD - cut V0
- AliESDv0 *esdV0 =
dynamic_cast<AliESDv0 *>(v0);
+ AliESDv0 *esdV0 =
static_cast<AliESDv0 *>(v0);
if(!fV0cuts->K0Cuts(esdV0)) return kFALSE;
v0id = esdV0->GetLabel();
pIndex = esdV0->GetPindex();
mPt = esdV0->Pt();
} else {
// AOD Analysis - not possible to cut
- AliAODv0 *aodV0 =
dynamic_cast<AliAODv0 *>(v0);
+ AliAODv0 *aodV0 =
static_cast<AliAODv0 *>(v0);
aodV0->GetID();
pIndex = aodV0->GetPosID();
nIndex = aodV0->GetNegID();
//const Double_t kPhiMass=TDatabasePDG::Instance()->GetParticle(333)->Mass(); // PDG phi mass
//AliVTrack* daughter[2];
//Double_t invMass = 0.;
+
+ if(!v0) return kFALSE;
Int_t pIndex = 0, nIndex = 0;
if(IsESDanalysis()){
// ESD - cut V0
- AliESDv0 *esdV0 =
dynamic_cast<AliESDv0 *>(v0);
+ AliESDv0 *esdV0 =
static_cast<AliESDv0 *>(v0);
pIndex = esdV0->GetPindex();
nIndex = esdV0->GetNindex();
} else {
//
// Identify Lambda
//
+
+ if(!v0) return kFALSE;
+
AliVTrack* daughter[2];
Int_t pIndex = 0, nIndex = 0;
Double_t invMass = 0.;
Int_t v0id = -1;
if(IsESDanalysis()){
// ESD - cut V0
- AliESDv0 *esdV0 =
dynamic_cast<AliESDv0 *>(v0);
+ AliESDv0 *esdV0 =
static_cast<AliESDv0 *>(v0);
v0id = esdV0->GetLabel();
if(!fV0cuts->LambdaCuts(esdV0,isLambda)) return kFALSE;
mPt = esdV0->Pt();
// AOD Analysis - not possible to cut
// again - two cases as above
- AliAODv0 *aodV0 =
dynamic_cast<AliAODv0 *>(v0);
+ AliAODv0 *aodV0 =
static_cast<AliAODv0 *>(v0);
v0id = aodV0->GetID();
pIndex = aodV0->GetPosID();
nIndex = aodV0->GetNegID();
fColl->Fill("h_QA_nParticles", 0);
// only ESD for now
AliESDtrack *track = dynamic_cast<AliESDtrack *>(recTrack);
+ if(!track) continue;
const AliExternalTrackParam *ext = track->GetOuterParam();
if(!ext) continue;
// MC label
//
// default constructor
//
-
- fList = new THashList();
- fList->SetOwner();
- if(!fList){
- AliError("Initialization of the list failed");
- }
- else{
- // list is owner of the objects. Once list is deleted, the objects
- // it contains will be deleted too
- //fList->SetOwner(kTRUE);
- }
- //Printf("%s:%d,%p",(char*)__FILE__,__LINE__,fInstance);
-
}
//___________________________________________________________________
AliHFEcollection::AliHFEcollection(const char* name, const char* title):
//
fList = new THashList();
- fList->SetOwner();
- fList->SetName(Form("list_%s", name));
- if(!fList){
- AliError("Initialization of the list failed");
- }
- else{
- // list is owner of the objects. Once list is deleted, the objects
- // it contains will be deleted too
- // fList->SetOwner(kTRUE);
+ if(fList){
+ fList->SetOwner();
+ fList->SetName(Form("list_%s", name));
}
}
//___________________________________________________________________
return kTRUE;
}
if(fList->FindObject(name)->InheritsFrom("TProfile")){
- (dynamic_cast<TProfile*>(fList->FindObject(name)))->Fill(v1, v2);
+ TProfile *pr = dynamic_cast<TProfile*>(fList->FindObject(name));
+ if(pr) pr->Fill(v1, v2);
return kTRUE;
}
if(list->IsEmpty())
return 1;
- TIterator *iter = list->MakeIterator();
+ TIter it(list);
TObject *o = NULL;
Int_t index = 0;
- while((o = iter->Next())){
+ TList templist; // Create temporary list containing all the lists to merge
+ while((o = it())){
AliHFEcollection *coll = dynamic_cast<AliHFEcollection *>(o);
if(!coll) continue;
- TList templist; // Create temporary list containing all the lists to merge
templist.Add(coll->fList);
- fList->Merge(&templist);
index++;
}
- delete iter;
+ fList->Merge(&templist);
return index + 1;
}
//____________________________________________________________________
virtual void Browse(TBrowser *b);
+ virtual Bool_t IsFolder() const { return kTRUE; }
// Set & Create functions
Bool_t CreateTH1F(const char* name, const char* title, Int_t nBin, Float_t nMin, Float_t nMax, Int_t logAxis = -1);
//
if(fNVars){
fVariables = new TObjArray(fNVars);
- for(UInt_t ivar = 0; ivar < fNVars; ivar++)
- fVariables->AddAt(new AliHFEvarInfo(*dynamic_cast<AliHFEvarInfo *>(ref.fVariables->UncheckedAt(ivar))), ivar);
+ AliHFEvarInfo *vtmp = NULL;
+ for(UInt_t ivar = 0; ivar < fNVars; ivar++){
+ vtmp = dynamic_cast<AliHFEvarInfo *>(ref.fVariables->UncheckedAt(ivar));
+ if(vtmp) fVariables->AddAt(new AliHFEvarInfo(*vtmp), ivar);
+ }
}
fContainers = new THashList;
fContainers->SetOwner();
AliCFContainer *ctmp = NULL;
for(Int_t ien = 0; ien < ref.fContainers->GetEntries(); ien++){
ctmp = dynamic_cast<AliCFContainer *>(ref.fContainers->At(ien));
- CreateContainer(ctmp->GetName(), ctmp->GetTitle(), ctmp->GetNStep());
+ if(ctmp) CreateContainer(ctmp->GetName(), ctmp->GetTitle(), ctmp->GetNStep());
}
// Copy also correlation matrices
if(ref.fCorrelationMatrices){
fCorrelationMatrices->SetOwner();
for(Int_t ien = 0; ien < ref.fCorrelationMatrices->GetEntries(); ien++){
htmp = dynamic_cast<THnSparseF *>(ref.fCorrelationMatrices->At(ien));
- CreateCorrelationMatrix(htmp->GetName(), htmp->GetTitle());
+ if(htmp) CreateCorrelationMatrix(htmp->GetName(), htmp->GetTitle());
}
}
}
TNamed::operator=(ref);
fContainers = new THashList();
fNVars = ref.fNVars;
- for(Int_t ien = 0; ien < ref.fContainers->GetEntries(); ien++)
- fContainers->Add(new AliCFContainer(*dynamic_cast<AliCFContainer *>(ref.fContainers->At(ien))));
+ AliCFContainer *ctmp = NULL;
+ for(Int_t ien = 0; ien < ref.fContainers->GetEntries(); ien++){
+ ctmp = dynamic_cast<AliCFContainer *>(ref.fContainers->At(ien));
+ fContainers->Add(new AliCFContainer(*ctmp));
+ }
if(fNVars){
fVariables = new TObjArray(fNVars);
- for(UInt_t ivar = 0; ivar < fNVars; ivar++)
- fVariables->AddAt(new AliHFEvarInfo(*dynamic_cast<AliHFEvarInfo *>(ref.fVariables->UncheckedAt(ivar))), ivar);
+ AliHFEvarInfo *vtmp = NULL;
+ for(UInt_t ivar = 0; ivar < fNVars; ivar++){
+ vtmp = dynamic_cast<AliHFEvarInfo *>(ref.fVariables->UncheckedAt(ivar));
+ if(vtmp) fVariables->AddAt(new AliHFEvarInfo(*vtmp), ivar);
+ }
} else {
fVariables = NULL;
}
fCorrelationMatrices->SetOwner();
for(Int_t ien = 0; ien < ref.fCorrelationMatrices->GetEntries(); ien++){
htmp = dynamic_cast<THnSparseF *>(ref.fCorrelationMatrices->At(ien));
- CreateCorrelationMatrix(htmp->GetName(), htmp->GetTitle());
+ if(htmp) CreateCorrelationMatrix(htmp->GetName(), htmp->GetTitle());
}
}
return *this;
}
Int_t *nBins = new Int_t[fNVars];
- for(UInt_t ivar = 0; ivar < fNVars; ivar++) nBins[ivar] = (dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar)))->GetNumberOfBins();
AliHFEvarInfo *var = NULL;
+ for(UInt_t ivar = 0; ivar < fNVars; ivar++){
+ var = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar));
+ nBins[ivar] = var ? var->GetNumberOfBins() : 0;
+ }
AliCFContainer *cont = new AliCFContainer(name, title, nStep, fNVars, nBins);
for(UInt_t ivar = 0; ivar < fNVars; ivar++){
var = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar));
- cont->SetBinLimits(ivar, var->GetBinning());
- cont->SetVarTitle(ivar, var->GetVarName()->Data());
+ if(var){
+ cont->SetBinLimits(ivar, var->GetBinning());
+ cont->SetVarTitle(ivar, var->GetVarName()->Data());
+ }
}
delete[] nBins;
fContainers->Add(cont);
AliHFEvarInfo *var = NULL;
for(UInt_t ivar = 0; ivar < fNVars; ivar++){
var = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar));
- nBins[ivar] = var->GetNumberOfBins();
- nBins[ivar+fNVars] = var->GetNumberOfBins();
+ if(var){
+ nBins[ivar] = var->GetNumberOfBins();
+ nBins[ivar+fNVars] = var->GetNumberOfBins();
+ }
}
THnSparseF * hTmp = new THnSparseF(name, title, 2*fNVars, nBins);
for(UInt_t ivar = 0; ivar < fNVars; ivar++){
var = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar));
- hTmp->SetBinEdges(ivar,var->GetBinning());
- //hTmp->GetAxis(ivar)->Set(var->GetNumberOfBins(), var->GetBinning());
- hTmp->GetAxis(ivar)->SetTitle(var->GetVarName()->Data());
- //hTmp->GetAxis(ivar + fNVars)->Set(var->GetNumberOfBins(), var->GetBinning());
- hTmp->GetAxis(ivar + fNVars)->SetTitle(Form("%s_{MC}", var->GetVarName()->Data()));
- hTmp->SetBinEdges(ivar+fNVars,var->GetBinning());
+ if(var){
+ hTmp->SetBinEdges(ivar,var->GetBinning());
+ //hTmp->GetAxis(ivar)->Set(var->GetNumberOfBins(), var->GetBinning());
+ hTmp->GetAxis(ivar)->SetTitle(var->GetVarName()->Data());
+ //hTmp->GetAxis(ivar + fNVars)->Set(var->GetNumberOfBins(), var->GetBinning());
+ hTmp->GetAxis(ivar + fNVars)->SetTitle(Form("%s_{MC}", var->GetVarName()->Data()));
+ hTmp->SetBinEdges(ivar+fNVars,var->GetBinning());
+ }
}
hTmp->Sumw2();
fCorrelationMatrices->AddLast(hTmp);
Int_t *dummyBinning = new Int_t[fNVars];
for(UInt_t ibin = 0; ibin < fNVars; ibin++) dummyBinning[ibin] = 1;
AliCFContainer *cmerged = new AliCFContainer(name, title, nStepMerged, fNVars, dummyBinning);
- delete dummyBinning;
+ delete[] dummyBinning;
// Fill container with content
AliInfo("Filling new container");
Int_t cstep = 0;
//
// Set Linear binning for the given container
//
- (dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(var)))->SetBinning(nBins, AliHFEtools::MakeLinearBinning(nBins, begin, end));
+ AliHFEvarInfo *myvar = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(var));
+ if(myvar) myvar->SetBinning(nBins, AliHFEtools::MakeLinearBinning(nBins, begin, end));
}
//__________________________________________________________________
//
// Set Logarithmic binning for the given container
//
- (dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(var)))->SetBinning(nBins, AliHFEtools::MakeLogarithmicBinning(nBins, begin, end));
+ AliHFEvarInfo *myvar = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(var));
+ if(myvar) myvar->SetBinning(nBins, AliHFEtools::MakeLogarithmicBinning(nBins, begin, end));
}
//__________________________________________________________________
//
// Variable name
//
- (dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(var)))->SetVarName(varname);
+ AliHFEvarInfo *myvar = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(var));
+ if(myvar) myvar->SetVarName(varname);
}
//__________________________________________________________________
if(nVars != fNVars)
std::cout << "Inconsistency in number of Variables [" << fNVars << "|" << nVars << "]" << std::endl;
AliHFEvarInfo *var = NULL;
- for(UInt_t ivar = 0; ivar < fNVars; ivar++){
- var = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar));
- std::cout << "Variable " << ivar << ": Name: " << var->GetVarName()->Data() << ", Number of Bins: " << var->GetNumberOfBins() << std::endl;
+ if(fVariables){
+ for(UInt_t ivar = 0; ivar < fNVars; ivar++){
+ var = dynamic_cast<AliHFEvarInfo *>(fVariables->UncheckedAt(ivar));
+ if(var)
+ std::cout << "Variable " << ivar << ": Name: " << var->GetVarName()->Data() << ", Number of Bins: " << var->GetNumberOfBins() << std::endl;
+ }
}
}
std::cout << std::endl;
// Print CF Containers:
- std::cout << "Containers[" << fContainers->GetEntries() << "]: "<< std::endl;
- std::cout << "=====================================================\n";
- for(Int_t icont = 0; icont < fContainers->GetEntries(); icont++){
- AliCFContainer *c = dynamic_cast<AliCFContainer *>(fContainers->At(icont));
- std::cout << "Name: " << c->GetName() << ", Title: " << c->GetTitle() << std::endl;
- for(Int_t istep = 0; istep < c->GetNStep(); istep++)
- std::cout << "Step " << istep << ": Title " << c->GetStepTitle(istep) << std::endl;
- std::cout << "------------------------------------------------------\n";
+ if(fContainers){
+ std::cout << "Containers[" << fContainers->GetEntries() << "]: "<< std::endl;
+ std::cout << "=====================================================\n";
+ for(Int_t icont = 0; icont < fContainers->GetEntries(); icont++){
+ AliCFContainer *c = dynamic_cast<AliCFContainer *>(fContainers->At(icont));
+ if(c){
+ std::cout << "Name: " << c->GetName() << ", Title: " << c->GetTitle() << std::endl;
+ for(Int_t istep = 0; istep < c->GetNStep(); istep++)
+ std::cout << "Step " << istep << ": Title " << c->GetStepTitle(istep) << std::endl;
+ }
+ std::cout << "------------------------------------------------------\n";
+ }
}
std::cout << "Number of Events: " << fNEvents << std::endl;
}
//
// Set MC information to the cuts in the cut step
//
+ AliCFCutBase *cfc = NULL;
for(Int_t icut = 0; icut < fCuts->GetEntriesFast(); icut++){
- if(fCuts->UncheckedAt(icut)->InheritsFrom("AliCFCutBase"))
- (dynamic_cast<AliCFCutBase *>(fCuts->UncheckedAt(icut)))->SetMCEventInfo(mc);
+ if((cfc = dynamic_cast<AliCFCutBase *>(fCuts->UncheckedAt(icut)))) cfc->SetMCEventInfo(mc);
}
}
//
// Publish rec event to the cut step
//
+ AliCFCutBase *cfc = NULL;
for(Int_t icut = 0; icut < fCuts->GetEntriesFast(); icut++){
- if(fCuts->UncheckedAt(icut)->InheritsFrom("AliCFCutBase"))
- (dynamic_cast<AliCFCutBase *>(fCuts->UncheckedAt(icut)))->SetRecEventInfo(rec);
+ if((cfc = dynamic_cast<AliCFCutBase *>(fCuts->UncheckedAt(icut)))) cfc->SetRecEventInfo(rec);
}
}
}
if(fCutList){
target.fCutList = dynamic_cast<TObjArray *>(fCutList->Clone());
- target.fCutList->SetOwner();
+ if(target.fCutList) target.fCutList->SetOwner(); // Coverity
}
if(target.fHistQA){
target.fHistQA->Clear();
class TString;
class TList;
-class TParticle;
class AliLog;
if(fCutTRD){
AliHFEcutStep *hfeTRD = fFilter->GetCutStep("HFETRD");
AliHFEextraCuts *hfecuts = dynamic_cast<AliHFEextraCuts *>(hfeTRD->GetCut("HFETRDCuts"));
- hfecuts->SetMinTrackletsTRD(4);
+ if(hfecuts) hfecuts->SetMinTrackletsTRD(4);
}
AliHFEcutStep *hfeITS = fFilter->GetCutStep("HFEITS");
AliHFEextraCuts *hfeitscuts = dynamic_cast<AliHFEextraCuts *>(hfeITS->GetCut("HFEPixelsCuts"));
- hfeitscuts->SetRequireITSpixel(AliHFEextraCuts::kFirst);
+ if(hfeitscuts)hfeitscuts->SetRequireITSpixel(AliHFEextraCuts::kFirst);
AliHFEcutStep *primary = fFilter->GetCutStep("Primary");
AliCFTrackIsPrimaryCuts *primcuts = dynamic_cast<AliCFTrackIsPrimaryCuts *>(primary->GetCut("PrimaryCuts"));
- primcuts->SetMaxDCAToVertexXY(3);
- primcuts->SetMaxDCAToVertexZ(5);
+ if(primcuts){
+ primcuts->SetMaxDCAToVertexXY(3);
+ primcuts->SetMaxDCAToVertexZ(5);
+ }
fAcceptanceCuts = new AliCFAcceptanceCuts("Acceptance", "MC Acceptance Cuts");
fAcceptanceCuts->SetMinNHitITS(3);
PostProcess();
TList *l = dynamic_cast<TList *>(GetOutputData(2));
- DrawPtResolution(l);
+ if(l) DrawPtResolution(l);
}
void AliHFEefficiency::FilterMC(){
//
TFile *input = TFile::Open(filename);
AliHFEcontainer *cin = dynamic_cast<AliHFEcontainer *>(input->Get("Efficiency"));
- fEfficiency = dynamic_cast<AliHFEcontainer *>(cin->Clone());
+
if(cin) fEfficiency = dynamic_cast<AliHFEcontainer *>(cin->Clone());
input->Close();
delete input;
}
// Draw pt resolution
//
TH2 *h2 = dynamic_cast<TH2 *>(l->FindObject("ptres"));
+ if(!h2) return;
TGraphErrors *mean = new TGraphErrors(h2->GetNbinsX());
TGraphErrors *rms = new TGraphErrors(h2->GetNbinsX());
Double_t norradius = TMath::Sqrt(fkVertex->GetX()*fkVertex->GetX()+fkVertex->GetY()*fkVertex->GetY());
- AliESDtrack *trackCopy = new AliESDtrack(*track);
- AliESDtrack *trackPartCopy = new AliESDtrack(*trackPart);
+ AliESDtrack trackCopy = AliESDtrack(*track);
+ AliESDtrack trackPartCopy = AliESDtrack(*trackPart);
Bool_t propagateok = kTRUE;
- if((!(trackPartCopy->PropagateTo(norradius,fBz))) || (!(trackCopy->PropagateTo(norradius,fBz)))) propagateok = kFALSE;
+ if((!(trackPartCopy.PropagateTo(norradius,fBz))) || (!(trackCopy.PropagateTo(norradius,fBz)))) propagateok = kFALSE;
if(!propagateok) {
- if(trackCopy) delete trackCopy;
- if(trackPartCopy) delete trackPartCopy;
+ //if(trackCopy) delete trackCopy;
+ //if(trackPartCopy) delete trackPartCopy;
return;
}
- CalculateMotherVariable(trackCopy,trackPartCopy,&results[0]);
- CalculateMotherVariableR(trackCopy,trackPartCopy,&resultsr[0]);
+ CalculateMotherVariable(&trackCopy,&trackPartCopy,&results[0]);
+ CalculateMotherVariableR(&trackCopy,&trackPartCopy,&resultsr[0]);
- if(trackCopy) delete trackCopy;
- if(trackPartCopy) delete trackPartCopy;
+ //if(trackCopy) delete trackCopy;
+ //if(trackPartCopy) delete trackPartCopy;
}
else {
class AliMCEvent;
class AliHFEpid;
-class TH2F;
//________________________________________________________________
if(IsQAOn()){
fIsQAOn = kTRUE;
fQAlist = dynamic_cast<TList *>(c.fQAlist->Clone());
- fQAlist->SetOwner();
+ if(fQAlist) fQAlist->SetOwner();
}
}
if(IsQAOn()){
fIsQAOn = kTRUE;
fQAlist = dynamic_cast<TList *>(c.fQAlist->Clone());
- fQAlist->SetOwner();
+ if(fQAlist) fQAlist->SetOwner();
}else fQAlist = 0x0;
}
return *this;
Float_t impactR, impactZ;
GetImpactParameters(track, impactR, impactZ);
Int_t nTPCf = GetTPCfindableClusters(track, fTPCiter1), nclsTPC = GetTPCncls(track, fTPCiter1);
- (dynamic_cast<TH1F *>(fQAlist->At(0 + when * kNhistos)))->Fill(impactR);
- (dynamic_cast<TH1F *>(fQAlist->At(1 + when * kNhistos)))->Fill(impactZ);
+ TH1 *htmp = NULL;
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(0 + when * kNhistos)))) htmp->Fill(impactR);
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(1 + when * kNhistos)))) htmp->Fill(impactZ);
// printf("TPC findable clusters: %d, found Clusters: %d\n", track->GetTPCNclsF(), track->GetTPCNcls());
- (dynamic_cast<TH1F *>(fQAlist->At(2 + when * kNhistos)))->Fill(nTPCf > 0. ? static_cast<Float_t>(nclsTPC)/static_cast<Float_t>(nTPCf) : 1.);
- (dynamic_cast<TH1F *>(fQAlist->At(3 + when * kNhistos)))->Fill(GetTRDnTrackletsPID(track));
- (dynamic_cast<TH1F *>(fQAlist->At(4 + when * kNhistos)))->Fill(nclsTPC);
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(2 + when * kNhistos)))) htmp->Fill(nTPCf > 0. ? static_cast<Float_t>(nclsTPC)/static_cast<Float_t>(nTPCf) : 1.);
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(3 + when * kNhistos)))) htmp->Fill(GetTRDnTrackletsPID(track));
+ if((htmp = dynamic_cast<TH1F *>(fQAlist->At(4 + when * kNhistos)))) htmp->Fill(nclsTPC);
UChar_t itsPixel = track->GetITSClusterMap();
TH1 *pixelHist = dynamic_cast<TH1F *>(fQAlist->At(5 + when * kNhistos));
//Int_t firstEntry = pixelHist->GetXaxis()->GetFirst();
- Double_t firstEntry = 0.5;
- if(!((itsPixel & BIT(0)) || (itsPixel & BIT(1))))
- pixelHist->Fill(firstEntry + 3);
- else{
- if(itsPixel & BIT(0)){
- pixelHist->Fill(firstEntry);
- if(itsPixel & BIT(1)) pixelHist->Fill(firstEntry + 2);
- else pixelHist->Fill(firstEntry + 4);
- }
- if(itsPixel & BIT(1)){
- pixelHist->Fill(firstEntry + 1);
- if(!(itsPixel & BIT(0))) pixelHist->Fill(firstEntry + 5);
+ if(pixelHist){
+ Double_t firstEntry = 0.5;
+ if(!((itsPixel & BIT(0)) || (itsPixel & BIT(1))))
+ pixelHist->Fill(firstEntry + 3);
+ else{
+ if(itsPixel & BIT(0)){
+ pixelHist->Fill(firstEntry);
+ if(itsPixel & BIT(1)) pixelHist->Fill(firstEntry + 2);
+ else pixelHist->Fill(firstEntry + 4);
+ }
+ if(itsPixel & BIT(1)){
+ pixelHist->Fill(firstEntry + 1);
+ if(!(itsPixel & BIT(0))) pixelHist->Fill(firstEntry + 5);
+ }
}
}
}
//
const Int_t kNhistos = 6;
TH2 *correlation = dynamic_cast<TH2F *>(fQAlist->At(2 * kNhistos));
- for(Int_t icut = 0; icut < kNcuts; icut++){
- if(!TESTBIT(fRequirements, icut)) continue;
- for(Int_t jcut = icut; jcut < kNcuts; jcut++){
- if(!TESTBIT(fRequirements, jcut)) continue;
- if(TESTBIT(survivedCut, icut) && TESTBIT(survivedCut, jcut))
- correlation->Fill(icut, jcut);
+ if(correlation){
+ for(Int_t icut = 0; icut < kNcuts; icut++){
+ if(!TESTBIT(fRequirements, icut)) continue;
+ for(Int_t jcut = icut; jcut < kNcuts; jcut++){
+ if(!TESTBIT(fRequirements, jcut)) continue;
+ if(TESTBIT(survivedCut, icut) && TESTBIT(survivedCut, jcut))
+ correlation->Fill(icut, jcut);
+ }
}
}
}
Int_t nTracklets = 0;
if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- nTracklets = esdtrack->GetTRDntrackletsPID();
+ if(esdtrack) nTracklets = esdtrack->GetTRDntrackletsPID();
} else if(!TString(track->IsA()->GetName()).CompareTo("AliAODTrack")){
AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
- AliAODPid *pidobject = aodtrack->GetDetPid();
+ AliAODPid *pidobject = NULL;
+ if(aodtrack) pidobject = aodtrack->GetDetPid();
// this is normally NOT the way to do this, but due to limitation in the
// AOD track it is not possible in a different way
if(pidobject){
Int_t det;
Float_t xloc, zloc;
AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- esdtrack->GetITSModuleIndexInfo(layer, det, status, xloc, zloc);
+ if(esdtrack) esdtrack->GetITSModuleIndexInfo(layer, det, status, xloc, zloc);
}
return status;
}
Int_t nClusters = 159; // in case no Information available consider all clusters findable
if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- nClusters = esdtrack->GetTPCNclsF();
+ if(esdtrack) nClusters = esdtrack->GetTPCNclsF();
}
return nClusters;
}
TString type = track->IsA()->GetName();
if(!type.CompareTo("AliESDtrack")){
AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- if(iter1)
- nClusters = esdtrack->GetTPCNclsIter1();
- else
- nClusters = esdtrack->GetTPCNcls();
+ if(esdtrack){ // coverity
+ if(iter1)
+ nClusters = esdtrack->GetTPCNclsIter1();
+ else
+ nClusters = esdtrack->GetTPCNcls();
+ }
}
else if(!type.CompareTo("AliAODTrack")){
AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
- const TBits &tpcmap = aodtrack->GetTPCClusterMap();
- for(UInt_t ibit = 0; ibit < tpcmap.GetNbits(); ibit++)
- if(tpcmap.TestBitNumber(ibit)) nClusters++;
-
+ const TBits &tpcmap = aodtrack->GetTPCClusterMap();
+ for(UInt_t ibit = 0; ibit < tpcmap.GetNbits(); ibit++)
+ if(tpcmap.TestBitNumber(ibit)) nClusters++;
}
return nClusters;
}
TString type = track->IsA()->GetName();
if(!type.CompareTo("AliESDtrack")){
AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
- esdtrack->GetImpactParameters(radial, z);
+ if(esdtrack) esdtrack->GetImpactParameters(radial, z);
}
else if(!type.CompareTo("AliAODTrack")){
AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
- Double_t xyz[3];
- aodtrack->XYZAtDCA(xyz);
- z = xyz[2];
- radial = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]+xyz[1]);
+ if(aodtrack){
+ Double_t xyz[3];
+ aodtrack->XYZAtDCA(xyz);
+ z = xyz[2];
+ radial = TMath::Sqrt(xyz[0]*xyz[0] + xyz[1]+xyz[1]);
+ }
}
}
AliESDtrack *esdtrack = NULL;
if(!TString(track->IsA()->GetName()).CompareTo("AliESDtrack")){
// Case ESD track: take copy constructor
- esdtrack = new AliESDtrack(*dynamic_cast<AliESDtrack *>(track));
+ AliESDtrack *tmptrack = dynamic_cast<AliESDtrack *>(track);
+ if(tmptrack) esdtrack = new AliESDtrack(*tmptrack);
} else {
// Case AOD track: take different constructor
esdtrack = new AliESDtrack(track);
// Get HFE impact parameter cut(pt dependent)
//
- TString type = track->IsA()->GetName();
- if(!type.CompareTo("AliESDtrack")){
-
AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
-
- Double_t pt = esdtrack->Pt();
- //hfeimpactRcut=0.0064+0.078*exp(-0.56*pt); // used Carlo's old parameter
- hfeimpactRcut=0.011+0.077*exp(-0.65*pt); // used Carlo's new parameter
- hfeimpactnsigmaRcut=3; // 3 sigma trail cut
+ TString type = track->IsA()->GetName();
+ if(!type.CompareTo("AliESDtrack")){
+ AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+ if(!esdtrack) return;
+ Double_t pt = esdtrack->Pt();
+ //hfeimpactRcut=0.0064+0.078*exp(-0.56*pt); // used Carlo's old parameter
+ hfeimpactRcut=0.011+0.077*exp(-0.65*pt); // used Carlo's new parameter
+ hfeimpactnsigmaRcut=3; // 3 sigma trail cut
}
}
//_______________________________________________________________________________________________
AliHFEmcQA::AliHFEmcQA() :
- fMCEvent(NULL)
+ fMCEvent(NULL)
,fMCHeader(NULL)
,fMCArray(NULL)
,fQAhistos(NULL)
{
// Default constructor
+ for(Int_t mom = 0; mom < 50; mom++){
+ fHeavyQuark[mom] = NULL;
+ }
+ for(Int_t mom = 0; mom < 2; mom++){
+ fIsHeavy[mom] = 0;
+ }
+
}
//_______________________________________________________________________________________________
,fNparents(p.fNparents)
{
// Copy constructor
+ for(Int_t mom = 0; mom < 50; mom++){
+ fHeavyQuark[mom] = NULL;
+ }
+ for(Int_t mom = 0; mom < 2; mom++){
+ fIsHeavy[mom] = 0;
+ }
}
//_______________________________________________________________________________________________
{
// decay particle's origin
- if ( abs(mcpart->GetPdgCode()) != AliHFEmcQA::kElectronPDG ) return kMisID;
-
- Int_t origin = -1;
- Bool_t isFinalOpenCharm = kFALSE;
-
if(!mcpart){
AliDebug(1, "no mcparticle, return\n");
return -1;
}
+ if ( abs(mcpart->GetPdgCode()) != AliHFEmcQA::kElectronPDG ) return kMisID;
+
+ Int_t origin = -1;
+ Bool_t isFinalOpenCharm = kFALSE;
+
Int_t iLabel = mcpart->GetFirstMother();
if (iLabel<0){
AliDebug(1, "Stack label is negative, return\n");
AliHistsComm fHistComm[2][6]; // struct of additional histograms of given particles
TList *fQAhistos; // Container for QA histos
- TParticle *fHeavyQuark[50]; // store pointer of heavy flavour quark
+ TParticle *fHeavyQuark[50]; //! store pointer of heavy flavour quark
Int_t fIsHeavy[2]; // count of heavy flavour
Int_t fNparents; // number of heavy hadrons to be considered
Int_t fParentSelect[2][7]; // heavy hadron species
}
AliDebug(2, "Particlae selected by detector");
if(fVarManager && cont){
- Char_t reccontname[256];
- sprintf(reccontname, "%sReco", contname);
- AliDebug(2, Form("Filling container %s", reccontname));
+ TString reccontname = contname; reccontname += "Reco";
+ AliDebug(2, Form("Filling container %s", reccontname.Data()));
if(fVarManager->IsSignalTrack())
- fVarManager->FillContainerStepname(cont, reccontname, SortedDetectorName(idet));
+ fVarManager->FillContainerStepname(cont, reccontname.Data(), SortedDetectorName(idet));
if(HasMCData()){
- Char_t mccontname[256];
- sprintf(mccontname, "%sMC", contname);
- AliDebug(2, Form("MC Information available, Filling container %s", mccontname));
- if(fVarManager->IsSignalTrack())
- fVarManager->FillContainerStepname(cont, mccontname, SortedDetectorName(idet), kTRUE);
- if(cont->GetCorrelationMatrix("correlationstepafterTOF")){
- TString tstept("TOFPID");
- if(!tstept.CompareTo(SortedDetectorName(idet))) {
- fVarManager->FillCorrelationMatrix(cont->GetCorrelationMatrix("correlationstepafterTOF"));
- //printf("Step %s\n",(const char*) SortedDetectorName(idet));
+ TString mccontname = contname; mccontname += "MC";
+ AliDebug(2, Form("MC Information available, Filling container %s", mccontname.Data()));
+ if(fVarManager->IsSignalTrack()) {
+ fVarManager->FillContainerStepname(cont, mccontname.Data(), SortedDetectorName(idet), kTRUE);
+ if(cont->GetCorrelationMatrix("correlationstepafterTOF")){
+ TString tstept("TOFPID");
+ if(!tstept.CompareTo(SortedDetectorName(idet))) {
+ fVarManager->FillCorrelationMatrix(cont->GetCorrelationMatrix("correlationstepafterTOF"));
+ //printf("Step %s\n",(const char*) SortedDetectorName(idet));
+ }
}
}
}
// TPC alone, symmetric 3 sigma cut and asymmetric sigma cut in the momentum region between 2GeV/c and 10 GeV/c and sigma between -1 and 100
//
AliHFEpidTPC *pid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
- pid->SetTPCnSigma(3);
- pid->SetAsymmetricTPCsigmaCut(pmin, pmax, sigmamin, sigmamax);
+ if(pid){
+ pid->SetTPCnSigma(3);
+ pid->SetAsymmetricTPCsigmaCut(pmin, pmax, sigmamin, sigmamax);
+ }
}
//____________________________________________________________
// TPC alone, symmetric 3 sigma cut and 2 - -100 sigma pion rejection
//
AliHFEpidTPC *pid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
- pid->SetTPCnSigma(3);
- pid->SetRejectParticle(AliPID::kPion, 0., -100., 10., 1.);
+ if(pid){
+ pid->SetTPCnSigma(3);
+ pid->SetRejectParticle(AliPID::kPion, 0., -100., 10., 1.);
+ }
}
//____________________________________________________________
if(HasMCData()) printf("Configuring TPC for MC\n");
AliHFEpidTPC *tpcpid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
AliHFEpidTOF *tofpid = dynamic_cast<AliHFEpidTOF *>(fDetectorPID[kTOFpid]);
- tpcpid->SetTPCnSigma(2);
- tofpid->SetTOFnSigma(3);
+ if(tofpid) tofpid->SetTOFnSigma(3);
//TF1 *upperCut = new TF1("upperCut", "[0] * TMath::Exp([1]*x)", 0, 20);
TF1 *upperCut = new TF1("upperCut", "[0]", 0, 20); // Use constant upper cut
//upperCut->SetParameter(1, -0.4357);
if(HasMCData()) lowerCut->SetParameter(0, -2.5);
- else lowerCut->SetParameter(0, -3.7);
+ else lowerCut->SetParameter(0, -3.71769);
+ //else lowerCut->SetParameter(0, -3.7);
- lowerCut->SetParameter(1, -0.8);
+ lowerCut->SetParameter(1, -0.40263);
+ //lowerCut->SetParameter(1, -0.8);
if(HasMCData()) lowerCut->SetParameter(2, -2.2);
- else lowerCut->SetParameter(2, -0.35);
+ else lowerCut->SetParameter(2, 0.267857);
+ //else lowerCut->SetParameter(2, -0.35);
- tpcpid->SetUpperSigmaCut(upperCut);
- tpcpid->SetLowerSigmaCut(lowerCut);
+ if(tpcpid){
+ tpcpid->SetTPCnSigma(2);
+ tpcpid->SetUpperSigmaCut(upperCut);
+ tpcpid->SetLowerSigmaCut(lowerCut);
+ }
AddCommonObject(upperCut);
AddCommonObject(lowerCut);
}
// TPC alone, symmetric 3 sigma cut and 2 - -100 sigma pion rejection
//
AliHFEpidTPC *pid = dynamic_cast<AliHFEpidTPC *>(fDetectorPID[kTPCpid]);
- pid->SetTPCnSigma(2);
- pid->SetRejectParticle(AliPID::kProton, 0., -3., 10., 3.);
- pid->SetRejectParticle(AliPID::kKaon, 0., -3., 10., 3.);
+ if(pid){
+ pid->SetTPCnSigma(2);
+ pid->SetRejectParticle(AliPID::kProton, 0., -3., 10., 3.);
+ pid->SetRejectParticle(AliPID::kKaon, 0., -3., 10., 3.);
+ }
}
//____________________________________________________________
return 0.;
}
AliESDtrack *track = dynamic_cast<AliESDtrack *>(vtrack);
+ if(!track) return 0.;
Double_t signal = 0.;
#ifdef TRUNK
Double_t dedx[4];
return 0.;
}
AliESDtrack *track = dynamic_cast<AliESDtrack *>(vtrack);
+ if(!track) return 0.;
Double_t dedx[4], tmp[4];
Int_t indices[4];
track->GetITSdEdxSamples(tmp);
if(!TString(o->IsA()->GetName()).CompareTo("AliESDtrack")){
// work on local copy in order to not spoil others
esdtrack = new AliESDtrack(*(dynamic_cast<AliESDtrack *>(o)));
+ if(!esdtrack) continue;
} else if(!TString(o->IsA()->GetName()).CompareTo("AliAODrack")){
// Bad hack: Fill ESD track with AOD information
esdtrack = new AliESDtrack(dynamic_cast<AliAODTrack *>(o));
+ if(!esdtrack) continue;
} else {
// Non usable
continue;
if(!TString(mcTrack->IsA()->GetName()).CompareTo("AliMCParticle")){
// case ESD
AliMCParticle *mcp = dynamic_cast<AliMCParticle *>(mcTrack);
+ if(!mcp) continue;
trackPdg = TMath::Abs(mcp->Particle()->GetPdgCode());
} else {
// case AOD
AliAODMCParticle *aodmcp = dynamic_cast<AliAODMCParticle *>(mcTrack);
+ if(!aodmcp) continue;
trackPdg = TMath::Abs(aodmcp->GetPdgCode());
}
if(trackPdg == pdg) // Correct identification
if(!TString(recTrack->IsA()->GetName()).CompareTo("AliESDtrack")){
// case ESD
AliESDtrack *esdTrack = dynamic_cast<AliESDtrack *>(recTrack);
+ if(!esdTrack) continue;
status = esdTrack->GetStatus();
//TPC momentum and likelihoods
// 8) TRD Ntrk, 9) TRD Ncls, 10) TRD dEdx,
Double_t data[12];
- memset(data, -99, sizeof(Double_t) *12);
+
Int_t run = fEvent->GetRunNumber();
AliVParticle *recTrack = NULL;
TIter trackIter(particles);
while((recTrack = dynamic_cast<AliVParticle *>(trackIter()))){
- memset(data, -99, sizeof(Double_t) *10);
+ for(Int_t i=0; i<12; ++i) data[i] = -99.;
// ESD
if(!TString(recTrack->IsA()->GetName()).CompareTo("AliESDtrack")){
// case ESD
//
if(status & AliESDtrack::kTOFpid){
Double_t p = esdTrack->GetOuterParam() ? esdTrack->GetOuterParam()->P() : esdTrack->P();
- Double_t t0 = fESDpid->GetTOFResponse().GetTimeZero();
+ Double_t t0 = fESDpid->GetTOFResponse().GetStartTime(esdTrack->P());
//TOF beta
sprintf(hname, "hTOF_beta_%s", typeName[species]);
if((esd = dynamic_cast<AliESDEvent *>(fEvent))){
AliESDtrack *track = NULL, *partnerTrack = NULL;
while((hfetrack = dynamic_cast<AliHFEV0info *>(candidates()))){
+ if(!hfetrack) continue;
track = dynamic_cast<AliESDtrack *>(hfetrack->GetTrack());
+ if(!track) continue;
partnerTrack = esd->GetTrack(hfetrack->GetPartnerID());
Double_t nSigmaTrack = TMath::Abs(fESDpid->NumberOfSigmasTPC(track, AliPID::kElectron) - shift);
Double_t nSigmaPartner = TMath::Abs(fESDpid->NumberOfSigmasTPC(partnerTrack, AliPID::kElectron) - shift);
}
} else {
aod = dynamic_cast<AliAODEvent *>(fEvent);
- AliAODTrack *track = NULL, *partnerTrack = NULL;
+ if(!aod) return NULL;
+ //AliAODTrack *track = NULL, *partnerTrack = NULL;
while((hfetrack = dynamic_cast<AliHFEV0info *>(candidates()))){
- track = dynamic_cast<AliAODTrack *>(hfetrack->GetTrack());
- partnerTrack = aod->GetTrack(hfetrack->GetPartnerID());
+ if(!hfetrack) continue;
+ //track = dynamic_cast<AliAODTrack *>(hfetrack->GetTrack());
+ //partnerTrack = aod->GetTrack(hfetrack->GetPartnerID());
// will be coming soon
}
}
if(p < 0) return kFALSE;
Int_t mombin = TRDmomBin(p); // momentum bin
+ if(mombin < 0) return -1.;
Float_t dEdxTRDsum = 0; // dEdxTRDsum for checking if tracklet is available
Float_t dEdxTRD[8]; // dEdx for a tracklet in the ESD slices
Double_t ddEdxTRD[8]; // dEdx as Double_t for TMultiLayerPerceptron::Evaluate()
fDetPID[idet] = pid->GetDetPID(static_cast<AliHFEpid::EDETtype_t>(idet));
if(pid->HasDetector(static_cast<AliHFEpid::EDETtype_t>(idet))){
CreateDetPIDqa(static_cast<AliHFEpid::EDETtype_t>(idet));
- fDetPIDqa[idet]->Initialize();
fDetPIDqa[idet]->SetPIDqaManager(this);
+ fDetPIDqa[idet]->Initialize();
}
}
}
AliHFEdetPIDqa *GetDetectorPIDqa(AliHFEpid::EDETtype_t detector) const { return fDetPIDqa[detector]; }
AliHFEpidBase *GetDetectorPID(AliHFEpid::EDETtype_t detector) const { return fDetPID[detector]; }
TList *MakeList(const Char_t *name);
+ void SetHighResolutionHistos() { SetBit(kHighResolutionHistos, kTRUE); };
+ Bool_t HasHighResolutionHistos() const { return TestBit(kHighResolutionHistos); }
protected:
enum{
- kIsOwner = BIT(14)
+ kIsOwner = BIT(14),
+ kHighResolutionHistos = BIT(15)
};
Bool_t IsOwner() const { return TestBit(kIsOwner); }
void SetOwner() { SetBit(kIsOwner, kTRUE); }
AliDebug(2, "PID object available");
AliHFEpidObject::AnalysisType_t anaType = track->IsESDanalysis() ? AliHFEpidObject::kESDanalysis : AliHFEpidObject::kAODanalysis;
- if(!((dynamic_cast<const AliVTrack *>(track->GetRecTrack()))->GetStatus() & AliESDtrack::kTOFpid)) return 0;
+ const AliVTrack *vtrack = dynamic_cast<const AliVTrack *>(track->GetRecTrack());
+ if(!(vtrack && vtrack->GetStatus() & AliESDtrack::kTOFpid)) return 0;
AliDebug(2, "Track Has TOF PID");
if(pidqa) pidqa->ProcessTrack(track, AliHFEpid::kTOFpid, AliHFEdetPIDqa::kBeforePID);
if(anaType == AliHFEpidObject::kESDanalysis){
// ESD analysis
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- if(esdtrack && fESDpid) nSigmas = fESDpid->NumberOfSigmasTOF(esdtrack, species, fESDpid->GetTOFResponse().GetTimeZero());
+ if(esdtrack && fESDpid) nSigmas = fESDpid->NumberOfSigmasTOF(esdtrack, species, fESDpid->GetTOFResponse().GetStartTime(track->P()));
} else {
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
if(aodtrack && fAODpid) nSigmas = fAODpid->NumberOfSigmasTOF(aodtrack, species);
if(anatype == AliHFEpidObject::kESDanalysis){
// ESD analysis
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- tofSignal = esdtrack->GetTOFsignal();
+ if(esdtrack) tofSignal = esdtrack->GetTOFsignal();
} else {
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
- tofSignal = aodtrack->GetDetPid() ? aodtrack->GetDetPid()->GetTOFsignal() : 0.;
+ if(aodtrack) tofSignal = aodtrack->GetDetPid() ? aodtrack->GetDetPid()->GetTOFsignal() : 0.;
}
return tofSignal;
}
if(anatype == AliHFEpidObject::kESDanalysis){
// ESD analysis
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- esdtrack->GetIntegratedTimes(times);
+ if(esdtrack) esdtrack->GetIntegratedTimes(times);
} else {
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
- aodtrack->GetDetPid()->GetIntegratedTimes(times);
+ if(aodtrack) aodtrack->GetDetPid()->GetIntegratedTimes(times);
}
}
#ifndef ALIHFEPIDBASE_H
#include "AliHFEpidBase.h"
#endif
-#include "AliPID.h"
class AliVParticle;
+class AliPID;
+
class AliHFEpidQAmanager;
class AliHFEpidTOF : public AliHFEpidBase{
, fLineCrossingsEnabled(0)
, fUpperSigmaCut(NULL)
, fLowerSigmaCut(NULL)
+ , fElectronMeanCorrection(NULL)
, fNsigmaTPC(3)
, fRejectionEnabled(0)
, fPID(NULL)
//
// default constructor
//
+ memset(fRejection, 0, sizeof(Float_t) * 4 * AliPID::kSPECIES);
+ memset(fLineCrossingSigma, 0, sizeof(Double_t) * AliPID::kSPECIES);
+ memset(fPAsigCut, 0, sizeof(Float_t) * 2);
+ memset(fNAsigmaTPC, 0, sizeof(Float_t) * 2);
+
}
//___________________________________________________________________
, fLineCrossingsEnabled(0)
, fUpperSigmaCut(NULL)
, fLowerSigmaCut(NULL)
+ , fElectronMeanCorrection(NULL)
, fNsigmaTPC(3)
, fRejectionEnabled(0)
, fPID(NULL)
//
// default constructor
//
+ //
+ memset(fRejection, 0, sizeof(Float_t) * 4 * AliPID::kSPECIES);
memset(fLineCrossingSigma, 0, sizeof(Double_t) * AliPID::kSPECIES);
memset(fPAsigCut, 0, sizeof(Float_t) * 2);
memset(fNAsigmaTPC, 0, sizeof(Float_t) * 2);
, fLineCrossingsEnabled(0)
, fUpperSigmaCut(NULL)
, fLowerSigmaCut(NULL)
+ , fElectronMeanCorrection(NULL)
, fNsigmaTPC(2)
, fRejectionEnabled(0)
, fPID(NULL)
target.fLineCrossingsEnabled = fLineCrossingsEnabled;
target.fUpperSigmaCut = fUpperSigmaCut;
target.fLowerSigmaCut = fLowerSigmaCut;
+ target.fElectronMeanCorrection = fElectronMeanCorrection;
target.fNsigmaTPC = fNsigmaTPC;
target.fRejectionEnabled = fRejectionEnabled;
target.fPID = new AliPID(*fPID);
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
if(aodtrack && fAODpid) nSigmas = fAODpid->NumberOfSigmasTPC(aodtrack, species);
}
+ // Correct for the mean o
+ if(fElectronMeanCorrection)
+ nSigmas -= fElectronMeanCorrection->Eval(GetP(track, anaType));
return nSigmas;
}
if(anatype == AliHFEpidObject::kESDanalysis){
// ESD analysis: Use Inner Params for the momentum estimate
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- p = esdtrack->GetInnerParam() ? esdtrack->GetInnerParam()->GetP() : esdtrack->P();
+ if(esdtrack) p = esdtrack->GetInnerParam() ? esdtrack->GetInnerParam()->GetP() : esdtrack->P();
} else {
// AOD analysis: Use TPC momentum stored in the AliAODpid object
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
- p = aodtrack->GetDetPid() ? aodtrack->GetDetPid()->GetTPCmomentum() : aodtrack->P();
+ if(aodtrack) p = aodtrack->GetDetPid() ? aodtrack->GetDetPid()->GetTPCmomentum() : aodtrack->P();
}
return p;
}
void AddTPCdEdxLineCrossing(Int_t species, Double_t sigma);
Bool_t HasAsymmetricSigmaCut() const { return TestBit(kAsymmetricSigmaCut);}
Bool_t HasParticleRejection() const { return TestBit(kRejection); }
+ void SetElectronMeanCorrection(TF1 *electronLineCorrection) { fElectronMeanCorrection = electronLineCorrection; }
void SetTPCnSigma(Short_t nSigma) { fNsigmaTPC = nSigma; };
inline void SetAsymmetricTPCsigmaCut(Float_t pmin, Float_t pmax, Float_t sigmaMin, Float_t sigmaMax);
inline void SetRejectParticle(Int_t species, Float_t pmin, Float_t sigmaMin, Float_t pmax, Float_t sigmaMax);
UChar_t fLineCrossingsEnabled; // Bitmap showing which line crossing is set
TF1 *fUpperSigmaCut; // Upper Sigma Cut
TF1 *fLowerSigmaCut; // Lower Sigma Cut
+ TF1 *fElectronMeanCorrection; // Correct the mean of the electron line position as function of the momentum
Float_t fPAsigCut[2]; // Momentum region where to perform asymmetric sigma cut
Float_t fNAsigmaTPC[2]; // Asymmetric TPC Sigma band
Short_t fNsigmaTPC; // TPC sigma band
//
// default constructor
//
+ memset(fThreshParams, 0, sizeof(Double_t) * kThreshParams);
+ SetUseDefaultParameters();
}
//___________________________________________________________________
// default constructor
//
memset(fThreshParams, 0, sizeof(Double_t) * kThreshParams);
+ SetUseDefaultParameters();
}
//___________________________________________________________________
// Performs the copying of the object
//
AliHFEpidTRD &target = dynamic_cast<AliHFEpidTRD &>(ref);
-
+
+ Bool_t defaultParameters = UseDefaultParameters();
+ target.SetUseDefaultParameters(defaultParameters);
target.fMinP = fMinP;
target.fPIDMethod = fPIDMethod;
target.fElectronEfficiency = fElectronEfficiency;
// InitializePID: Load TRD thresholds and create the electron efficiency axis
// to navigate
//
- if(fPIDMethod == kLQ)
- InitParameters1DLQ();
- else
- InitParameters();
+ if(UseDefaultParameters()){
+ if(fPIDMethod == kLQ)
+ InitParameters1DLQ();
+ else
+ InitParameters();
+ }
return kTRUE;
}
return TMath::Max(TMath::Min(threshold, 0.99), 0.2); // truncate the threshold upperwards to 0.999 and lowerwards to 0.2 and exclude unphysical values
}
+//___________________________________________________________________
+void AliHFEpidTRD::SetThresholdParameters(Double_t electronEff, Double_t *params){
+ //
+ // Set threshold parameters for the given bin
+ //
+ if(electronEff >= 1. || electronEff < 0.7) return;
+ Int_t effbin = static_cast<Int_t>((electronEff - 0.7)/0.05);
+ memcpy(&fThreshParams[effbin * 4], params, sizeof(Double_t) * 4);
+ SetUseDefaultParameters(kFALSE);
+}
+
//___________________________________________________________________
void AliHFEpidTRD::InitParameters(){
//
Double_t pidProbs[AliPID::kSPECIES]; memset(pidProbs, 0, sizeof(Double_t) * AliPID::kSPECIES);
if(anaType == AliHFEpidObject::kESDanalysis){
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- esdtrack->GetTRDpid(pidProbs);
+ if(esdtrack) esdtrack->GetTRDpid(pidProbs);
} else {
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
- fAODpid->MakeTRDPID(const_cast<AliAODTrack *>(aodtrack), pidProbs);
+
if(aodtrack)fAODpid->MakeTRDPID(const_cast<AliAODTrack *>(aodtrack), pidProbs);
}
return pidProbs[AliPID::kElectron];
}
Double_t p = 0.;
if(anaType == AliHFEpidObject::kESDanalysis){
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- p = esdtrack->GetOuterParam() ? esdtrack->GetOuterParam()->P() : esdtrack->P();
+ if(esdtrack) p = esdtrack->GetOuterParam() ? esdtrack->GetOuterParam()->P() : esdtrack->P();
} else {
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
- p = aodtrack->P();
+ if(aodtrack) p = aodtrack->P();
}
return p;
}
Double_t charge = 0.;
if(anaType == AliHFEpidObject::kESDanalysis){
const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
- for(Int_t islice = 0; islice < esdtrack->GetNumberOfTRDslices(); islice++) charge += esdtrack->GetTRDslice(static_cast<UInt_t>(layer), islice);
+ if(esdtrack)
+ for(Int_t islice = 0; islice < esdtrack->GetNumberOfTRDslices(); islice++) charge += esdtrack->GetTRDslice(static_cast<UInt_t>(layer), islice);
} else {
const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
- AliAODPid *aoddetpid = aodtrack->GetDetPid();
- for(Int_t islice = 0; islice < aoddetpid->GetTRDnSlices(); islice++) charge += aoddetpid->GetTRDsignal()[layer * aoddetpid->GetTRDnSlices() + islice];
+ AliAODPid *aoddetpid = aodtrack ? aodtrack->GetDetPid() : NULL;
+ if(aoddetpid)
+ for(Int_t islice = 0; islice < aoddetpid->GetTRDnSlices(); islice++) charge += aoddetpid->GetTRDsignal()[layer * aoddetpid->GetTRDnSlices() + islice];
}
return charge;
}
//___________________________________________________________________
-Double_t AliHFEpidTRD::GetTRDSignalV1(const AliESDtrack *track) const {
+void AliHFEpidTRD::GetTRDmomenta(const AliVParticle *track, AliHFEpidObject::AnalysisType_t anaType, Double_t *mom) const {
+ //
+ // Fill Array with momentum information at the TRD tracklet
+ //
+ if(anaType == AliHFEpidObject::kESDanalysis){
+ const AliESDtrack *esdtrack = dynamic_cast<const AliESDtrack *>(track);
+ if(esdtrack)
+ for(Int_t itl = 0; itl < 6; itl++)
+ mom[itl] = esdtrack->GetTRDmomentum(itl);
+ } else {
+ const AliAODTrack *aodtrack = dynamic_cast<const AliAODTrack *>(track);
+ AliAODPid *aoddetpid = aodtrack ? aodtrack->GetDetPid() : NULL;
+ if(aoddetpid){
+ Float_t *trdmom = aoddetpid->GetTRDmomentum();
+ for(Int_t itl = 0; itl < 6; itl++){
+ mom[itl] = trdmom[itl];
+ }
+ }
+ }
+}
+
+//___________________________________________________________________
+Double_t AliHFEpidTRD::GetTRDSignalV1(const AliESDtrack *track, Float_t truncation) const {
//
// Calculation of the TRD Signal via truncated mean
// Method 1: Take all Slices available
// Calculate mean over the last 2/3 slices
//
const Int_t kNSlices = 48;
- AliDebug(3, Form("Number of Tracklets: %d\n", track->GetTRDpidQuality()));
+ const Int_t kSlicePerLayer = 7;
+ // Weight the slice to equalize the MPV of the dQ/dl-distribution per slice to the one in the first slice
+ // Pions are used as reference for the equalization
+ const Double_t kWeightSlice[8] = {1., 2.122, 1.8, 1.635, 1.595, 1.614, 1.16, 7.0};
+ AliDebug(3, Form("Number of Tracklets: %d\n", track->GetTRDntrackletsPID()));
Double_t trdSlices[kNSlices], tmp[kNSlices];
Int_t indices[48];
Int_t icnt = 0;
for(Int_t idet = 0; idet < 6; idet++)
- for(Int_t islice = 0; islice < 8; islice++){
+ for(Int_t islice = 0; islice < kSlicePerLayer; islice++){
AliDebug(2, Form("Chamber[%d], Slice[%d]: TRDSlice = %f", idet, islice, track->GetTRDslice(idet, islice)));
if(TMath::Abs(track->GetTRDslice(idet, islice)) < fgkVerySmall) continue;;
- trdSlices[icnt++] = track->GetTRDslice(idet, islice);
+ trdSlices[icnt++] = track->GetTRDslice(idet, islice) * kWeightSlice[islice];
}
AliDebug(1, Form("Number of Slices: %d\n", icnt));
if(icnt < 6) return 0.; // We need at least 6 Slices for the truncated mean
memcpy(tmp, trdSlices, sizeof(Double_t) * icnt);
for(Int_t ien = 0; ien < icnt; ien++)
trdSlices[ien] = tmp[indices[ien]];
- Double_t trdSignal = TMath::Mean(static_cast<Int_t>(static_cast<Double_t>(icnt) * 2./3.), trdSlices);
+ Double_t trdSignal = TMath::Mean(static_cast<Int_t>(static_cast<Float_t>(icnt) * truncation), trdSlices);
Double_t mom = track->GetOuterParam() ? track->GetOuterParam()->P() : -1;
AliDebug(3, Form("PID Meth. 1: p[%f], TRDSignal[%f]", mom, trdSignal));
return trdSignal;
}
//___________________________________________________________________
-Double_t AliHFEpidTRD::GetTRDSignalV2(const AliESDtrack *track) const {
+Double_t AliHFEpidTRD::GetTRDSignalV2(const AliESDtrack *track, Float_t truncation) const {
//
// Calculation of the TRD Signal via truncated mean
// Method 2: Take only first 5 slices per chamber
// Cut out upper half
// Now do mean with the reamining 3 slices per chamber
//
- Double_t trdSlicesLowTime[30], trdSlicesRemaining[32];
- Int_t indices[30];
+ const Double_t kWeightSlice[8] = {1., 2.122, 1.8, 1.635, 1.595, 1.614, 1.16, 7.0};
+ const Int_t kLayers = 6;
+ const Int_t kSlicesLow = 6;
+ const Int_t kSlicesHigh = 1;
+ Double_t trdSlicesLowTime[kLayers*kSlicesLow], trdSlicesRemaining[kLayers*(kSlicesHigh + kSlicesLow)];
+ Int_t indices[kLayers*kSlicesLow];
Int_t cntLowTime=0, cntRemaining = 0;
for(Int_t idet = 0; idet < 6; idet++)
- for(Int_t islice = 0; islice < 8; islice++){
+ for(Int_t islice = 0; islice < kSlicesLow+kSlicesHigh; islice++){
if(TMath::Abs(track->GetTRDslice(idet, islice)) < fgkVerySmall) continue;;
- if(islice < 5){
+ if(islice < kSlicesLow){
AliDebug(3, Form("Part 1, Det[%d], Slice[%d], TRDSlice: %f", idet, islice, track->GetTRDslice(idet, islice)));
- trdSlicesLowTime[cntLowTime++] = track->GetTRDslice(idet, islice);
+ trdSlicesLowTime[cntLowTime++] = track->GetTRDslice(idet, islice) * kWeightSlice[islice];
} else{
AliDebug(3, Form("Part 1, Det[%d], Slice[%d], TRDSlice: %f", idet, islice, track->GetTRDslice(idet, islice)));
- trdSlicesRemaining[cntRemaining++] = track->GetTRDslice(idet, islice);
+ trdSlicesRemaining[cntRemaining++] = track->GetTRDslice(idet, islice) * kWeightSlice[islice];
}
}
if(cntLowTime < 4 || cntRemaining < 2) return 0.; // Min. Number of Slices at high time is 2 (matches with 1 layer), for the truncated mean we need at least 4 Slices
TMath::Sort(cntLowTime, trdSlicesLowTime, indices, kFALSE);
// Fill the second array with the lower half of the first time bins
- for(Int_t ien = 0; ien < static_cast<Int_t>(static_cast<Double_t>(cntLowTime) * 0.5); ien++)
+ for(Int_t ien = 0; ien < static_cast<Int_t>(static_cast<Float_t>(cntLowTime) * truncation); ien++)
trdSlicesRemaining[cntRemaining++] = trdSlicesLowTime[indices[ien]];
Double_t trdSignal = TMath::Mean(cntRemaining, trdSlicesRemaining);
Double_t mom = track->GetOuterParam() ? track->GetOuterParam()->P() : -1;
AliHFEpidTRD(const AliHFEpidTRD &ref);
AliHFEpidTRD& operator=(const AliHFEpidTRD &ref);
virtual ~AliHFEpidTRD();
-
+
virtual Bool_t InitializePID();
virtual Int_t IsSelected(const AliHFEpidObject *track, AliHFEpidQAmanager *pidqa) const;
- Double_t GetTRDSignalV1(const AliESDtrack *track) const;
- Double_t GetTRDSignalV2(const AliESDtrack *track) const;
+ Double_t GetTRDSignalV1(const AliESDtrack *track, Float_t truncation = 0.7) const;
+ Double_t GetTRDSignalV2(const AliESDtrack *track, Float_t trucation = 0.7) const;
Bool_t IsCalculateTRDSignals() const { return TestBit(kTRDsignals); }
void SetPIDMethod(PIDMethodTRD_t method) { fPIDMethod = method; };
void SetElectronEfficiency(Double_t electronEfficiency) { fElectronEfficiency = electronEfficiency; }
+ void SetThresholdParameters(Double_t electronEff, Double_t *params);
void SetMinP(Double_t p) { fMinP = p; }
void CalculateTRDSignals(Bool_t docalc) { SetBit(kTRDsignals, docalc); }
Double_t GetElectronLikelihood(const AliVParticle *track, AliHFEpidObject::AnalysisType_t anaType) const;
+ void GetTRDmomenta(const AliVParticle *track, AliHFEpidObject::AnalysisType_t anaType, Double_t *mom) const;
Double_t GetP(const AliVParticle *track, AliHFEpidObject::AnalysisType_t anaType) const;
Double_t GetTRDthresholds(Double_t electronEff, Double_t p) const;
Double_t GetChargeLayer(const AliVParticle *track, UInt_t layer, AliHFEpidObject::AnalysisType_t anatype) const;
protected:
enum{
- kTRDsignals = BIT(16)
+ kTRDsignals = BIT(16),
+ kTRDdefaultThresholds = BIT(17)
};
void Copy(TObject &ref) const;
void InitParameters();
void InitParameters1DLQ();
void GetParameters(Double_t electronEff, Double_t *parameters) const;
+ void SetUseDefaultParameters(Bool_t useDefault = kTRUE) { SetBit(kTRDdefaultThresholds, useDefault); }
+ Bool_t UseDefaultParameters() const { return TestBit(kTRDdefaultThresholds); }
private:
static const Double_t fgkVerySmall; // Check for 0
if(charge) fSignalToBackgroundMC->GetAxis(3)->SetRange(0, fSignalToBackgroundMC->GetAxis(3)->GetLast() + 1);
TH1 *hEff = dynamic_cast<TH1D *>(hNom->Clone());
- char hname[256];
- sprintf(hname, mode ? "sigToBack" : "sigEff");
- char cname[256];
- Color_t mycolor = kBlack;
- switch(charge){
- case 0: mycolor = kBlue; sprintf(cname, "All"); break;
- case 1: mycolor = kBlack; sprintf(cname, "Neg"); break;
- case 2: mycolor = kRed; sprintf(cname, "Pos"); break;
- default: break;
- }
- sprintf(hname, "%s%s", hname, cname);
- hEff->SetName(hname);
- hEff->SetTitle(mode ? "Signal/Background" : "Signal/(Signal+Background)");
- hEff->Divide(hDenom);
+ if(hEff){
+ TString hname, cname;
+ hname = mode ? "sigToBack" : "sigEff";
+ Color_t mycolor = kBlack;
+ switch(charge){
+ case 0: mycolor = kBlue; cname = "All"; break;
+ case 1: mycolor = kBlack; cname = "Neg"; break;
+ case 2: mycolor = kRed; cname ="Pos"; break;
+ default: break;
+ }
+ hname += cname;
+ hEff->SetName(hname);
+ hEff->SetTitle(mode ? "Signal/Background" : "Signal/(Signal+Background)");
+ hEff->Divide(hDenom);
- // Make nice plots
- hEff->GetXaxis()->SetTitle("p_{T} / GeV/c");
- hEff->GetYaxis()->SetTitle("Efficiency");
- hEff->SetStats(kFALSE);
- hEff->SetLineColor(kBlack);
- hEff->SetLineWidth(1);
- hEff->SetMarkerStyle(22);
- hEff->SetMarkerColor(mycolor);
+ // Make nice plots
+ hEff->GetXaxis()->SetTitle("p_{T} / GeV/c");
+ hEff->GetYaxis()->SetTitle("Efficiency");
+ hEff->SetStats(kFALSE);
+ hEff->SetLineColor(kBlack);
+ hEff->SetLineWidth(1);
+ hEff->SetMarkerStyle(22);
+ hEff->SetMarkerColor(mycolor);
+ }
delete hNom; delete hDenom;
return hEff;
fPIDperformance->GetAxis(3)->SetRange(0, fPIDperformance->GetAxis(3)->GetNbins() + 1);
TH1 *hNom = NULL, *hDenom = NULL;
- char hname[256], htitle[256], cname[256];
+ TString hname, htitle, cname;
Color_t mycolor = kBlack;
if(charge) fPIDperformance->GetAxis(3)->SetRange(charge, charge);
// Normalisation by all candidates - no restriction in axis 4 - only for mode == 1
switch(mode){
case 0: // Electron purity
fPIDperformance->GetAxis(4)->SetRange(2,3);
- sprintf(hname, "electronPurity");
- sprintf(htitle, "Electron Purity");
+ hname = "electronPurity";
+ htitle = "Electron Purity";
break;
case 1: // Signal purity
fPIDperformance->GetAxis(4)->SetRange(3,3); // here signal not divided into charm and beauty
- sprintf(hname, "signalPurity");
- sprintf(htitle, "Signal Purity");
+ hname = "signalPurity";
+ htitle = "Signal Purity";
break;
case 2: // Fake contamination
fPIDperformance->GetAxis(4)->SetRange(1,1);
- sprintf(hname, "fakeContamination");
- sprintf(htitle, "Contamination of misidentified hadrons");
+ hname = "fakeContamination";
+ htitle = "Contamination of misidentified hadrons";
break;
default: break;
}
switch(charge){
case 0:
- sprintf(cname, "All");
+ cname = "All";
mycolor = kBlue;
break;
case 1:
- sprintf(cname, "Neg");
+ cname = "Neg";
mycolor = kBlack;
break;
case 2:
- sprintf(cname, "Pos");
+ cname = "Pos";
mycolor = kRed;
break;
}
- sprintf(hname, "%s%s", hname, cname);
+ hname += cname;
hNom = fPIDperformance->Projection(0);
hNom->Sumw2();
hNom->SetName("hNom");
// Create Efficiency histogram
TH1 *hEff = dynamic_cast<TH1D *>(hNom->Clone());
- hEff->SetName(hname);
- hEff->SetTitle(htitle);
- hEff->Divide(hDenom);
- hEff->Scale(100.);
- hEff->GetXaxis()->SetTitle("p_{T} / GeV/c");
- hEff->GetYaxis()->SetTitle(mode < 2 ? "Purity / %" : "Contamination / %");
- hEff->GetYaxis()->SetRangeUser(0., 100.);
- hEff->SetStats(kFALSE);
- hEff->SetLineColor(kBlack);
- hEff->SetLineWidth(1);
- hEff->SetMarkerColor(mycolor);
- hEff->SetMarkerStyle(22);
- delete hNom; delete hDenom;
+ if(hEff){
+ hEff->SetName(hname.Data());
+ hEff->SetTitle(htitle.Data());
+ hEff->Divide(hDenom);
+ hEff->Scale(100.);
+ hEff->GetXaxis()->SetTitle("p_{T} / GeV/c");
+ hEff->GetYaxis()->SetTitle(mode < 2 ? "Purity / %" : "Contamination / %");
+ hEff->GetYaxis()->SetRangeUser(0., 100.);
+ hEff->SetStats(kFALSE);
+ hEff->SetLineColor(kBlack);
+ hEff->SetLineWidth(1);
+ hEff->SetMarkerColor(mycolor);
+ hEff->SetMarkerStyle(22);
+ delete hNom; delete hDenom;
+ }
return hEff;
}
AliESDtrack *track = dynamic_cast<AliESDtrack *>(signalTrack);
+ if(!track){
+ AliDebug(1, "no esd track pointer, return\n");
+ return;
+ }
+
FillHistos(0,track); // wo any cuts
InitHFEpairs();
//
// make container
//
-
+
step = step*7;
AliMCParticle *mctrack = NULL;
Int_t esource=fMCQA->GetElecSource(mcpart);
if(esource==1) {
+ if(!(dynamic_cast<TH1F *>(fSecVtxList->At(step+1)))) return;
(dynamic_cast<TH1F *>(fSecVtxList->At(step+1)))->Fill(mcpart->Pt()); //charm
}
else if(esource==2 || esource==3) {
+ if(!(dynamic_cast<TH1F *>(fSecVtxList->At(step+2)))) return;
(dynamic_cast<TH1F *>(fSecVtxList->At(step+2)))->Fill(mcpart->Pt()); //beauty
}
else if(esource==4) {
+ if(!(dynamic_cast<TH1F *>(fSecVtxList->At(step+3)))) return;
(dynamic_cast<TH1F *>(fSecVtxList->At(step+3)))->Fill(mcpart->Pt()); //conversion
}
else if(esource==7) {
+ if(!(dynamic_cast<TH1F *>(fSecVtxList->At(step+5)))) return;
(dynamic_cast<TH1F *>(fSecVtxList->At(step+5)))->Fill(mcpart->Pt()); //contamination
}
else if(!(esource<0)) {
+ if(!(dynamic_cast<TH1F *>(fSecVtxList->At(step+4)))) return;
(dynamic_cast<TH1F *>(fSecVtxList->At(step+4)))->Fill(mcpart->Pt()); //e backgrounds
}
else {
+ if(!(dynamic_cast<TH1F *>(fSecVtxList->At(step+6)))) return;
(dynamic_cast<TH1F *>(fSecVtxList->At(step+6)))->Fill(mcpart->Pt()); //something else?
}
}
//
// Check if mother is coming from Charm
//
+ if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
if(esources == AliHFEmcQA::kDirectCharm) // 1: direct D->e
return kTRUE;
//
// Check if mother is coming from Beauty
//
+ if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
if(esources == AliHFEmcQA::kDirectBeauty || esources == AliHFEmcQA::kBeautyCharm) // 2: B->e 3: B->D->e
return kTRUE;
//
// Check for MC if the electron is coming from Gamma
//
+ if(!dynamic_cast<const AliVParticle *>(o)) return kFALSE;
Int_t esources = GetElecSource(dynamic_cast<const AliVParticle *>(o));
if(esources == AliHFEmcQA::kGamma) // 4: conversion electrons
return kTRUE;
Int_t motherPDG = 0;
if(TString(mctrack->IsA()->GetName()).CompareTo("AliMCParticle") == 0){
// case MC Particle
- motherParticle = fMC->GetTrack((dynamic_cast<const AliMCParticle *>(mctrack)->Particle()->GetFirstMother()));
- if(motherParticle)
- motherPDG = TMath::Abs((dynamic_cast<const AliMCParticle *>(motherParticle))->Particle()->GetPdgCode());
+ const AliMCParticle *esdmctrack = dynamic_cast<const AliMCParticle *>(mctrack);
+ if(esdmctrack) motherParticle = fMC->GetTrack(esdmctrack->Particle()->GetFirstMother());
+ if(motherParticle){
+ const AliMCParticle *esdmcmother = dynamic_cast<const AliMCParticle *>(motherParticle);
+ if(esdmcmother) motherPDG = TMath::Abs(esdmcmother->Particle()->GetPdgCode());
+ }
} else {
// case AODMCParticle
- motherParticle = fMC->GetTrack((dynamic_cast<const AliAODMCParticle *>(mctrack))->GetMother());
- if(motherParticle)
- motherPDG = TMath::Abs((dynamic_cast<const AliAODMCParticle *>(motherParticle))->GetPdgCode());
+ const AliAODMCParticle *aodmctrack = dynamic_cast<const AliAODMCParticle *>(mctrack);
+ if(aodmctrack) motherParticle = fMC->GetTrack(aodmctrack->GetMother());
+ if(motherParticle){
+ const AliAODMCParticle *aodmcmother = dynamic_cast<const AliAODMCParticle *>(motherParticle);
+ if(aodmcmother) motherPDG = TMath::Abs(aodmcmother->GetPdgCode());
+ }
}
return motherPDG;
}
Int_t trackPdg = 0;
if(!mctype.CompareTo("AliMCParticle")){
const AliMCParticle *esdmc = dynamic_cast<const AliMCParticle *>(mctrack);
- trackPdg = esdmc->Particle()->GetPdgCode();
+ if(esdmc) trackPdg = esdmc->Particle()->GetPdgCode();
} else {
const AliAODMCParticle *aodmc = dynamic_cast< const AliAODMCParticle *>(mctrack);
- trackPdg = aodmc->GetPdgCode();
+ if(aodmc) trackPdg = aodmc->GetPdgCode();
}
return trackPdg;
}
AliDebug(1, "No MCQA Available\n");
return 0;
}
+ if(!track){
+ AliDebug(1, "Track not Available\n");
+ return 0;
+ }
TString sourcetype = track->IsA()->GetName();
const AliVParticle *mctrack = NULL;
// MC container: correlation matrix
THnSparseF *mccorrelation = 0x0;
- if(fInclusiveSpectrum) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterPID");
+ if(fInclusiveSpectrum) {
+ if((fStepMC==(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 2))) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterPID");
+ if((fStepMC==(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 1))) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterTOF");
+ if((fStepMC==(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD))) mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepbeforePID");
+ }
else mccorrelation = mchfecontainer->GetCorrelationMatrix("correlationstepafterDE");
if(!mccorrelation) return kFALSE;
THnSparseF *mccorrelationD = GetSlicedCorrelation(mccorrelation, fNbDimensions, dims);
if(fDebugLevel > 0) {
+ printf("Step MC: %d\n",fStepMC);
+ printf("Step tracking: %d\n",AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack);
+ printf("Step MC true: %d\n",fStepTrue);
AliCFEffGrid *efficiencymcPID = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerMC),fStepMC,AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack);
AliCFEffGrid *efficiencymctrackinggeo = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerMC),AliHFEcuts::kStepHFEcutsTRD + AliHFEcuts::kNcutStepsMCTrack,fStepTrue);
AliCFEffGrid *efficiencymcall = (AliCFEffGrid*) GetEfficiency(GetContainer(kMCContainerMC),fStepMC,fStepTrue);
void SetStepToCorrect(Int_t step) { fStepData = step; };
void SetStepBeforeCutsV0(Int_t step) { fStepBeforeCutsV0 = step; };
void SetStepAfterCutsV0(Int_t step) { fStepAfterCutsV0 = step; };
+ void SetNbDimensions(Int_t nbDimensions) { fNbDimensions = nbDimensions; };
void SetStepGuessedUnfolding(Int_t stepGuessedUnfolding) { fStepGuessedUnfolding = stepGuessedUnfolding; };
void SetNumberOfIteration(Int_t numberOfIteration) { fNumberOfIterations = numberOfIteration; };
#include "AliCFCutBase.h"
#include "AliCFContainer.h"
#include "AliCFManager.h"
+#include "AliESDtrack.h"
#include "AliLog.h"
#include "AliPID.h"
+#include "AliHFEcollection.h"
#include "AliHFEcontainer.h"
#include "AliHFEcuts.h"
#include "AliHFEpid.h"
, fPIDqa(NULL)
, fCuts(NULL)
, fCFM(NULL)
+ , fQAhistos(NULL)
, fClean(kFALSE)
{
//
, fPIDqa(ref.fPIDqa)
, fCuts(ref.fCuts)
, fCFM(ref.fCFM)
+ , fQAhistos(ref.fQAhistos)
, fClean(ref.fClean)
{
//
fPIDqa = ref.fPIDqa;
fCuts = ref.fCuts;
fCFM = ref.fCFM;
+ fQAhistos = ref.fQAhistos;
fClean = ref.fClean;
if(ref.fContainer) InitContainer();
fContainer->SetStepTitle("taggedTrackContainerReco", fPID->SortedDetectorName(ipid), ipid + AliHFEcuts::kNcutStepsRecTrack);
}
fCFM->SetParticleContainer(fContainer->GetCFContainer("taggedTrackContainerReco"));
+
+ // temporarily special QA
+ fQAhistos = new AliHFEcollection("taggedTrackQA", "Special QA for the TaggedTrackAnalysis");
+ fQAhistos->CreateTH2F("TPCclusters2_1", "TPCclusterInfo for findable clusters for 2 neighbors", 30, 0.1, 10., 162, 0., 161.);
+ fQAhistos->CreateTH2F("TPCclusters2_0", "TPCclusterInfo for the ratio for 2 neighbors", 30, 0.1, 10., 100, 0., 1.);
+ fQAhistos->BinLogAxis("TPCclusters2_1", 0); // pt axis in logarithmic binning
+ fQAhistos->BinLogAxis("TPCclusters2_0", 0); // pt axis in logarithmic binning
}
//____________________________________________________________
AliDebug(2, Form("Cut passed, filling container %d", icut - offset + 1));
fVarManager->FillContainer(fCFM->GetParticleContainer(), icut - offset + 1);
}
- if(survived){
+
+ if(survived){
AliDebug(2, "Use track in the PID");
// Try a loose cut to reject pion contamination
if(fClean) {
if(abinitioPID == AliPID::kElectron){
- AliHFEpidTPC *pidTPC = (AliHFEpidTPC *) fPID->GetDetPID(AliHFEpid::kTPCpid);
- if(pidTPC) {
- Double_t numberOfSigmaTPC = pidTPC->NumberOfSigmas(track,AliPID::kElectron,AliHFEpidObject::kESDanalysis);
- if(numberOfSigmaTPC < -5) return;
- }
+ AliHFEpidTPC *pidTPC = (AliHFEpidTPC *) fPID->GetDetPID(AliHFEpid::kTPCpid);
+ if(pidTPC) {
+ Double_t numberOfSigmaTPC = pidTPC->NumberOfSigmas(track,AliPID::kElectron,AliHFEpidObject::kESDanalysis);
+ if(numberOfSigmaTPC < -5) return;
+ }
}
}
+ // temporarily monitoring of the number of TPC clusters
+ AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
+ if(esdtrack && abinitioPID == AliPID::kElectron){
+ fQAhistos->Fill("TPCclusters2_1", track->Pt(), esdtrack->GetTPCClusterInfo(2,1));
+ fQAhistos->Fill("TPCclusters2_0", track->Pt(), esdtrack->GetTPCClusterInfo(2,0));
+ }
// Apply PID
AliHFEpidObject hfetrack;
hfetrack.SetAnalysisType(AliHFEpidObject::kESDanalysis);
fCuts->CreateStandardCuts();
fCuts->SetQAOn();
SetBit(kIsOwnerCuts);
+ } else {
+ AliDebug(1, "Setting single track cuts");
+ fCuts = cuts;
}
- AliDebug(1, "Setting single track cuts");
- fCuts = cuts;
const Int_t kNcutSteps = AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kNcutStepsRecTrack + AliHFEcuts::kNcutStepsDETrack;
printf("Setting Number of cut steps %d\n", kNcutSteps);
fCFM->SetNStepParticle(kNcutSteps);
#include <TObject.h>
#endif
+class AliHFEcollection;
class AliHFEcontainer;
class AliHFEcuts;
class AliHFEpid;
AliHFEpidQAmanager *GetPIDqa() const { return fPIDqa; }
TList * GetPIDQA() const;
TList * GetCutQA() const;
+ AliHFEcollection * GetQAcollection() const { return fQAhistos; }
Bool_t GetClean() const { return fClean; };
void SetCuts(AliHFEcuts *cuts);
AliHFEpidQAmanager *fPIDqa; // PID monitoring
AliHFEcuts *fCuts; // Single track cuts
AliCFManager *fCFM; // CF Manager used for the track filtering
+ AliHFEcollection *fQAhistos; // QA histos
Bool_t fClean; // Clean
ClassDef(AliHFEtaggedTrackAnalysis, 0)
// Make common binning
const Int_t kPIDbins = AliPID::kSPECIES + 1;
- const Int_t kPbins = 100;
const Int_t kSteps = 2;
const Double_t kMinPID = -1;
const Double_t kMinP = 0.;
const Double_t kMaxPID = (Double_t)AliPID::kSPECIES;
const Double_t kMaxP = 20.;
-
+ // Quantities where one can switch between low and high resolution
+ Int_t kPbins = fQAmanager->HasHighResolutionHistos() ? 1000 : 100;
+ Int_t kSigmaBins = fQAmanager->HasHighResolutionHistos() ? 1400 : 240;
+
// 1st histogram: TOF sigmas: (species, p nsigma, step)
- const Int_t kSigmaBins = 240;
Int_t nBinsSigma[4] = {kPIDbins, kPbins, kSigmaBins, kSteps};
Double_t minSigma[4] = {kMinPID, kMinP, -12., 0};
Double_t maxSigma[4] = {kMaxPID, kMaxP, 12., 2.};
Double_t contentSignal[4];
contentSignal[0] = species;
contentSignal[1] = track->GetRecTrack()->P();
- contentSignal[2] = tofpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype);
+ contentSignal[2] = tofpid ? tofpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype): 0.;
contentSignal[3] = step;
- (dynamic_cast<THnSparseF *>(fHistos->Get("tofnSigma")))->Fill(contentSignal);
- Double_t timeTof = tofpid->GetTOFsignal(track->GetRecTrack(), anatype);
- Double_t time0 = tofpid->GetTime0(anatype);
- Double_t tof = timeTof - time0;
- Double_t times[AliPID::kSPECIES]; tofpid->GetIntegratedTimes(track->GetRecTrack(), times, anatype);
- fHistos->Fill("tofTimeRes",contentSignal[1], tof - times[AliPID::kPion]);
+ fHistos->Fill("tofnSigma", contentSignal);
+ if(tofpid){
+ Double_t timeTof = tofpid->GetTOFsignal(track->GetRecTrack(), anatype);
+ Double_t time0 = tofpid->GetTime0(anatype);
+ Double_t tof = timeTof - time0;
+ Double_t times[AliPID::kSPECIES]; tofpid->GetIntegratedTimes(track->GetRecTrack(), times, anatype);
+ fHistos->Fill("tofTimeRes",contentSignal[1], tof - times[AliPID::kPion]);
+ }
if(species > -1){
- contentSignal[2] = tpcpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype);
+ contentSignal[2] = tpcpid ? tpcpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype) : 0.;
fHistos->Fill("tofMonitorTPC", contentSignal);
}
}
// Plot the Spectrum
//
THnSparseF *hSignal = dynamic_cast<THnSparseF *>(fHistos->Get("tofnSigma"));
+ if(!hSignal) return NULL;
hSignal->GetAxis(3)->SetRange(istep + 1, istep + 1);
if(species > 0 && species < AliPID::kSPECIES)
hSignal->GetAxis(0)->SetRange(2 + species, 2 + species);
TH2 *hTmp = hSignal->Projection(2,1);
- Char_t hname[256], htitle[256];
- sprintf(hname, "hTPCsigma%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
- sprintf(htitle, "TPC dE/dx Spectrum[#sigma] %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
+ TString hname = Form("hTPCsigma%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after"),
+ htitle = Form("TPC dE/dx Spectrum[#sigma] %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
if(species > -1){
- strncat(hname, AliPID::ParticleName(species), strlen(AliPID::ParticleName(species)));
- Char_t speciesname[256];
- sprintf(speciesname, " for %ss", AliPID::ParticleName(species));
- strncat(htitle, speciesname, strlen(speciesname));
+ hname += AliPID::ParticleName(species);
+ htitle += Form(" for %ss", AliPID::ParticleName(species));
}
- hTmp->SetName(hname);
- hTmp->SetTitle(htitle);
+ hTmp->SetName(hname.Data());
+ hTmp->SetTitle(htitle.Data());
hTmp->SetStats(kFALSE);
hTmp->GetXaxis()->SetTitle("p [GeV/c]");
hTmp->GetYaxis()->SetTitle("TOF time|_{el} - expected time|_{el} [#sigma]");
TAxis *axis = 0x0;
if(o->InheritsFrom("TH1")){
- axis = (dynamic_cast<TH1F*>(o))->GetXaxis();
+ TH1 *h1 = dynamic_cast<TH1F*>(o);
+ if(h1) axis = h1->GetXaxis();
}
else if(o->InheritsFrom("TH2")){
- if(0 == dim){
- axis = (dynamic_cast<TH2F*>(o))->GetXaxis();
+ TH2 *h2 = dynamic_cast<TH2F*>(o);
+ if(h2 && 0 == dim){
+ axis = h2->GetXaxis();
}
- else if(1 == dim){
- axis = (dynamic_cast<TH2F*>(o))->GetYaxis();
+ else if(h2 && 1 == dim){
+ axis = h2->GetYaxis();
}
else{
AliError("Only dim = 0 or 1 possible for TH2F");
}
}
else if(o->InheritsFrom("THnSparse")){
- axis = (dynamic_cast<THnSparse*>(o))->GetAxis(dim);
+ THnSparse *hs = dynamic_cast<THnSparse*>(o);
+ if(hs) axis = hs->GetAxis(dim);
}
else{
AliError("Type of input object not recognized, please check your code or update this finction");
newBins[i] = factor * newBins[i-1];
}
axis->Set(bins, newBins);
- delete newBins;
+ delete[] newBins;
return kTRUE;
}
// Author:
// Markus Fasel <M.Fasel@gsi.de>
//
+#include <TBrowser.h>
#include <TClass.h>
#include <TH2.h>
#include <THnSparse.h>
return count + 1;
}
+//_________________________________________________________
+void AliHFEtpcPIDqa::Browse(TBrowser *b){
+ //
+ // Browse the PID QA
+ //
+ if(b){
+ if(fHistos){
+ b->Add(fHistos, fHistos->GetName());
+
+ // Make Projections of the dE/dx Spectra and add them to a new Folder
+ TString specnames[4] = {"All", "Electrons", "Pions", "Protons"};
+ Int_t specind[4] = {-1, AliPID::kElectron, AliPID::kPion, AliPID::kProton};
+ TList *listdEdx = new TList;
+ listdEdx->SetOwner();
+ TList *listNsigma = new TList;
+ listNsigma->SetOwner();
+
+ TH2 *hptr = NULL;
+ for(Int_t ispec = 0; ispec < 4; ispec++){
+ for(Int_t istep = 0; istep < 2; istep++){
+ hptr = MakeSpectrumdEdx(static_cast<AliHFEdetPIDqa::EStep_t>(istep), specind[ispec]);
+ hptr->SetName(Form("hTPCdEdx%s%s", specnames[ispec].Data(), istep == 0 ? "Before" : "After"));
+ listdEdx->Add(hptr);
+ hptr = MakeSpectrumNSigma(static_cast<AliHFEdetPIDqa::EStep_t>(istep), specind[ispec]);
+ hptr->SetName(Form("hTPCnsigma%s%s", specnames[ispec].Data(), istep == 0 ? "Before" : "After"));
+ listNsigma->Add(hptr);
+ }
+ }
+
+ b->Add(listdEdx, "Projections dE/dx");
+ b->Add(listNsigma, "Projections NSigma");
+ }
+ }
+}
+
//_________________________________________________________
void AliHFEtpcPIDqa::Initialize(){
//
// Make common binning
const Int_t kNdim = 5;
const Int_t kPIDbins = AliPID::kSPECIES + 1;
- const Int_t kPbins = 100;
const Int_t kSteps = 2;
const Int_t kCentralityBins = 11;
const Double_t kMinPID = -1;
const Double_t kMinP = 0.;
const Double_t kMaxPID = (Double_t)AliPID::kSPECIES;
const Double_t kMaxP = 20.;
-
+ // Quantities where one can switch between low and high resolution
+ Int_t kPbins = fQAmanager->HasHighResolutionHistos() ? 1000 : 100;
+ Int_t kDedxbins = fQAmanager->HasHighResolutionHistos() ? 400 : 200;
+ Int_t kSigmaBins = fQAmanager->HasHighResolutionHistos() ? 1400 : 240;
+
// 1st histogram: TPC dEdx: (species, p, dEdx, step)
- const Int_t kDedxbins = 200;
Int_t nBinsdEdx[kNdim] = {kPIDbins, kPbins, kDedxbins, kSteps, kCentralityBins};
Double_t mindEdx[kNdim] = {kMinPID, kMinP, 0., 0., 0.};
Double_t maxdEdx[kNdim] = {kMaxPID, kMaxP, 200, 2., 11.};
fHistos->CreateTHnSparse("tpcDedx", "TPC signal; species; p [GeV/c]; TPC signal [a.u.]; Centrality; Selection Step", kNdim, nBinsdEdx, mindEdx, maxdEdx);
// 2nd histogram: TPC sigmas: (species, p nsigma, step)
- const Int_t kSigmaBins = 240;
Int_t nBinsSigma[kNdim] = {kPIDbins, kPbins, kSigmaBins, kSteps, kCentralityBins};
Double_t minSigma[kNdim] = {kMinPID, kMinP, -12., 0., 0.};
Double_t maxSigma[kNdim] = {kMaxPID, kMaxP, 12., 2., 100.};
AliHFEpidTPC *tpcpid = dynamic_cast<AliHFEpidTPC *>(fQAmanager->GetDetectorPID(AliHFEpid::kTPCpid));
Double_t contentSignal[5];
contentSignal[0] = species;
- contentSignal[1] = tpcpid->GetP(track->GetRecTrack(), anatype);
+ contentSignal[1] = tpcpid ? tpcpid->GetP(track->GetRecTrack(), anatype) : 0.;
contentSignal[2] = GetTPCsignal(track->GetRecTrack(), anatype);
contentSignal[3] = step;
contentSignal[4] = centrality;
- (dynamic_cast<THnSparseF *>(fHistos->Get("tpcDedx")))->Fill(contentSignal);
+ fHistos->Fill("tpcDedx", contentSignal);
- contentSignal[2] = tpcpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype);
- (dynamic_cast<THnSparseF *>(fHistos->Get("tpcnSigma")))->Fill(contentSignal);
+ contentSignal[2] = tpcpid ? tpcpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype) : 0.;
+ fHistos->Fill("tpcnSigma", contentSignal);
}
//_________________________________________________________
// Plot the Spectrum
//
THnSparseF *hSignal = dynamic_cast<THnSparseF *>(fHistos->Get("tpcDedx"));
+ if(!hSignal) return NULL;
hSignal->GetAxis(3)->SetRange(istep + 1, istep + 1);
- if(species > 0 && species < AliPID::kSPECIES)
+ if(species >= 0 && species < AliPID::kSPECIES)
hSignal->GetAxis(0)->SetRange(2 + species, 2 + species);
TH2 *hTmp = hSignal->Projection(2,1);
- Char_t hname[256], htitle[256];
- sprintf(hname, "hTPCsignal%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
- sprintf(htitle, "TPC dE/dx Spectrum %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
+ TString hname = Form("hTPCsignal%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after"),
+ htitle = Form("TPC dE/dx Spectrum %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
if(species > -1){
- strncat(hname, AliPID::ParticleName(species), strlen(AliPID::ParticleName(species)));
- Char_t speciesname[256];
- sprintf(speciesname, " for %ss", AliPID::ParticleName(species));
- strncat(htitle, speciesname, strlen(speciesname));
+ hname += AliPID::ParticleName(species);
+ htitle += Form(" for %ss", AliPID::ParticleName(species));
}
- hTmp->SetName(hname);
- hTmp->SetTitle(htitle);
+ hTmp->SetName(hname.Data());
+ hTmp->SetTitle(htitle.Data());
hTmp->SetStats(kFALSE);
hTmp->GetXaxis()->SetTitle("p [GeV/c]");
hTmp->GetYaxis()->SetTitle("TPC signal [a.u.]");
// Plot the Spectrum
//
THnSparseF *hSignal = dynamic_cast<THnSparseF *>(fHistos->Get("tpcnSigma"));
+ if(!hSignal) return NULL;
hSignal->GetAxis(3)->SetRange(istep + 1, istep + 1);
if(species >= 0 && species < AliPID::kSPECIES)
hSignal->GetAxis(0)->SetRange(2 + species, 2 + species);
TH2 *hTmp = hSignal->Projection(2,1);
- Char_t hname[256], htitle[256];
- sprintf(hname, "hTPCsigma%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
- sprintf(htitle, "TPC dE/dx Spectrum[#sigma] %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
+ TString hname = Form("hTPCsigma%s", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after"),
+ htitle = Form("TPC dE/dx Spectrum[#sigma] %s selection", istep == AliHFEdetPIDqa::kBeforePID ? "before" : "after");
if(species > -1){
- strncat(hname, AliPID::ParticleName(species), strlen(AliPID::ParticleName(species)));
- Char_t speciesname[256];
- sprintf(speciesname, " for %ss", AliPID::ParticleName(species));
- strncat(htitle, speciesname, strlen(speciesname));
+ hname += AliPID::ParticleName(species);
+ htitle += Form(" for %ss", AliPID::ParticleName(species));
}
- hTmp->SetName(hname);
- hTmp->SetTitle(htitle);
+ hTmp->SetName(hname.Data());
+ hTmp->SetTitle(htitle.Data());
hTmp->SetStats(kFALSE);
hTmp->GetXaxis()->SetTitle("p [GeV/c]");
hTmp->GetYaxis()->SetTitle("TPC dE/dx - <dE/dx>|_{el} [#sigma]");
#include "AliHFEpidBase.h"
#endif
+class TBrowser;
class TH2;
class AliHFEcollection;
class AliVParticle;
~AliHFEtpcPIDqa();
void Copy(TObject &o) const;
virtual Long64_t Merge(TCollection *col);
+ virtual void Browse(TBrowser *b);
+ virtual Bool_t IsFolder() const { return kTRUE; };
virtual void Initialize();
virtual void ProcessTrack(const AliHFEpidObject *track, AliHFEdetPIDqa::EStep_t step);
AliHFEcutStep *cutStep = NULL;
for(Int_t istep = 0; istep < fCutSteps->GetEntriesFast(); istep++){
cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(istep));
- ctmp->SetStepTitle(istep + 1, cutStep->GetName());
+ if(cutStep) ctmp->SetStepTitle(istep + 1, cutStep->GetName());
}
fEfficiencyContainers->AddAt(ctmp, 0);
}
AliHFEcutStep *cutStep = NULL;
for(Int_t istep = 0; istep < fCutSteps->GetEntriesFast(); istep++){
cutStep = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(istep));
- ctmp->SetStepTitle(istep + 1, cutStep->GetName());
+ if(cutStep) ctmp->SetStepTitle(istep + 1, cutStep->GetName());
}
}
ReleaseContainers();
//
fMC = mc;
fMCsignal->SetMC(fMC);
+ AliHFEcutStep *cs = NULL;
for(Int_t icut = 0; icut < fCutSteps->GetEntriesFast(); icut++)
- (dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))->SetMC(fMC);
+ if((cs = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))) cs->SetMC(fMC);
}
//__________________________________________________________________
//
// Publish MC event to the single cut steps
//
+ AliHFEcutStep *cs = NULL;
for(Int_t icut = 0; icut < fCutSteps->GetEntriesFast(); icut++)
- (dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))->SetRecEvent(rec);
+ if((cs = dynamic_cast<AliHFEcutStep *>(fCutSteps->UncheckedAt(icut)))) cs->SetRecEvent(rec);
}
// Markus Fasel <M.Fasel@gsi.de>
//
#include <TAxis.h>
+#include <TBrowser.h>
#include <TClass.h>
#include <TCanvas.h>
#include <TF1.h>
AliHFEtrdPIDqa::AliHFEtrdPIDqa():
TNamed("trdPIDqa", ""),
fTRDpid(NULL),
- fLikeTRD(NULL),
- fQAtrack(NULL),
- fQAdEdx(NULL),
- fTRDtruncMean(NULL),
+ fHistos(NULL),
fPionEfficiencies(NULL),
fProtonEfficiencies(NULL),
fKaonEfficiencies(NULL),
AliHFEtrdPIDqa::AliHFEtrdPIDqa(const Char_t *name):
TNamed(name, ""),
fTRDpid(NULL),
- fLikeTRD(NULL),
- fQAtrack(NULL),
- fQAdEdx(NULL),
- fTRDtruncMean(NULL),
+ fHistos(NULL),
fPionEfficiencies(NULL),
fProtonEfficiencies(NULL),
fKaonEfficiencies(NULL),
AliHFEtrdPIDqa::AliHFEtrdPIDqa(const AliHFEtrdPIDqa &ref):
TNamed(ref),
fTRDpid(NULL),
- fLikeTRD(NULL),
- fQAtrack(NULL),
- fQAdEdx(NULL),
- fTRDtruncMean(NULL),
+ fHistos(NULL),
fPionEfficiencies(NULL),
fProtonEfficiencies(NULL),
fKaonEfficiencies(NULL),
// Destructor
//
if(fTRDpid) delete fTRDpid;
- if(fLikeTRD) delete fLikeTRD;
- if(fQAtrack) delete fQAtrack;
- if(fQAdEdx) delete fQAdEdx;
- if(fTRDtruncMean) delete fTRDtruncMean;
-
+ if(fHistos) delete fHistos;
if(fPionEfficiencies) delete fPionEfficiencies;
if(fProtonEfficiencies) delete fProtonEfficiencies;
if(fKaonEfficiencies) delete fKaonEfficiencies;
AliHFEtrdPIDqa &target = dynamic_cast<AliHFEtrdPIDqa &>(ref);
target.fTRDpid = fTRDpid;
- target.fLikeTRD = dynamic_cast<THnSparseF *>(fLikeTRD->Clone());
- target.fQAtrack = dynamic_cast<THnSparseF *>(fQAtrack->Clone());
- target.fQAdEdx = dynamic_cast<THnSparseF *>(fQAdEdx->Clone());
- target.fTRDtruncMean = dynamic_cast<THnSparseF *>(fTRDtruncMean->Clone());
+ target.fHistos = dynamic_cast<AliHFEcollection *>(fHistos->Clone());
}
//__________________________________________________________________
if(!coll) return 0;
if(coll->IsEmpty()) return 1;
- TIterator *it = coll->MakeIterator();
+ AliHFEtrdPIDqa *refQA = NULL;
+ TIter it(coll);
TObject *o = NULL;
Long64_t count = 0;
- while((o = it->Next())){
- AliHFEtrdPIDqa *refQA = dynamic_cast<AliHFEtrdPIDqa *>(o);
+ TList listHistos;
+ while((o = it())){
+ refQA = dynamic_cast<AliHFEtrdPIDqa *>(o);
if(!refQA) continue;
- if(fLikeTRD) fLikeTRD->Add(refQA->fLikeTRD);
- if(fQAtrack) fQAtrack->Add(refQA->fQAtrack);
- if(fQAdEdx) fQAdEdx->Add(refQA->fQAdEdx);
- if(fTRDtruncMean) fTRDtruncMean->Add(refQA->fTRDtruncMean);
- count++;
+ listHistos.Add(refQA->fHistos);
+ count++;
}
- delete it;
+ fHistos->Merge(&listHistos);
return count+1;
}
+//__________________________________________________________________
+void AliHFEtrdPIDqa::Browse(TBrowser *b){
+ //
+ // Enable Browser functionality
+ //
+ if(b){
+ // Add objects to the browser
+ if(fHistos) b->Add(fHistos, fHistos->GetName());
+ if(fPionEfficiencies) b->Add(fPionEfficiencies, "Pion Efficiencies");
+ if(fProtonEfficiencies) b->Add(fProtonEfficiencies, "Proton Efficiencies");
+ if(fKaonEfficiencies) b->Add(fKaonEfficiencies, "Kaon Efficiencies");
+ if(fThresholds) b->Add(fThresholds, "Thresholds");
+ }
+}
+
//__________________________________________________________________
void AliHFEtrdPIDqa::Init(){
//
// Initialize Object
//
+
+ fHistos = new AliHFEcollection("TRDqa", "Histos for TRD QA");
CreateLikelihoodHistogram();
CreateQAHistogram();
binMin[kElectronLike] = 0.;
binMax[kElectronLike] = 1.;
- fLikeTRD = new THnSparseF("fLikeTRD","TRD Likelihood Studies", kQuantitiesLike, nbins, binMin, binMax);
- Double_t *binLog = AliHFEtools::MakeLogarithmicBinning(nbins[kP], binMin[kP], binMax[kP]);
- fLikeTRD->GetAxis(kP)->Set(nbins[kP], binLog);
- delete[] binLog;
+ fHistos->CreateTHnSparse("fLikeTRD","TRD Likelihood Studies", kQuantitiesLike, nbins, binMin, binMax);
+ fHistos->BinLogAxis("fLikeTRD", kP);
}
//__________________________________________________________________
binMax[kNonZeroTrackletCharge] = AliESDtrack::kTRDnPlanes + 1.;
binMax[kNClusters] = 200.;
- fQAtrack = new THnSparseF("fQAtrack","TRD QA Histogram", kQuantitiesQA, nbins, binMin, binMax);
- Double_t *binLog = AliHFEtools::MakeLogarithmicBinning(nbins[kP], binMin[kP], binMax[kP]);
- fQAtrack->GetAxis(kP)->Set(nbins[kP], binLog);
- delete[] binLog;
+ fHistos->CreateTHnSparse("fQAtrack","TRD QA Histogram", kQuantitiesQA, nbins, binMin, binMax);
+ fHistos->BinLogAxis("fQAtrack", kP);
}
//__________________________________________________________________
binMax[kNclusters] = 260.;
binMax[kNonZeroSlices] = 8.;
- fQAdEdx = new THnSparseF("fQAdEdx","TRD summed dEdx", kQuantitiesdEdx, nbins, binMin, binMax);
- Double_t *binLog = AliHFEtools::MakeLogarithmicBinning(nbins[kP], binMin[kP], binMax[kP]);
- fQAdEdx->GetAxis(kP)->Set(nbins[kP], binLog);
- delete[] binLog;
+ fHistos->CreateTHnSparse("fQAdEdx","TRD summed dEdx", kQuantitiesdEdx, nbins, binMin, binMax);
+ fHistos->BinLogAxis("fQAdEdx", kP);
+ fHistos->Sumw2("fQAdEdx");
}
//__________________________________________________________________
binMax[kTRDdEdxMethod1] = 20000.;
binMax[kTRDdEdxMethod2] = 20000.;
- fTRDtruncMean = new THnSparseF("fTRDtruncMean","TRD TruncatedMean studies", kQuantitiesTruncMean, nbins, binMin, binMax);
- Double_t *binLog = AliHFEtools::MakeLogarithmicBinning(nbins[kP], binMin[kP], binMax[kP]);
- fTRDtruncMean->GetAxis(kP)->Set(nbins[kP], binLog);
- delete[] binLog;
+ fHistos->CreateTHnSparse("fTRDtruncMean","TRD TruncatedMean studies", kQuantitiesTruncMean, nbins, binMin, binMax);
+ fHistos->BinLogAxis("fTRDtruncMean", kP);
+ fHistos->CreateTH2F("fTRDslicesPions","TRD dEdx per slice for Pions", 8, 0, 8, 500, 0, 2000);
+ fHistos->CreateTH2F("fTRDslicesElectrons","TRD dEdx per slice for Electrons", 8, 0, 8, 500, 0, 2000);
}
//
// Process single ESD track
//
+ if(!track) return;
if((track->GetStatus() & AliESDtrack::kTRDout) == 0) return; // require TRD track
FillTRDLikelihoods(track, species);
FillTRDQAplots(track, species);
// Process single AOD track
// AOD PID object required
//
+ if(!track) return;
AliAODPid *trackPID = track->GetDetPid();
if(!trackPID) return;
quantities[kNTracklets] = track->GetTRDntrackletsPID();
quantities[kElectronLike] = trdLike[AliPID::kElectron];
- fLikeTRD->Fill(quantities);
+ fHistos->Fill("fLikeTRD", quantities);
}
//__________________________________________________________________
quantitiesQA[kSpecies] = quantitiesdEdx[kSpecies]
= quantitiesTruncMean[kSpecies]
= species;
- quantitiesQA[kP] = quantitiesdEdx[kP]
- = quantitiesTruncMean[kP]
+ quantitiesQA[kP] = quantitiesTruncMean[kP]
= outerPars ? outerPars->P() : track->P();
quantitiesQA[kNTracklets] = quantitiesdEdx[kNTracklets]
= quantitiesTruncMean[kNTracklets]
Double_t dEdxSum = 0., qSlice = 0.;
// remove the last slice from the histogram
- Int_t ntrackletsNonZero = 0, nSlices = track->GetNumberOfTRDslices() - 1, nSlicesNonZero = 0;
+ Int_t ntrackletsNonZero = 0, nSlices = track->GetNumberOfTRDslices(), nSlicesNonZero = 0;
+ TString speciesname = "pions";
+ Bool_t selectedForSlicemon = kFALSE;
+
+ switch(species){
+ case AliPID::kElectron: speciesname = "Electrons"; selectedForSlicemon = kTRUE; break;
+ case AliPID::kPion: speciesname = "Pions"; selectedForSlicemon = kTRUE; break;
+ default: speciesname = "undefined"; selectedForSlicemon = kFALSE; break;
+ };
+ AliDebug(1, Form("species %d, speciesname %s, momentum %f, selected %s", species, speciesname.Data(), track->P(), selectedForSlicemon ? "yes" : "no"));
for(Int_t iplane = 0; iplane < AliESDtrack::kTRDnPlanes; iplane++){
+ quantitiesdEdx[kP] = track->GetTRDmomentum(iplane);
dEdxSum = 0.;
for(Int_t islice = 0; islice < nSlices; islice++){
qSlice = track->GetTRDslice(iplane, islice);
// cut out 0 slices
nSlicesNonZero++;
dEdxSum += qSlice;
+ // Reweighting of the slices for the truncated mean: select all pion tracks above
+ // 1.5 GeV and monitor the dEdx as function of slice
+ if(selectedForSlicemon && track->P() > 1.5){
+ AliDebug(2, Form("Filling Histogram fTRDslices%s", speciesname.Data()));
+ fHistos->Fill(Form("fTRDslices%s", speciesname.Data()), static_cast<Double_t>(islice), qSlice);
+ }
}
}
quantitiesdEdx[kNonZeroSlices] = nSlicesNonZero;
quantitiesdEdx[kdEdx] = dEdxSum;
if(dEdxSum) ntrackletsNonZero++;
// Fill dEdx histogram
- if(dEdxSum > 1e-1) fQAdEdx->Fill(quantitiesdEdx); // Cut out 0 entries
+ if(dEdxSum > 1e-1) fHistos->Fill("fQAdEdx", quantitiesdEdx); // Cut out 0 entries
}
quantitiesQA[kNonZeroTrackletCharge] = ntrackletsNonZero;
- fQAtrack->Fill(quantitiesQA);
+ fHistos->Fill("fQAtrack", quantitiesQA);
quantitiesTruncMean[kTPCdEdx] = track->GetTPCsignal();
- quantitiesTruncMean[kTRDdEdxMethod1] = fTRDpid->GetTRDSignalV1(track);
- quantitiesTruncMean[kTRDdEdxMethod2] = fTRDpid->GetTRDSignalV2(track);
- fTRDtruncMean->Fill(quantitiesTruncMean);
+ quantitiesTruncMean[kTRDdEdxMethod1] = fTRDpid->GetTRDSignalV1(track, 0.6);
+ quantitiesTruncMean[kTRDdEdxMethod2] = fTRDpid->GetTRDSignalV2(track, 0.6);
+ fHistos->Fill("fTRDtruncMean", quantitiesTruncMean);
}
/////////////////////////////////////////////////////////
printf(">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>\n");
threshhist = dynamic_cast<TGraph *>(lHistos->FindObject(Form("eff%d", static_cast<Int_t>(fgkElectronEff[ieff] * 100))));
+ if(!threshhist) continue;
threshparam = MakeThresholds(threshhist);
threshparam->SetName(Form("thresh_%d_%d", itracklet, static_cast<Int_t>(fgkElectronEff[ieff] * 100)));
+ lFormulas->Add(threshparam);
}
}
// elPion Efficiency, Proton Efficiency and Kaon Efficiency at discrete
// electron efficiencies
//
- Int_t binTracklets = fLikeTRD->GetAxis(kNTracklets)->FindBin(nTracklets);
- fLikeTRD->GetAxis(kNTracklets)->SetRange(binTracklets, binTracklets);
+ THnSparse *hLikeTRD = dynamic_cast<THnSparseF *>(fHistos->Get("fLikeTRD"));
+ if(!hLikeTRD){
+ AliError("Likelihood Histogram not available");
+ return;
+ }
+ Int_t binTracklets = hLikeTRD->GetAxis(kNTracklets)->FindBin(nTracklets);
+ hLikeTRD->GetAxis(kNTracklets)->SetRange(binTracklets, binTracklets);
- Int_t binElectrons = fLikeTRD->GetAxis(kSpecies)->FindBin(AliPID::kElectron);
+ Int_t binElectrons = hLikeTRD->GetAxis(kSpecies)->FindBin(AliPID::kElectron);
AliDebug(1, Form("BinElectrons %d", binElectrons));
- Int_t binPions = fLikeTRD->GetAxis(kSpecies)->FindBin(AliPID::kPion);
+ Int_t binPions = hLikeTRD->GetAxis(kSpecies)->FindBin(AliPID::kPion);
AliDebug(1, Form("BinPions %d", binPions));
- Int_t binProtons = fLikeTRD->GetAxis(kSpecies)->FindBin(AliPID::kProton);
+ Int_t binProtons = hLikeTRD->GetAxis(kSpecies)->FindBin(AliPID::kProton);
AliDebug(1, Form("BinProtons %d", binProtons));
- fLikeTRD->GetAxis(kSpecies)->SetRange(binElectrons, binElectrons);
- TH2 *likeElectron = fLikeTRD->Projection(kElectronLike, kP);
+ hLikeTRD->GetAxis(kSpecies)->SetRange(binElectrons, binElectrons);
+ TH2 *likeElectron = hLikeTRD->Projection(kElectronLike, kP);
likeElectron->SetName("likeElectron");
- fLikeTRD->GetAxis(kSpecies)->SetRange(binPions, binPions);
- TH2 *likePion = fLikeTRD->Projection(kElectronLike, kP);
+ hLikeTRD->GetAxis(kSpecies)->SetRange(binPions, binPions);
+ TH2 *likePion = hLikeTRD->Projection(kElectronLike, kP);
likePion->SetName("likePion");
- fLikeTRD->GetAxis(kSpecies)->SetRange(binProtons, binProtons);
- TH2 *likeProton = fLikeTRD->Projection(kElectronLike, kP);
+ hLikeTRD->GetAxis(kSpecies)->SetRange(binProtons, binProtons);
+ TH2 *likeProton = hLikeTRD->Projection(kElectronLike, kP);
likeProton->SetName("likeProton");
// Undo ranges
- fLikeTRD->GetAxis(kSpecies)->SetRange(0, fLikeTRD->GetAxis(kSpecies)->GetNbins());
- fLikeTRD->GetAxis(kNTracklets)->SetRange(0, fLikeTRD->GetAxis(kNTracklets)->GetNbins());
+ hLikeTRD->GetAxis(kSpecies)->SetRange(0, hLikeTRD->GetAxis(kSpecies)->GetNbins());
+ hLikeTRD->GetAxis(kNTracklets)->SetRange(0, hLikeTRD->GetAxis(kNTracklets)->GetNbins());
// Prepare List for output
TList *listPions = new TList; listPions->SetName(Form("%dTracklets", nTracklets));
}
//__________________________________________________________________
-void AliHFEtrdPIDqa::DrawTracklet(Int_t itracklet){
+void AliHFEtrdPIDqa::DrawTracklet(Int_t itracklet, Bool_t doFit){
//
// Draw efficiencies and threshold as function of p
//
TList *lprotons = dynamic_cast<TList *>(fProtonEfficiencies->FindObject(Form("%dTracklets", itracklet)));
TList *lthresholds = dynamic_cast<TList *>(fThresholds->FindObject(Form("%dTracklets", itracklet)));
+ if(!(lpions && lprotons && lthresholds)){
+ AliDebug(1, "Relevant lists not found. Did you forget to run FinishAnalysis()?");
+ return;
+ }
TGraphErrors *pi, *pr;
TGraph *tr;
TLegend *leg;
TCanvas *c1 = new TCanvas(Form("tracklet%d", itracklet), Form("Tracklet %d", itracklet), 1024, 768);
c1->Divide(3,2);
+ TF1 *threshfit = NULL;
for(Int_t ieff = 0; ieff < kNElectronEffs; ieff++){
c1->cd(ieff + 1);
leg = new TLegend(0.6, 0.7, 0.89, 0.89);
pi = dynamic_cast<TGraphErrors *>(lpions->FindObject(Form("eff%d", static_cast<Int_t>(fgkElectronEff[ieff] * 100))));
pr = dynamic_cast<TGraphErrors *>(lprotons->FindObject(Form("eff%d", static_cast<Int_t>(fgkElectronEff[ieff] * 100))));
tr = dynamic_cast<TGraph *>(lthresholds->FindObject(Form("eff%d", static_cast<Int_t>(fgkElectronEff[ieff] * 100))));
+ if(!(pi && pr && tr)) continue;
// Define nice plot, draw
// Axis Title
// Draw
pi->Draw("ape"); pr->Draw("pesame"); tr->Draw("psame");
+ // Optionally do Fit
+ if(doFit){
+ threshfit = MakeThresholds(tr);
+ threshfit->SetLineColor(kBlack);
+ threshfit->Draw("same");
+ }
+
// Add entries to legend
leg->AddEntry(pi, "Pion Efficiency", "lp");
leg->AddEntry(pr, "Proton Efficiency", "lp");
//
TF1 *threshparam = new TF1("thresh", "1-[0]-[1]*x-[2]*TMath::Exp(-[3]*x)", 0.1, 10);
- threshist->Fit(threshparam, "NE", "", 0, 10);
+ threshist->Fit(threshparam, "NE", "", 0.1, 3.5);
return threshparam;
}
#include <THnSparse.h>
#endif
+#ifndef ROOT_TH2
+#include <TH2.h>
+#endif
+
+#ifndef ALIHFECOLLECTION_H
+#include "AliHFEcollection.h"
+#endif
+
+class TBrowser;
class TCollection;
class TF1;
class TGraph;
class TH1;
+//class TH2;
class TList;
class TObjArray;
virtual void Copy(TObject &o) const;
virtual Long64_t Merge(TCollection *coll);
virtual ~AliHFEtrdPIDqa();
+
+ virtual void Browse(TBrowser *b);
+ virtual Bool_t IsFolder() const { return kTRUE; }
void ProcessTracks(TObjArray * const l, Int_t species);
void ProcessTrack(AliVTrack *track, Int_t species);
void StoreResults(const Char_t *filename = "HFEtrdPIDqa.root");
void SaveThresholdParameters(const Char_t * filename = "TRD.Thresholds.root");
- void DrawTracklet(Int_t tracklet);
+ void DrawTracklet(Int_t tracklet, Bool_t doFit = kFALSE);
void ClearLists();
//---------------------------------------------------
// Getters for Histograms
- THnSparseF *GetLikelihoodHistogram() const { return fLikeTRD; }
- THnSparseF *GetQAHistogram() const { return fQAtrack; }
- THnSparseF *GetdEdxHistogram() const { return fQAdEdx; }
- THnSparseF *GetHistoTruncMean() const { return fTRDtruncMean; }
+ THnSparseF *GetLikelihoodHistogram() const { return fHistos ? dynamic_cast<THnSparseF *>(fHistos->Get("fLikeTRD")) : NULL; }
+ THnSparseF *GetQAHistogram() const { return fHistos ? dynamic_cast<THnSparseF *>(fHistos->Get("fQAtrack")) : NULL; }
+ THnSparseF *GetdEdxHistogram() const { return fHistos ? dynamic_cast<THnSparseF *>(fHistos->Get("fQAdEdx")) : NULL; }
+ THnSparseF *GetHistoTruncMean() const { return fHistos ? dynamic_cast<THnSparseF *>(fHistos->Get("fTRDtruncMean")) : NULL; }
+ TH2 *GetSliceChargePions(Bool_t pions = kTRUE) const {
+ const Char_t * species = pions ? "Pions" : "Electrons";
+ return fHistos ? dynamic_cast<TH2 *>(fHistos->Get(Form("fTRDslices%s", species))) : NULL;
+ }
+ AliHFEcollection *GetHistos() const { return fHistos; }
//---------------------------------------------------
protected:
// Description of the containers we use to store basic information
static const Double_t fgkMinBinCommon[kQuantitiesCommon]; // Bin Limits for common quantities (lower limit)
static const Double_t fgkMaxBinCommon[kQuantitiesCommon]; // Bin Limits for common quantities (upper limit)
AliHFEpidTRD *fTRDpid; // HFE PID for TRD
- THnSparseF *fLikeTRD; // Histo for Likelihoods
- THnSparseF *fQAtrack; // QA histo for quantities based on track level
- THnSparseF *fQAdEdx; // QA for tracklet charge
- THnSparseF *fTRDtruncMean; // QA for truncated mean
+ AliHFEcollection *fHistos; // Histogram collection
// List for Histograms:
TList *fPionEfficiencies; //! List for Pion efficiencies
TList *fThresholds; //! List for Threshold Graphs
- ClassDef(AliHFEtrdPIDqa, 2) // QA class for TRD PID
+ ClassDef(AliHFEtrdPIDqa, 3) // QA class for TRD PID
};
#endif
Double_t minTPCSigma[4] = {kMinPID, kMinP, -12., 0};
Double_t maxTPCSigma[4] = {kMaxPID, kMaxP, 12., 2.};
fHistos->CreateTHnSparse("hTPCsigma", "TPC sigma; species p [GeV/c]; TPC dEdx - <dE/dx>|_{el} [#sigma]; selection step", 4, nBinsTPCSigma, minTPCSigma, maxTPCSigma);
+ fHistos->Sumw2("hTPCsigma");
// Create Monitoring histogram for the Likelihood distribution
const Int_t kTRDLikelihoodBins = 100;
Int_t nBinsTRDlike[4] = {kPIDbins, kPbins, kTRDLikelihoodBins, kSteps};
Double_t minTRDlike[4] = {kMinPID, kMinP, 0., 0.};
Double_t maxTRDlike[4] = {kMaxPID, kMaxP, 1., 2.};
fHistos->CreateTHnSparse("hTRDlikelihood", "TRD Likelihood Distribution; species; p [GeV/c]; TRD electron Likelihood; selection step", 4, nBinsTRDlike, minTRDlike, maxTRDlike);
+ fHistos->Sumw2("hTRDlikelihood");
// Create Monitoring histogram for the TRD total charge
const Int_t kTRDchargeBins = 1000;
Int_t nBinsTRDcharge[4] = {kPIDbins, kPbins, kTRDchargeBins, kSteps};
Double_t minTRDcharge[4] = {kMinPID, kMinP, 0., 0.};
Double_t maxTRDcharge[4] = {kMaxPID, kMaxP, 100000., 2.};
fHistos->CreateTHnSparse("hTRDcharge", "Total TRD charge; species; p [GeV/c]; TRD charge [a.u.]; selection step", 4, nBinsTRDcharge, minTRDcharge, maxTRDcharge);
+ fHistos->Sumw2("hTRDcharge");
}
//____________________________________________________________
Double_t container[4];
container[0] = species;
- container[1] = trdpid->GetP(track->GetRecTrack(), anatype);
- container[2] = tpcpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype);
+ container[1] = trdpid ? trdpid->GetP(track->GetRecTrack(), anatype) : 0.;
+ container[2] = tpcpid ? tpcpid->NumberOfSigmas(track->GetRecTrack(), AliPID::kElectron, anatype) : 0.;
container[3] = step;
fHistos->Fill("hTPCsigma", container);
git://git.uio.no / u / mrichter / AliRoot.git / commitdiff