PWG2/RESONANCES/AliRsnReader.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 //
  17 // ==== Class AliRsnReader ========
  18 //
  19 // This object reads a 'standard' event and converts it into the internal
  20 // format used for resonance analysis (AliRsnEvent).
  21 // 'Standard' event means ESD, standard AOD and MC event.
  22 //
  23 // The input-2-AliRsnEvent conversion is done through a class which reads
  24 // from AliAnalysisTaskSE, which is the standard analysis object.
  25 // This class creates the AliRsnEvent's before the input event is read,
  26 // so this class has not to 'create' a new outpu event, but instead it has
  27 // to 'fill' one which has already been created elsewhere.
  28 // Then, the methods provided here accept an AliRsnEvent as argument passed
  29 // by reference, and they 'fill' this object using the data from the inputs
  30 // passed to them.
  31 //
  32 // author: A. Pulvirenti
  33 // email : alberto.pulvirenti@ct.infn.it
  34 //
  35
  36 #include <TString.h>
  37
  38 #include "AliLog.h"
  39
  40 #include "AliVEvent.h"
  41
  42 #include "AliESDEvent.h"
  43 #include "AliESDtrack.h"
  44 #include "AliESDVertex.h"
  45
  46 #include "AliAODEvent.h"
  47 #include "AliAODTrack.h"
  48 #include "AliAODVertex.h"
  49
  50 #include "AliStack.h"
  51 #include "AliMCEvent.h"
  52 #include "AliMCParticle.h"
  53 #include "AliGenEventHeader.h"
  54
  55 #include "AliRsnMCInfo.h"
  56 #include "AliRsnDaughter.h"
  57 #include "AliRsnEvent.h"
  58 #include "AliRsnPIDWeightsMgr.h"
  59 #include "AliRsnReader.h"
  60
  61 ClassImp(AliRsnReader)
  62
  63 //_____________________________________________________________________________
  64 AliRsnReader::AliRsnReader(ESource source, AliRsnPIDWeightsMgr *mgr) :
  65   TObject(),
  66   fSource(source),
  67   fCheckSplit(kFALSE),
  68   fRejectFakes(kFALSE),
  69   fWeightsMgr(mgr)
  70 {
  71 //
  72 // Constructor.
  73 // Initializes the base-type data members:
  74 //   - management of fake tracks
  75 //
  76 }
  77
  78 //_____________________________________________________________________________
  79 AliRsnReader::AliRsnReader(const AliRsnReader &copy) :
  80   TObject(copy),
  81   fSource(copy.fSource),
  82   fCheckSplit(copy.fCheckSplit),
  83   fRejectFakes(copy.fRejectFakes),
  84   fWeightsMgr(copy.fWeightsMgr)
  85 {
  86 //
  87 // Copy constructor.
  88 //
  89 }
  90
  91 //_____________________________________________________________________________
  92 AliRsnReader& AliRsnReader::operator=(const AliRsnReader &copy)
  93 {
  94 //
  95 // Assignment operator.
  96 //
  97
  98     fSource = copy.fSource;
  99     fCheckSplit = copy.fCheckSplit;
 100     fRejectFakes = copy.fRejectFakes;
 101     fWeightsMgr = copy.fWeightsMgr;
 102     return (*this);
 103 }
 104
 105 //_____________________________________________________________________________
 106 Bool_t AliRsnReader::Fill(AliRsnEvent *rsn, AliVEvent *event, AliMCEvent *mc)
 107 {
 108 //
 109 // According to the class type of event and the selected source
 110 // recalls one of the private reading methods to fill the RsnEvent
 111 // passed by reference as first argument.
 112 //
 113
 114     TString str(event->ClassName());
 115
 116     switch (fSource) {
 117         case kESD:
 118             if (!str.Contains("AliESDEvent")) {
 119                 AliError("A reader set to 'kESD' or 'kESDTPC' can read only ESD events");
 120                 return kFALSE;
 121             }
 122             return FillFromESD(rsn, (AliESDEvent*)event, mc);
 123         case kESDTPC:
 124             if (!str.Contains("AliESDEvent")) {
 125                 AliError("A reader set to 'kESD' or 'kESDTPC' can read only ESD events");
 126                 return kFALSE;
 127             }
 128             return FillFromESD(rsn, (AliESDEvent*)event, mc, kTRUE);
 129         case kAOD:
 130             if (!str.Contains("AliAODEvent")) {
 131                 AliError("A reader set to 'kAOD' can read only AOD events");
 132                 return kFALSE;
 133             }
 134             return FillFromAOD(rsn, (AliAODEvent*)event, mc);
 135         case kMC:
 136             if (!str.Contains("AliMCEvent")) {
 137                 AliError("A reader set to 'kMC' can read only MC events");
 138                 return kFALSE;
 139             }
 140             return FillFromMC(rsn, (AliMCEvent*)event);
 141         default:
 142             return kFALSE;
 143     }
 144 }
 145
 146 //_____________________________________________________________________________
 147 Bool_t AliRsnReader::FillFromESD
 148 (AliRsnEvent *rsn, AliESDEvent *esd, AliMCEvent *mc, Bool_t useTPCOnly)
 149 {
 150 //
 151 // Filler from an ESD event.
 152 // Stores all tracks (if a filter is defined, it will store
 153 // only the ones which survive the cuts).
 154 // If a reference MC event is provided, it is used to store
 155 // the MC informations for each track (true PDG code,
 156 // GEANT label of mother, PDG code of mother, if any).
 157 // When this is used, the 'source' flag of the output
 158 // AliRsnEvent object will be set to 'kESD'.
 159 //
 160
 161     // retrieve stack (if possible)
 162     AliStack *stack = 0x0;
 163     if (mc) stack = mc->Stack();
 164
 165     // get number of tracks
 166     Int_t ntracks = esd->GetNumberOfTracks();
 167     /*
 168     if (!ntracks) {
 169        AliWarning("No tracks in this event");
 170        return kFALSE;
 171     }
 172     */
 173
 174     // if required with the flag, scans the event
 175     // and searches all split tracks (= 2 tracks with the same label);
 176     // for each pair of split tracks, only the better (best chi2) is kept
 177     // and the other is rejected: this info is stored into a Boolean array
 178     Int_t i1, i2, lab1, lab2;
 179     Bool_t *accept = new Bool_t[ntracks];
 180     for (i1 = 0; i1 < ntracks; i1++) accept[i1] = kTRUE;
 181     if (fCheckSplit) {
 182         for (i1 = 0; i1 < ntracks; i1++) {
 183             AliESDtrack *trk1 = esd->GetTrack(i1);
 184             lab1 = TMath::Abs(trk1->GetLabel());
 185             for (i2 = i1+1; i2 < ntracks; i2++) {
 186                 AliESDtrack *trk2 = esd->GetTrack(i2);
 187                 lab2 = TMath::Abs(trk2->GetLabel());
 188                 // check if labels are equal
 189                 if (lab1 == lab2) {
 190                     if (trk1->GetConstrainedChi2() < trk2->GetConstrainedChi2()) {
 191                         accept[i1] = kTRUE;
 192                         accept[i2] = kFALSE;
 193                     }
 194                     else {
 195                         accept[i1] = kFALSE;
 196                         accept[i2] = kTRUE;
 197                     }
 198                 }
 199             }
 200         }
 201     }
 202
 203     // get primary vertex
 204     Double_t vertex[3];
 205     if (!useTPCOnly) {
 206         vertex[0] = esd->GetVertex()->GetXv();
 207         vertex[1] = esd->GetVertex()->GetYv();
 208         vertex[2] = esd->GetVertex()->GetZv();
 209     }
 210     else {
 211         vertex[0] = esd->GetPrimaryVertexTPC()->GetXv();
 212         vertex[1] = esd->GetPrimaryVertexTPC()->GetYv();
 213         vertex[2] = esd->GetPrimaryVertexTPC()->GetZv();
 214     }
 215     rsn->SetPrimaryVertex(vertex[0], vertex[1], vertex[2]);
 216
 217     // store tracks from ESD
 218     Int_t  i, index, label, labmum;
 219     Bool_t check;
 220     AliRsnDaughter temp;
 221     for (index = 0; index < ntracks; index++) {
 222         // skip track recognized as the worse one in a splitted pair
 223         if (!accept[index]) {
 224             AliInfo(Form("Rejecting split track #%d in this event", index));
 225             continue;
 226         }
 227         // get and check ESD track
 228         AliESDtrack *esdTrack = esd->GetTrack(index);
 229         label = esdTrack->GetLabel();
 230         if (fRejectFakes && (label < 0)) continue;
 231         // copy ESD track data into RsnDaughter
 232         // if unsuccessful, this track is skipped
 233         check = temp.Adopt(esdTrack, useTPCOnly);
 234         if (!check) continue;
 235         // if the AliRsnWeightsMgr object is initialized
 236         // this means that ESD PID weights are not used
 237         // and they are computed according to the Weight manager
 238         if (fWeightsMgr) {
 239             //AliInfo("Using customized weights");
 240             //AliInfo(Form("ESD weights = %f %f %f %f %f", temp.PID()[0], temp.PID()[1], temp.PID()[2], temp.PID()[3], temp.PID()[4], temp.PID()[5]));
 241             esdTrack->GetITSpid(fWeightsMgr->GetWeightArray(AliRsnPIDWeightsMgr::kITS));
 242             esdTrack->GetTPCpid(fWeightsMgr->GetWeightArray(AliRsnPIDWeightsMgr::kTPC));
 243             esdTrack->GetTRDpid(fWeightsMgr->GetWeightArray(AliRsnPIDWeightsMgr::kTRD));
 244             esdTrack->GetTOFpid(fWeightsMgr->GetWeightArray(AliRsnPIDWeightsMgr::kTOF));
 245             esdTrack->GetHMPIDpid(fWeightsMgr->GetWeightArray(AliRsnPIDWeightsMgr::kHMPID));
 246             for (i = 0; i < AliRsnPID::kSpecies; i++) {
 247                 temp.SetPIDWeight(i, fWeightsMgr->GetWeight((AliRsnPID::EType)i, temp.Pt()));
 248             }
 249             //AliInfo(Form("Used weights = %f %f %f %f %f", temp.PID()[0], temp.PID()[1], temp.PID()[2], temp.PID()[3], temp.PID()[4], temp.PID()[5]));
 250         }
 251         else {
 252             //AliInfo("Using standard ESD weights");
 253         }
 254
 255         // if stack is present, copy MC info
 256         if (stack) {
 257             TParticle *part = stack->Particle(TMath::Abs(label));
 258             if (part) {
 259                 temp.InitMCInfo(part);
 260                 labmum = part->GetFirstMother();
 261                 if (labmum >= 0) {
 262                     TParticle *mum = stack->Particle(labmum);
 263                     temp.GetMCInfo()->SetMotherPDG(mum->GetPdgCode());
 264                 }
 265             }
 266         }
 267         // set index and label and add this object to the output container
 268         temp.SetIndex(index);
 269         temp.SetLabel(label);
 270         AliRsnDaughter *ptr = rsn->AddTrack(temp);
 271         // if problems occurred while storing, that pointer is NULL
 272         if (!ptr) AliWarning(Form("Failed storing track#%d", index));
 273     }
 274
 275     // compute total multiplicity
 276     if (rsn->GetMultiplicity() <= 0) {
 277         AliWarning("Zero Multiplicity in this event");
 278         return kFALSE;
 279     }
 280
 281     return kTRUE;
 282 }
 283
 284 //_____________________________________________________________________________
 285 Bool_t AliRsnReader::FillFromAOD(AliRsnEvent *rsn, AliAODEvent *aod, AliMCEvent *mc)
 286 {
 287 //
 288 // Filler from an AOD event.
 289 // Stores all tracks (if a filter is defined, it will store
 290 // only the ones which survive the cuts).
 291 // If a reference MC event is provided, it is used to store
 292 // the MC informations for each track (true PDG code,
 293 // GEANT label of mother, PDG code of mother, if any).
 294 // When this is used, the 'source' flag of the output
 295 // AliRsnEvent object will be set to 'kAOD'.
 296 //
 297
 298     // retrieve stack (if possible)
 299     AliStack *stack = 0x0;
 300     if (mc) stack = mc->Stack();
 301
 302     // get number of tracks
 303     Int_t ntracks = aod->GetNTracks();
 304     if (!ntracks) {
 305        AliWarning("No tracks in this event");
 306        return kFALSE;
 307     }
 308
 309     // get primary vertex
 310     Double_t vertex[3];
 311     vertex[0] = aod->GetPrimaryVertex()->GetX();
 312     vertex[1] = aod->GetPrimaryVertex()->GetY();
 313     vertex[2] = aod->GetPrimaryVertex()->GetZ();
 314     rsn->SetPrimaryVertex(vertex[0], vertex[1], vertex[2]);
 315
 316     // store tracks from ESD
 317     Int_t  index, label, labmum;
 318     Bool_t check;
 319     AliAODTrack *aodTrack = 0;
 320     AliRsnDaughter temp;
 321     TObjArrayIter iter(aod->GetTracks());
 322     while ( (aodTrack = (AliAODTrack*)iter.Next()) ) {
 323         // retrieve index
 324         index = aod->GetTracks()->IndexOf(aodTrack);
 325         label = aodTrack->GetLabel();
 326         if (fRejectFakes && (label < 0)) continue;
 327         // copy ESD track data into RsnDaughter
 328         // if unsuccessful, this track is skipped
 329         check = temp.Adopt(aodTrack);
 330         if (!check) continue;
 331         // if stack is present, copy MC info
 332         if (stack) {
 333             TParticle *part = stack->Particle(TMath::Abs(label));
 334             if (part) {
 335                 temp.InitMCInfo(part);
 336                 labmum = part->GetFirstMother();
 337                 if (labmum >= 0) {
 338                     TParticle *mum = stack->Particle(labmum);
 339                     temp.GetMCInfo()->SetMotherPDG(mum->GetPdgCode());
 340                 }
 341             }
 342         }
 343         // set index and label and add this object to the output container
 344         temp.SetIndex(index);
 345         temp.SetLabel(label);
 346         AliRsnDaughter *ptr = rsn->AddTrack(temp);
 347         // if problems occurred while storin, that pointer is NULL
 348         if (!ptr) AliWarning(Form("Failed storing track#%d"));
 349     }
 350
 351     // compute total multiplicity
 352     if (rsn->GetMultiplicity() <= 0) {
 353         AliWarning("Zero multiplicity in this event");
 354         return kFALSE;
 355     }
 356
 357     return kTRUE;
 358 }
 359
 360 //_____________________________________________________________________________
 361 Bool_t AliRsnReader::FillFromMC(AliRsnEvent *rsn, AliMCEvent *mc)
 362 {
 363 //
 364 // Filler from an ESD event.
 365 // Stores all tracks which generate at least one
 366 // TrackReference (a point in a sensitive volume).
 367 // In this case, the MC info is stored by default and
 368 // perfect particle identification is the unique available.
 369 // When this is used, the 'source' flag of the output
 370 // AliRsnEvent object will be set to 'kMC'.
 371 //
 372
 373     // get number of tracks
 374     Int_t ntracks = mc->GetNumberOfTracks();
 375     if (!ntracks) {
 376        AliWarning("No tracks in this event");
 377        return kFALSE;
 378     }
 379
 380     AliStack *stack = mc->Stack();
 381
 382     // get primary vertex
 383     TArrayF fvertex(3);
 384     Double_t vertex[3];
 385     mc->GenEventHeader()->PrimaryVertex(fvertex);
 386     vertex[0] = (Double_t)fvertex[0];
 387     vertex[1] = (Double_t)fvertex[1];
 388     vertex[2] = (Double_t)fvertex[2];
 389     rsn->SetPrimaryVertex(vertex[0], vertex[1], vertex[2]);
 390
 391     // store tracks from MC
 392     Int_t  index, labmum;
 393     Bool_t check;
 394     AliRsnDaughter temp;
 395     for (index = 0; index < ntracks; index++) {
 396         // get and check MC track
 397         AliMCParticle *mcTrack = mc->GetTrack(index);
 398         // if particle has no track references, it is rejected
 399         if (mcTrack->GetNumberOfTrackReferences() <= 0) continue;
 400         // try to insert in the RsnDaughter its data
 401         check = temp.Adopt(mcTrack);
 402         if (!check) continue;
 403         labmum = temp.GetMCInfo()->Mother();
 404         if (labmum >= 0) {
 405             TParticle *mum = stack->Particle(labmum);
 406             temp.GetMCInfo()->SetMotherPDG(mum->GetPdgCode());
 407         }
 408         // if successful, set other data and stores it
 409         temp.SetIndex(index);
 410         temp.SetLabel(mcTrack->Label());
 411         AliRsnDaughter *ptr = rsn->AddTrack(temp);
 412         // if problems occurred while storin, that pointer is NULL
 413         if (!ptr) AliWarning(Form("Failed storing track#%d", index));
 414     }
 415
 416     // compute total multiplicity
 417     if (rsn->GetMultiplicity() <= 0) {
 418         AliWarning("Zero multiplicity in this event");
 419         return kFALSE;
 420     }
 421
 422     return kTRUE;
 423 }