PWG3/hfe/AliAnalysisTaskCheckV0tender.cxx

   1 /**************************************************************************
   2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3 *                                                                        *
   4 * Author: The ALICE Off-line Project.                                    *
   5 * Contributors are mentioned in the code where appropriate.              *
   6 *                                                                        *
   7 * Permission to use, copy, modify and distribute this software and its   *
   8 * documentation strictly for non-commercial purposes is hereby granted   *
   9 * without fee, provided that the above copyright notice appears in all   *
  10 * copies and that both the copyright notice and this permission notice   *
  11 * appear in the supporting documentation. The authors make no claims     *
  12 * about the suitability of this software for any purpose. It is          *
  13 * provided "as is" without express or implied warranty.                  *
  14 **************************************************************************/
  15 //
  16 // Task fir checking the performance of the V0 tender
  17 //
  18 //
  19 // Authors
  20 //   Matus Kalisky <matus.kalisky@cern.ch>
  21 //
  22
  23 #include <TH1F.h>
  24 #include <TList.h>
  25
  26 #include "AliAnalysisManager.h"
  27 #include "AliMCEventHandler.h"
  28 #include "AliESDInputHandler.h"
  29 #include "AliESDv0.h"
  30 #include "AliESDtrack.h"
  31 #include "AliMCParticle.h"
  32 #include "AliMCEvent.h"
  33
  34 #include "AliHFEtools.h"
  35 #include "AliHFEcollection.h"
  36
  37 #include "AliAnalysisTaskCheckV0tender.h"
  38
  39 ClassImp(AliAnalysisTaskCheckV0tender)
  40
  41 //__________________________________________________________
  42 AliAnalysisTaskCheckV0tender::AliAnalysisTaskCheckV0tender():
  43   AliAnalysisTaskSE("CheckV0tenderTask")
  44   , fOutput(0x0)
  45   , fColl(0x0)
  46   , fCollMC(0x0)
  47   , fEvents(0x0)
  48 {
  49   //
  50   // Default Constructor
  51   //
  52 }
  53 //__________________________________________________________
  54 AliAnalysisTaskCheckV0tender::AliAnalysisTaskCheckV0tender(const Char_t *name):
  55   AliAnalysisTaskSE(name)
  56   , fOutput(0x0)
  57   , fColl(0x0)
  58   , fCollMC(0x0)
  59   , fEvents(0x0)
  60 {
  61   //
  62   // Default Constructor
  63   //
  64   DefineOutput(1, TH1F::Class());
  65   DefineOutput(2, TList::Class());
  66
  67 }
  68 //__________________________________________________________
  69 AliAnalysisTaskCheckV0tender::~AliAnalysisTaskCheckV0tender(){
  70   //
  71   // Destructor
  72   //
  73
  74   if (fOutput) delete fOutput;
  75   if (fColl) delete fColl;
  76   if (fCollMC) delete fCollMC;
  77   if (fEvents) delete fEvents;
  78 }
  79 //__________________________________________________________
  80 void AliAnalysisTaskCheckV0tender::UserCreateOutputObjects(){
  81   //
  82   // prepare output objects
  83   //
  84
  85   fOutput = new TList();
  86   fOutput->SetOwner();
  87   // Counter for number of events
  88   fEvents = new TH1I("nEvents", "NumberOfEvents", 1, 1, 2);
  89
  90   fColl = new AliHFEcollection("V0_QA", "tender V0s for data");
  91   fCollMC = new AliHFEcollection("V0_MC_QA", "tender V0s for MC");
  92
  93   //
  94   // Data histos
  95   //
  96
  97   // total number of tagged V0s
  98   fColl->CreateTH1F("h_NumberOf_V0s", "Number of tagged V0s; type; counts", 4, -0.5, 3.5);
  99   // pT spectra of the tagged V0s
 100   fColl->CreateTH1F("h_Gamma_pt", "p_{T} spectrum of tagged gammas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 101   fColl->CreateTH1F("h_K0_pt", "p_{T} spectrum of tagged K0s; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 102   fColl->CreateTH1F("h_Lambda_pt", "p_{T} spectrum of tagged Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 103   fColl->CreateTH1F("h_ALambda_pt", "p_{T} spectrum of tagged A-Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 104   // invariant mass of the V0s
 105   fColl->CreateTH1F("h_Gamma_mass", "Inv. Mass of the ", 100, 0., 0.1);
 106   fColl->CreateTH1F("h_K0_mass", "Inv. Mass of the ", 100, 0.45, 0.55);
 107   fColl->CreateTH1F("h_Lambda_mass", "Inv. Mass of the ", 100, 1.08, 1.14);
 108   fColl->CreateTH1F("h_ALambda_mass", "Inv. Mass of the ", 100, 1.08, 1.14);
 109
 110   // total number of tagged daughter particles (should correlate with number of V0s !
 111   fColl->CreateTH1F("h_NumberOfDaughters", "Number of tagged daughters; type; counts", 3, -0.5, 2.5);
 112   // pT spectra of tagged daughter particles
 113   fColl->CreateTH1F("h_Electron_pt", "p_{T} spectrum of tagged; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 114   fColl->CreateTH1F("h_Pion_pt", "p_{T} spectrum of tagged; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 115   fColl->CreateTH1F("h_Proton_pt", "p_{T} spectrum of tagged; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 116
 117   //
 118   // MC histos
 119   //
 120
 121   // pT spectra of the tagged V0s
 122   fCollMC->CreateTH1F("h_Gamma_pt_S", "MC-S: p_{T} spectrum of tagged gammas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 123   fCollMC->CreateTH1F("h_K0_pt_S", "MC-S: p_{T} spectrum of tagged K0s; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 124   fCollMC->CreateTH1F("h_Lambda_pt_S", "MC-S: p_{T} spectrum of tagged Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 125   fCollMC->CreateTH1F("h_ALambda_pt_S", "MC-S: p_{T} spectrum of tagged A-Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 126   fCollMC->CreateTH1F("h_Gamma_pt_B", "MC-B: p_{T} spectrum of tagged gammas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 127   fCollMC->CreateTH1F("h_K0_pt_B", "MC-B: p_{T} spectrum of tagged K0s; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 128   fCollMC->CreateTH1F("h_Lambda_pt_B", "MC-B: p_{T} spectrum of tagged Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 129   fCollMC->CreateTH1F("h_ALambda_pt_B", "MC-B: p_{T} spectrum of tagged A-Lambdas; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 130   // invariant mass of the V0s
 131   fCollMC->CreateTH1F("h_Gamma_mass_S", "MC-S: Inv. Mass of the gamma; m (GeV/c^{2}); counts", 100, 0., 0.1);
 132   fCollMC->CreateTH1F("h_K0_mass_S", "MC-S: Inv. Mass of the K0; m (GeV/c^{2}); counts", 100, 0.45, 0.55);
 133   fCollMC->CreateTH1F("h_Lambda_mass_S", "MC-S: Inv. Mass of the Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
 134   fCollMC->CreateTH1F("h_ALambda_mass_S", "MC-S: Inv. Mass of the A-Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
 135   fCollMC->CreateTH1F("h_Gamma_mass_B", "MC-B: Inv. Mass of the gamma; m (GeV/c^{2}); counts", 100, 0., 0.1);
 136   fCollMC->CreateTH1F("h_K0_mass_B", "MC-B: Inv. Mass of the K0; m (GeV/c^{2}); counts", 100, 0.45, 0.55);
 137   fCollMC->CreateTH1F("h_Lambda_mass_B", "MC-B: Inv. Mass of the Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
 138   fCollMC->CreateTH1F("h_ALambda_mass_B", "MC-B: Inv. Mass of the A-Lambda; m (GeV/c^{2}); counts", 100, 1.08, 1.14);
 139   // pT spectra of tagged daughter particles
 140   fCollMC->CreateTH1F("h_Electron_pt_S", "MC-S: p_{T} spectrum of tagged electrons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 141   fCollMC->CreateTH1F("h_Pion_pt_S", "MC-S: p_{T} spectrum of tagged pions; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 142   fCollMC->CreateTH1F("h_Proton_pt_S", "MC-S: p_{T} spectrum of tagged protons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 143   fCollMC->CreateTH1F("h_Electron_pt_B", "MC-B: p_{T} spectrum of tagged electrons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 144   fCollMC->CreateTH1F("h_Pion_pt_B", "MC-B: p_{T} spectrum of tagged pions; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 145   fCollMC->CreateTH1F("h_Proton_pt_B", "MC-B: p_{T} spectrum of tagged protons; p_{T} (GeV/c); counts", 20, 0.1, 20, 0);
 146
 147
 148   TList *tmp = fColl->GetList();
 149   tmp->SetName(fColl->GetName());
 150   fOutput->Add(tmp);
 151   tmp = 0x0;
 152   tmp = fCollMC->GetList();
 153   tmp->SetName(fCollMC->GetName());
 154   fOutput->Add(tmp);
 155
 156
 157 }
 158 //__________________________________________________________
 159 void AliAnalysisTaskCheckV0tender::UserExec(Option_t *){
 160   //
 161   // Event Loop
 162   //
 163   AliMCEventHandler* mcHandler = (dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()));
 164   AliESDInputHandler *inh = dynamic_cast<AliESDInputHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 165   AliESDpid *workingPID = NULL;
 166   if(inh && (workingPID = inh->GetESDpid())) workingPID = inh->GetESDpid();
 167   else workingPID = AliHFEtools::GetDefaultPID(mcHandler ? kTRUE : kFALSE);
 168
 169   // check the MC data
 170   if(fMCEvent && !mcHandler ) return;
 171   if(fMCEvent &&  !mcHandler->InitOk() ) return;
 172   if(fMCEvent &&  !mcHandler->TreeK() ) return;
 173   if(fMCEvent &&  !mcHandler->TreeTR() ) return;
 174
 175   ProcessV0s();
 176   ProcessDaughters();
 177
 178   if(fMCEvent){
 179     ProcessV0sMC();
 180     ProcessDaughtersMC();
 181   }
 182
 183   fEvents->Fill(1.1);
 184   PostData(1, fEvents);
 185   PostData(2, fOutput);
 186 }
 187 //__________________________________________________________
 188 void AliAnalysisTaskCheckV0tender::Terminate(Option_t *){
 189   //
 190   // Do Post Processing
 191   //
 192 }
 193 //__________________________________________________________
 194 void AliAnalysisTaskCheckV0tender::ProcessV0s(){
 195   //
 196   // loop over the V0s and extract the information about
 197   // the V0s tagged by the V0 tender
 198   //
 199
 200   const char * type[4] = {"Gamma", "K0", "Lambda", "ALambda"};
 201   char name[256] = "";
 202
 203   Int_t nV0s = fInputEvent->GetNumberOfV0s();
 204   for(Int_t i=0; i<nV0s; ++i){
 205     AliESDv0 *esdV0 = (dynamic_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
 206     if(!esdV0) continue;
 207     if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
 208     Int_t pid = GetTenderPidV0(esdV0);
 209     if(pid < 0) continue;
 210     fColl->Fill("h_NumberOf_V0s", pid);
 211     sprintf(name, "h_%s_pt", type[pid]);
 212     Float_t pT = esdV0->Pt();
 213     fColl->Fill(name, pT);
 214     Float_t mass = MassV0(esdV0, pid);
 215     sprintf(name, "h_%s_mass", type[pid]);
 216     fColl->Fill(name, mass);
 217   }
 218
 219 }
 220 //__________________________________________________________
 221 void AliAnalysisTaskCheckV0tender::ProcessDaughters(){
 222   //
 223   // produce some check plots for V0 tender tagged single tracks
 224   //
 225
 226   const char * type[3] = {"Electron", "Pion", "Proton"};
 227   char name[256] = "";
 228
 229   Int_t nTracks = fInputEvent->GetNumberOfTracks();
 230   for(Int_t i=0; i<nTracks; ++i){
 231     AliESDtrack *track = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(i));
 232     if(!track) continue;
 233     Int_t  pid = GetTenderPidDaughter(track);
 234     if(pid < 0) continue;
 235     fColl->Fill("h_NumberOfDaughters", pid*1.0);
 236     Float_t pT = track->Pt();
 237     sprintf(name, "h_%s_pt", type[pid]);
 238     fColl->Fill(name, pT);
 239
 240   }
 241 }
 242 //__________________________________________________________
 243 void AliAnalysisTaskCheckV0tender::ProcessV0sMC(){
 244   //
 245   // check all V0tender selected V0 on their true identity
 246   //
 247
 248   const Int_t pid2pdg[4] = {22, 310, 3122, -3122};
 249   const char * type[4] = {"Gamma", "K0", "Lambda", "ALambda"};
 250   char name[256] = "";
 251   Int_t nV0s = fInputEvent->GetNumberOfV0s();
 252
 253   Int_t nTracks = fInputEvent->GetNumberOfTracks();
 254
 255   // V0 loop
 256   for(Int_t i=0; i<nV0s; ++i){
 257     Bool_t id = kFALSE;
 258     AliESDv0 *esdV0 = (dynamic_cast<AliESDEvent *>(fInputEvent))->GetV0(i);
 259     if(!esdV0) continue;
 260     if(!esdV0->GetOnFlyStatus()) continue; // Take only V0s from the On-the-fly v0 finder
 261     Int_t pid = GetTenderPidV0(esdV0);
 262     if(pid < 0) continue;
 263
 264     // both ESD daughtr tracks
 265     Int_t iN, iP;
 266     iN = iP = -1;
 267     iP = esdV0->GetPindex();
 268     iN = esdV0->GetNindex();
 269     if(iN < 0 || iP < 0) continue;
 270     if(iN >= nTracks || iP >= nTracks) continue;
 271     AliESDtrack *dP, *dN;
 272     dP = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iP));
 273     dN = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(iN));
 274     if(!dN || !dP) continue;
 275
 276     // MC labels of the daughter tracks
 277     Int_t lN, lP;
 278     lN = dN->GetLabel();
 279     lP = dP->GetLabel();
 280     if(lN < 0 || lP < 0) continue;
 281
 282     // MC daughter particles
 283     AliMCParticle *mcP, *mcN;
 284     mcP = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lP));
 285     mcN = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lN));
 286     if(!mcP || !mcN) continue;
 287
 288     // labels of the mother particles
 289     Int_t lPm, lNm;
 290     lPm = mcP->GetMother();
 291     lNm = mcN->GetMother();
 292     if(lPm < 0) continue;
 293     AliMCParticle *m = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(lPm));
 294     if(!m) continue;
 295     Int_t pdg = m->PdgCode();
 296     if((lPm == lNm) && (pdg == pid2pdg[pid])) id = kTRUE;
 297
 298     if(id) sprintf(name, "h_%s_pt_S", type[pid]);
 299     else sprintf(name, "h_%s_pt_B", type[pid]);
 300     Float_t pT = esdV0->Pt();
 301     fCollMC->Fill(name, pT);
 302
 303     if(id) sprintf(name, "h_%s_mass_S", type[pid]);
 304     else sprintf(name, "h_%s_mass_B", type[pid]);
 305     Float_t mass = MassV0(esdV0, pid);
 306     fCollMC->Fill(name, mass);
 307
 308    }
 309 }
 310 //__________________________________________________________
 311 void AliAnalysisTaskCheckV0tender::ProcessDaughtersMC(){
 312   //
 313   // check the identity of the V0tender selected V0 daughters
 314   // !!! for positive check only the true identity plays a role here,
 315   // not the true V0 mother identity (e.g. selected electron could come
 316   // from primary vertex or pion dalitz deca or true gamma conversion) !!!
 317   //
 318
 319   const Int_t pid2pdg [3] = {11, 211, 2212};
 320   const char * type[3] = {"Electron", "Pion", "Proton"};
 321   char name[256] = "";
 322
 323
 324   Int_t nTracks = fInputEvent->GetNumberOfTracks();
 325   for(Int_t i=0; i<nTracks; ++i){
 326     Bool_t id = kFALSE;
 327     AliESDtrack *track = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(i));
 328     if(!track) continue;
 329     Int_t  pid = GetTenderPidDaughter(track);
 330     if(pid < 0) continue;
 331     Float_t pT = track->Pt();
 332     Int_t label = track->GetLabel();
 333     if(label < 0) continue;
 334     AliMCParticle *mcp = 0x0;
 335     mcp = dynamic_cast<AliMCParticle*>(fMCEvent->GetTrack(label));
 336     if(!mcp) continue;
 337     Int_t pdg = TMath::Abs(mcp->PdgCode());
 338     if(pdg == pid2pdg[pid]) id = kTRUE;
 339     if(id) sprintf(name, "h_%s_pt_S", type[pid]);
 340     else sprintf(name, "h_%s_pt_B", type[pid]);
 341     fCollMC->Fill(name, pT);
 342   }
 343 }
 344 //__________________________________________________________
 345 Int_t AliAnalysisTaskCheckV0tender::GetTenderPidV0(AliESDv0 * const v0){
 346   //
 347   // retrieve the PID nformation stored in the status flags by the train tender
 348   //
 349   Int_t pid = -1;
 350   if(!v0){
 351     return pid;
 352   }
 353   Int_t nTimes = 0;
 354   if(v0->TestBit(BIT(14))){
 355     pid = 0;
 356     nTimes++;
 357   }
 358   if(v0->TestBit(BIT(15))){
 359     pid = 1;
 360     nTimes++;
 361   }
 362   if(v0->TestBit(BIT(16))){
 363     pid = 2;
 364     nTimes++;
 365   }
 366   if(v0->TestBit(BIT(17))){
 367     pid = 3;
 368     nTimes++;
 369   }
 370   if(nTimes > 1){
 371     AliWarning("V0 track labeled multiple times by the V0 tender");
 372     pid = -1;
 373   }
 374
 375   //printf(" -D: pid: %i \n", pid);
 376
 377   return pid;
 378 }
 379 //__________________________________________________________
 380 Int_t AliAnalysisTaskCheckV0tender::GetTenderPidDaughter(AliESDtrack * const track){
 381   //
 382   // retrieve the PID nformation stored in the status flags by the train tender
 383   //
 384
 385   Int_t pid = -1;
 386   if(!track){
 387     return pid;
 388   }
 389   Int_t nTimes = 0;
 390   if(track->TestBit(BIT(14))){
 391     pid = 0;
 392     nTimes++;
 393   }
 394   if(track->TestBit(BIT(15))){
 395     pid = 1;
 396     nTimes++;
 397   }
 398   if(track->TestBit(BIT(16))){
 399     pid = 2;
 400     nTimes++;
 401   }
 402   if(nTimes > 1){
 403     AliWarning("V0 track labeled multiple times by the V0 tender");
 404     pid = -1;
 405   }
 406   return pid;
 407 }
 408 //__________________________________________________________
 409 Float_t AliAnalysisTaskCheckV0tender::MassV0(AliESDv0 * const v0, Int_t id){
 410   //
 411   // Get the V0 effective mass
 412   //
 413
 414   Float_t mass = -0.1;
 415   Bool_t sign = CheckSigns(v0);
 416   if(0 == id){
 417     mass = v0->GetEffMass(0, 0);
 418   }
 419   else if(1 == id){
 420     mass = v0->GetEffMass(2, 2);
 421   }
 422   else if(2 == id){
 423     mass = (sign) ? v0->GetEffMass(4, 2) : v0->GetEffMass(2, 4);
 424   }
 425   else if(3 == id){
 426     mass = (sign) ? v0->GetEffMass(2, 4) : v0->GetEffMass(4, 2);
 427   }
 428   else{
 429     AliWarning(Form("Unrecognized V0 id: %i", id));
 430   }
 431
 432   return mass;
 433
 434 }
 435 //__________________________________________________________
 436 Bool_t AliAnalysisTaskCheckV0tender::CheckSigns(AliESDv0 * const v0){
 437   //
 438   // check wheter the sign was correctly applied to
 439   // V0 daughter tracks
 440   // This function should become obsolete once the V0 finder will be updated (fixed)
 441   //
 442
 443   Bool_t correct = kFALSE;
 444
 445   Int_t pIndex = 0, nIndex = 0;
 446   pIndex = v0->GetPindex();
 447   nIndex = v0->GetNindex();
 448
 449   AliESDtrack* d[2];
 450   d[0] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(pIndex));
 451   d[1] = dynamic_cast<AliESDtrack*>(fInputEvent->GetTrack(nIndex));
 452
 453   Int_t sign[2];
 454   sign[0] = (int)d[0]->GetSign();
 455   sign[1] = (int)d[1]->GetSign();
 456
 457   if(-1 == sign[0] && 1 == sign[1]){
 458     correct = kFALSE;
 459     //v0->SetIndex(0, pIndex);  // set the index of the negative v0 track
 460     //v0->SetIndex(1, nIndex);  // set the index of the positive v0 track
 461   }
 462   else{
 463     correct = kTRUE;
 464   }
 465
 466   //pIndex = v0->GetPindex();
 467   //nIndex = v0->GetNindex();
 468   //printf("-D2: P: %i, N: %i\n", pIndex, nIndex);
 469
 470   return correct;
 471 }