#include "AliHFEV0info.h"
#include "AliHFEV0pid.h"
#include "AliHFEV0pidMC.h"
+#include "AliHFEV0cuts.h"
#include "AliHFEtrdPIDqa.h"
// Prepare task output
//
+ // load networks
if(fNNref){
- for(Int_t mom = 0; mom < 11; mom++){ // load networks
+ for(Int_t mom = 0; mom < 11; mom++){
fNet[mom] = (TMultiLayerPerceptron*) fNNref->Get(Form("NN_Mom%d", mom));
if(!fNet[mom]){
AliError(Form("No reference network for momentum bin %d!", mom));
}
}
- fV0pid = new AliHFEV0pid;
+ fV0pid = new AliHFEV0pid();
if(HasV0pidQA()) fV0pid->InitQA();
fV0pidMC = new AliHFEV0pidMC();
fV0pidMC->Init();
//
// event based histograms
//
- Int_t n_T0[2] = {10000, 100};
- Double_t min_T0[2] = {114500, -1000};
- Double_t max_T0[2] = {124500, 1000};
- fOutput->CreateTHnSparse("hEvent_T0", "T0 as a function of run number; run number; T0 (ns)", 2, n_T0, min_T0, max_T0);
+ Int_t nT0[2] = {10000, 100};
+ Double_t minT0[2] = {114500, -1000};
+ Double_t maxT0[2] = {124500, 1000};
+ fOutput->CreateTHnSparse("hEvent_T0", "T0 as a function of run number; run number; T0 (ns)", 2, nT0, minT0, maxT0);
//
// test the tender V0 supply
}
// TRD THnSparse - entries per tracklet
// species, p, tracklet position, number of PID tracklets, number of slices (non 0), number of clusters, electron likelihood,
+
{
Int_t nbins[7] = {5, 20, 6, 7, 9, 45, 100};
Double_t min[7] = {-0.5, 0.5, -0.5, -0.5, -0.5, -0.5, 0.};
AliError("AliVEvent not available, returning");
}
+
if(fMC) fV0pidMC->SetMCEvent(fMC);
+ if(fMC) fV0pid->SetMCEvent(fMC);
fV0pid->Process(fEvent);
TObjArray *hfeelectrons = fV0pid->GetListOfElectrons();
TObjArray *cleanElectrons = MakeCleanListElectrons(hfeelectrons);
if(fMC){
- fV0pidMC->Process(electrons, AliHFEV0pid::kRecoElectron);
- fV0pidMC->Process(pionsK0, AliHFEV0pid::kRecoPionK0);
- fV0pidMC->Process(pionsL, AliHFEV0pid::kRecoPionL);
- fV0pidMC->Process(protons, AliHFEV0pid::kRecoProton);
+ fV0pidMC->Process(electrons, AliHFEV0cuts::kRecoElectron);
+ fV0pidMC->Process(pionsK0, AliHFEV0cuts::kRecoPionK0);
+ fV0pidMC->Process(pionsL, AliHFEV0cuts::kRecoPionL);
+ fV0pidMC->Process(protons, AliHFEV0cuts::kRecoProton);
}
AliDebug(2, Form("Number of Electrons : %d", electrons->GetEntries()));
if(fMC){
AliDebug(2, "MC Information available. Doing Purity checks...");
// Calculate the purity of the clean samples using MC
- MakePurity(electrons, AliHFEV0pid::kRecoElectron);
- MakePurity(pionsK0, AliHFEV0pid::kRecoPionK0);
- MakePurity(pionsL, AliHFEV0pid::kRecoPionL);
- MakePurity(protons, AliHFEV0pid::kRecoProton);
+ MakePurity(electrons, AliHFEV0cuts::kRecoElectron);
+ MakePurity(pionsK0, AliHFEV0cuts::kRecoPionK0);
+ MakePurity(pionsL, AliHFEV0cuts::kRecoPionL);
+ MakePurity(protons, AliHFEV0cuts::kRecoProton);
}
// some event wise checks
CheckEvent();
// Make Illumination Plot
- FillIllumination(electrons, AliHFEV0pid::kRecoElectron);
- FillIllumination(pionsK0, AliHFEV0pid::kRecoPionK0);
- FillIllumination(pionsL, AliHFEV0pid::kRecoPionL);
- FillIllumination(protons, AliHFEV0pid::kRecoProton);
+ FillIllumination(electrons, AliHFEV0cuts::kRecoElectron);
+ FillIllumination(pionsK0, AliHFEV0cuts::kRecoPionK0);
+ FillIllumination(pionsL, AliHFEV0cuts::kRecoPionL);
+ FillIllumination(protons, AliHFEV0cuts::kRecoProton);
// Now we can do studies on the PID itself
// For TRD use the TRD PID QA object
fTRDpidQA->ProcessTracks(pionsL, AliPID::kPion);
fTRDpidQA->ProcessTracks(protons, AliPID::kProton);
- FillElectronLikelihoods(electrons, AliHFEV0pid::kRecoElectron);
- FillElectronLikelihoods(pionsK0, AliHFEV0pid::kRecoPionK0);
- FillElectronLikelihoods(pionsL, AliHFEV0pid::kRecoPionL);
- FillElectronLikelihoods(protons, AliHFEV0pid::kRecoProton);
+ FillElectronLikelihoods(electrons, AliHFEV0cuts::kRecoElectron);
+ FillElectronLikelihoods(pionsK0, AliHFEV0cuts::kRecoPionK0);
+ FillElectronLikelihoods(pionsL, AliHFEV0cuts::kRecoPionL);
+ FillElectronLikelihoods(protons, AliHFEV0cuts::kRecoProton);
- FillPIDresponse(electrons, AliHFEV0pid::kRecoElectron);
- FillPIDresponse(pionsK0, AliHFEV0pid::kRecoPionK0);
- FillPIDresponse(pionsL, AliHFEV0pid::kRecoPionL);
- FillPIDresponse(protons, AliHFEV0pid::kRecoProton);
+ FillPIDresponse(electrons, AliHFEV0cuts::kRecoElectron);
+ FillPIDresponse(pionsK0, AliHFEV0cuts::kRecoPionK0);
+ FillPIDresponse(pionsL, AliHFEV0cuts::kRecoPionL);
+ FillPIDresponse(protons, AliHFEV0cuts::kRecoProton);
// check the tender V0s
- CheckTenderV0pid(electrons, AliHFEV0pid::kRecoElectron);
- CheckTenderV0pid(pionsK0, AliHFEV0pid::kRecoPionK0);
- CheckTenderV0pid(pionsL, AliHFEV0pid::kRecoPionL);
- CheckTenderV0pid(protons, AliHFEV0pid::kRecoProton);
+ CheckTenderV0pid(electrons, AliHFEV0cuts::kRecoElectron);
+ CheckTenderV0pid(pionsK0, AliHFEV0cuts::kRecoPionK0);
+ CheckTenderV0pid(pionsL, AliHFEV0cuts::kRecoPionL);
+ CheckTenderV0pid(protons, AliHFEV0cuts::kRecoProton);
// Analysis done, flush the containers
fV0pid->Flush();
}
//__________________________________________
-void AliHFEpidQA::FillIllumination(TObjArray * const tracks, Int_t species){
+void AliHFEpidQA::FillIllumination(const TObjArray * const tracks, Int_t species){
//
// Fill Illumination Plot
//
while((o = trackIter())){
if(!TString(o->IsA()->GetName()).CompareTo("AliESDtrack")){
// work on local copy in order to not spoil others
- esdtrack = new AliESDtrack(*(dynamic_cast<AliESDtrack *>(o)));
+ esdtrack = new AliESDtrack(*(static_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));
+ esdtrack = new AliESDtrack(static_cast<AliAODTrack *>(o));
+ if(!esdtrack) continue;
} else {
// Non usable
continue;
}
//__________________________________________
-void AliHFEpidQA::MakePurity(TObjArray *tracks, Int_t species){
+void AliHFEpidQA::MakePurity(const TObjArray *tracks, Int_t species){
//
// Fill the QA histos for a given species
//
if(!fMC) return;
AliDebug(3, Form("Doing Purity checks for species %d", species));
Int_t pdg = GetPDG(species);
- Char_t hname[256];
+ TString hname;
switch(species){
- case AliHFEV0pid::kRecoElectron:
- sprintf(hname, "purityElectron");
+ case AliHFEV0cuts::kRecoElectron:
+ hname = "purityElectron";
break;
- case AliHFEV0pid::kRecoPionK0:
- sprintf(hname, "purityPionK0");
+ case AliHFEV0cuts::kRecoPionK0:
+ hname = "purityPionK0";
break;
- case AliHFEV0pid::kRecoPionL:
- sprintf(hname, "purityPionL");
+ case AliHFEV0cuts::kRecoPionL:
+ hname = "purityPionL";
break;
- case AliHFEV0pid::kRecoProton:
- sprintf(hname, "purityProton");
+ case AliHFEV0cuts::kRecoProton:
+ hname = "purityProton";
break;
default: // non investigated species
AliDebug(3, "Species not investigated");
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
}
//__________________________________________
-void AliHFEpidQA::FillElectronLikelihoods(TObjArray * const particles, Int_t species){
+void AliHFEpidQA::FillElectronLikelihoods(const TObjArray * const particles, Int_t species){
//
// Fill electron Likelihoods for the ITS, TPC and TOF
// Required for the calculation of the electron efficiency,
// pion and proton efficiency and the thresholds
//
Long_t status = 0;
- const Char_t *detname[4] = {"ITS", "TPC", "TRD", "TOF"};
- Char_t specname[256];
+ const TString detname[4] = {"ITS", "TPC", "TRD", "TOF"};
+ TString specname;
switch(species){
- case AliHFEV0pid::kRecoElectron:
- sprintf(specname, "Electron");
+ case AliHFEV0cuts::kRecoElectron:
+ specname = "Electron";
break;
- case AliHFEV0pid::kRecoPionK0:
- sprintf(specname, "PionK0");
+ case AliHFEV0cuts::kRecoPionK0:
+ specname = "PionK0";
break;
- case AliHFEV0pid::kRecoPionL:
- sprintf(specname, "PionL");
+ case AliHFEV0cuts::kRecoPionL:
+ specname = "PionL";
break;
- case AliHFEV0pid::kRecoProton:
- sprintf(specname, "Proton");
+ case AliHFEV0cuts::kRecoProton:
+ specname = "Proton";
break;
default:
AliDebug(2, Form("Species %d not investigated", species));
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
quantities[0] = pTPC;
Bool_t detFlagSet = kFALSE;
for(Int_t idet = 0; idet < 4; idet++){
- Char_t histname[256], histnameMC[256];
- sprintf(histname, "h%s_El_like_%s", detname[idet], specname);
- sprintf(histnameMC, "h%s_El_like_MC_%s", detname[idet], specname);
+ TString histname, histnameMC;
+ histname = "h" + detname[idet] + "_El_like_" + specname;
+ histnameMC = "h" + detname[idet] + "_El_like_MC_" + specname;
+
switch(idet){
case kITS: esdTrack->GetITSpid(pidProbs);
detFlagSet = status & AliESDtrack::kITSpid;
}//.. while tracks
}
//__________________________________________
-void AliHFEpidQA::FillPIDresponse(TObjArray * const particles, Int_t species){
+void AliHFEpidQA::FillPIDresponse(const TObjArray * const particles, Int_t species){
//
// Fill the PID response of different detectors to V0 daughter particles
//
- Char_t hname[256] = "";
- Char_t hname2[256] = "";
- Char_t hname3[256] = "";
+ TString hname, hname2, hname3;
- const Char_t *typeName[5] = {"Electron", "PionK0", "PionL", "Kaon", "Proton"};
+ const TString typeName[5] = {"Electron", "PionK0", "PionL", "Kaon", "Proton"};
const Int_t typePID[5] = {0, 2, 2, 3, 4};
// PID THnSparse
// 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
esdTrack->GetITSdEdxSamples(dEdxSamples);
Int_t nSamples = 0;
Double_t dEdxSum = 0.;
- sprintf(hname, "hITS_dEdx_%s", typeName[species]);
+ hname = "hITS_dEdx_" + typeName[species];
for(Int_t i=0; i<4; ++i){
if(dEdxSamples[i] > 0){
nSamples++;
fOutput->Fill("hITS_dEdx_nSamples", nSamples);
Double_t signal = esdTrack->GetITSsignal();
- sprintf(hname, "hITS_Signal_%s", typeName[species]);
+ hname = "hITS_Signal_" + typeName[species];
fOutput->Fill(hname, p, signal);
data[4] = signal;
// ITS number of signas
Double_t nsigma = fESDpid->NumberOfSigmasITS(esdTrack,(AliPID::EParticleType)typePID[species]);
- sprintf(hname, "hITS_nSigma_%s", typeName[species]);
+ hname = "hITS_nSigma_" + typeName[species];
fOutput->Fill(hname, p, nsigma);
-
// ITS PID response
Double_t itsPID[5] = {-1, -1, -1, -1, -1};
esdTrack->GetITSpid(itsPID);
Int_t ix = GetMaxPID(itsPID);
- sprintf(hname, "hITS_PID_p_%s", typeName[species]);
+ hname = "hITS_PID_p_" + typeName[species];
fOutput->Fill(hname, p, ix);
}//.. kITSpid
Double_t p = esdTrack->GetInnerParam() ? esdTrack->GetInnerParam()->P() : esdTrack->P();
// TPC dEdx
Double_t dEdx = esdTrack->GetTPCsignal();
- sprintf(hname, "hTPC_dEdx_%s", typeName[species]);
+ hname = "hTPC_dEdx_" + typeName[species];
fOutput->Fill(hname, p, dEdx);
data[6] = dEdx;
//TPC number of sigmas
Double_t nsigma = fESDpid->NumberOfSigmasTPC(esdTrack,(AliPID::EParticleType)typePID[species]);
- sprintf(hname, "hTPC_nSigma_%s", typeName[species]);
+ hname = "hTPC_nSigma_" + typeName[species];
fOutput->Fill(hname, p, nsigma);
data[7] = nsigma;
// TPC PID response
- sprintf(hname, "hTPC_PID_p_%s", typeName[species]);
+ hname = "hTPC_PID_p_" + typeName[species];
Double_t tpcPID[5] = {-1, -1, -1, -1, -1};
esdTrack->GetTPCpid(tpcPID);
Int_t ix = GetMaxPID(tpcPID);
// TRD number of tracklets
Int_t ntrk = esdTrack->GetTRDntrackletsPID();
- sprintf(hname, "hTRD_trk_%s", typeName[species]);
+ hname = "hTRD_trk_" + typeName[species];
fOutput->Fill(hname, p, ntrk);
data[8] = ntrk;
// TRD PID response
- sprintf(hname, "hTRD_PID_p_%s", typeName[species]);
+ hname = "hTRD_PID_p_" + typeName[species];
Double_t trdPID[5] = {-1., -1., -1., -1., -1};
esdTrack->GetTRDpid(trdPID);
Int_t ix = GetMaxPID(trdPID);
fOutput->Fill(hname, p, ix);
// TRD n clusters
Int_t ncls = esdTrack->GetTRDncls();
- sprintf(hname, "hTRD_cls_%s", typeName[species]);
+ hname = "hTRD_cls_" + typeName[species];
fOutput->Fill(hname, p, ncls);
data[9] = ncls;
// TRD - per tracklet - dEdx per, likelihood
- sprintf(hname, "hTRD_Nslc_%s", typeName[species]);
- sprintf(hname2, "hTRD_slc_%s", typeName[species]);
- sprintf(hname3, "hTRD_dEdx_%s", typeName[species]);
+ hname = "hTRD_Nslc_" + typeName[species];
+ hname2 = "hTRD_slc_" + typeName[species];
+ hname3 = "hTRD_dEdx_" + typeName[species];
Int_t nSlices = esdTrack->GetNumberOfTRDslices();
Double_t sumTot = 0.;
Int_t not0Tot = 0;
//
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]);
+ hname = "hTOF_beta_" + typeName[species];
Float_t beta = TOFbeta(esdTrack);
fOutput->Fill(hname, p, beta);
fOutput->Fill("hTOF_beta_all", p, beta);
// TOF nSigma
- Double_t nsigma = fESDpid->NumberOfSigmasTOF(esdTrack,(AliPID::EParticleType)typePID[species], t0);
- sprintf(hname, "hTOF_nSigma_%s", typeName[species]);
+ Double_t nsigma = fESDpid->NumberOfSigmasTOF(esdTrack,(AliPID::EParticleType)typePID[species], 0.);
+ hname = "hTOF_nSigma_" + typeName[species];
fOutput->Fill(hname, p, nsigma);
if(beta > 0.97 && beta < 1.03){
data[11] = 1;
}
// TOF PID response
- sprintf(hname, "hTOF_PID_p_%s", typeName[species]);
+ hname = "hTOF_PID_p_" + typeName[species];
Double_t tofPID[5] = {-1., -1., -1., -1., -1};
esdTrack->GetTOFpid(tofPID);
Int_t ix = GetMaxPID(tofPID);
return beta;
}
//____________________________________________
-Int_t AliHFEpidQA::GetMaxPID(Double_t *pidProbs) const {
+Int_t AliHFEpidQA::GetMaxPID(const Double_t *pidProbs) const {
//
// return the index of maximal PID probability
//
Int_t pdg = 0;
switch(species){
- case AliHFEV0pid::kRecoElectron:
+ case AliHFEV0cuts::kRecoElectron:
pdg = TMath::Abs(kElectron);
break;
- case AliHFEV0pid::kRecoPionK0:
+ case AliHFEV0cuts::kRecoPionK0:
pdg = TMath::Abs(kPiPlus);
break;
- case AliHFEV0pid::kRecoPionL:
+ case AliHFEV0cuts::kRecoPionL:
pdg = TMath::Abs(kPiPlus);
break;
- case AliHFEV0pid::kRecoProton:
+ case AliHFEV0cuts::kRecoProton:
pdg = TMath::Abs(kProton);
break;
default: // non investigated species
}
//_____________________________________________
-TObjArray * AliHFEpidQA::MakeTrackList(TObjArray *tracks) const {
+TObjArray * AliHFEpidQA::MakeTrackList(const TObjArray *tracks) const {
//
// convert list of AliHFEV0Info into a list of AliVParticle
//
}
//_____________________________________________
-TObjArray * AliHFEpidQA::MakeCleanListElectrons(TObjArray *electrons) const {
+TObjArray * AliHFEpidQA::MakeCleanListElectrons(const TObjArray *electrons) const {
//
// Cleanup electron sample using TPC PID
// PID requirement will allways be implemented to the pair
TIter candidates(electrons);
AliESDEvent *esd; AliAODEvent *aod;
AliHFEV0info *hfetrack;
+ // const Double_t kSigmaTight = 1.;
+ // const Double_t kSigmaLoose = 4.;
+ const Double_t kSigmaTight = 2.;
+ const Double_t kSigmaLoose = 2.;
+ const Double_t shift = 0.57;
if((esd = dynamic_cast<AliESDEvent *>(fEvent))){
AliESDtrack *track = NULL, *partnerTrack = NULL;
while((hfetrack = dynamic_cast<AliHFEV0info *>(candidates()))){
track = dynamic_cast<AliESDtrack *>(hfetrack->GetTrack());
+ if(!track) continue;
partnerTrack = esd->GetTrack(hfetrack->GetPartnerID());
- Double_t nSigmaTrack = TMath::Abs(fESDpid->NumberOfSigmasTPC(track, AliPID::kElectron));
- Double_t nSigmaPartner = TMath::Abs(fESDpid->NumberOfSigmasTPC(partnerTrack, AliPID::kElectron));
- if((nSigmaTrack < 1 && nSigmaPartner < 4) || (nSigmaTrack < 4 && nSigmaPartner < 1))
+ Double_t nSigmaTrack = TMath::Abs(fESDpid->NumberOfSigmasTPC(track, AliPID::kElectron) - shift);
+ Double_t nSigmaPartner = TMath::Abs(fESDpid->NumberOfSigmasTPC(partnerTrack, AliPID::kElectron) - shift);
+ if((nSigmaTrack < kSigmaTight && nSigmaPartner < kSigmaLoose) || (nSigmaTrack < kSigmaLoose && nSigmaPartner < kSigmaTight))
tracks->Add(track);
}
} 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
}
}
return tracks;
}
//___________________________________________________________
-void AliHFEpidQA::CheckTenderV0pid(TObjArray * const particles, Int_t species){
+void AliHFEpidQA::CheckTenderV0pid(const TObjArray * const particles, Int_t species){
//
// retrieve the status bits from the TObject used to flag
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()
return sum;
}
//__________________________________________________________________________
-Int_t AliHFEpidQA::TRDmomBin(Double_t p){
+Int_t AliHFEpidQA::TRDmomBin(Double_t p) const {
//
// compute the momentum bin position
//
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff