PWGHF/hfe/AliHFEreducedEventCreatorESD.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 // Debug tree task
  17 // the tree is represented as reduced events
  18 //
  19 // Authors:
  20 //   M.Fasel <M.Fasel@gsi.de>
  21 //
  22 //
  23 #include <TArrayI.h>
  24 #include <TBits.h>
  25 #include <TFile.h>
  26 #include <TParticle.h>
  27 #include <TString.h>
  28 #include <TTree.h>
  29
  30 #include "AliAnalysisManager.h"
  31 #include "AliAODMCHeader.h"
  32 #include "AliAnalysisUtils.h"
  33 #include "AliCentrality.h"
  34 #include "AliESDEvent.h"
  35 #include "AliESDtrack.h"
  36 #include "AliESDVertex.h"
  37 #include "AliInputEventHandler.h"
  38 #include "AliHFEcuts.h"
  39 #include "AliHFEextraCuts.h"
  40 #include "AliHFEpidTPC.h"
  41 #include "AliHFEreducedEvent.h"
  42 #include "AliHFEreducedTrack.h"
  43 #include "AliHFEreducedMCParticle.h"
  44 #include "AliHFEsignalCuts.h"
  45 #include "AliLog.h"
  46 #include "AliMCEvent.h"
  47 #include "AliMCParticle.h"
  48 #include "AliMCEventHandler.h"
  49 #include "AliMultiplicity.h"
  50 #include "AliPIDResponse.h"
  51 #include "AliVEvent.h"
  52 #include "AliVTrack.h"
  53 #include "AliVCluster.h"
  54 #include "AliVVZERO.h"
  55 #include "AliVZDC.h"
  56 #include "AliTRDTriggerAnalysis.h"
  57 #include "TTreeStream.h"
  58
  59 #include "AliHFEreducedEventCreatorESD.h"
  60
  61 ClassImp(AliHFEreducedEventCreatorESD)
  62
  63 AliHFEreducedEventCreatorESD::AliHFEreducedEventCreatorESD():
  64   AliAnalysisTaskSE(),
  65   fHFEtree(NULL),
  66   fAnalysisUtils(NULL),
  67   fHFEevent(NULL),
  68   fTrackCuts(NULL),
  69   fExtraCuts(NULL),
  70   fSignalCuts(NULL),
  71   fTPCpid(NULL),
  72   fTRDTriggerAnalysis(NULL),
  73   fEventNumber(0),
  74   fNclustersTPC(70),
  75   fNclustersTPCPID(0),
  76   fNclustersITS(2),
  77   fRemoveFirstEvent(kFALSE)
  78 {
  79   //
  80   // Default constructor
  81   //
  82 }
  83
  84 AliHFEreducedEventCreatorESD::AliHFEreducedEventCreatorESD(const char *name):
  85   AliAnalysisTaskSE(name),
  86   fHFEtree(NULL),
  87   fAnalysisUtils(NULL),
  88   fHFEevent(NULL),
  89   fTrackCuts(NULL),
  90   fExtraCuts(NULL),
  91   fSignalCuts(NULL),
  92   fTPCpid(NULL),
  93   fTRDTriggerAnalysis(NULL),
  94   fEventNumber(0),
  95   fNclustersTPC(70),
  96   fNclustersTPCPID(0),
  97   fNclustersITS(2),
  98   fRemoveFirstEvent(kFALSE)
  99 {
 100   //
 101   // Default constructor
 102   //
 103   fTPCpid = new AliHFEpidTPC("QAtpcPID");
 104   fAnalysisUtils = new AliAnalysisUtils;
 105   fTRDTriggerAnalysis = new AliTRDTriggerAnalysis();
 106   fTRDTriggerAnalysis->SetRequireMatchElectron(kTRUE);
 107   DefineOutput(1, TTree::Class());
 108 }
 109
 110 AliHFEreducedEventCreatorESD::~AliHFEreducedEventCreatorESD(){
 111   //
 112   // Default destructor
 113   //
 114   if(fAnalysisUtils) delete fAnalysisUtils;
 115   if(fTPCpid) delete fTPCpid;
 116   if(fTRDTriggerAnalysis) delete fTRDTriggerAnalysis;
 117   if(fHFEevent) delete fHFEevent;
 118   if(fSignalCuts) delete fSignalCuts;
 119   if(fTrackCuts) delete fTrackCuts;
 120
 121 }
 122
 123 void AliHFEreducedEventCreatorESD::UserCreateOutputObjects(){
 124   //
 125   // Create debug tree, signal cuts and track cuts
 126   //
 127
 128   //printf("test\n");
 129
 130  // printf("testa\n");
 131   fSignalCuts = new AliHFEsignalCuts("HFEsignalCuts", "HFE MC Signal definition");
 132   //printf("testb\n");
 133
 134   fTrackCuts = new AliHFEcuts("fTrackCuts", "Basic HFE track cuts");
 135   fTrackCuts->CreateStandardCuts();
 136   // Track cuts
 137   fTrackCuts->SetMaxChi2perClusterTPC(1000.);   // switch off this cut (for the moment);
 138   fTrackCuts->SetMinNClustersTPC(fNclustersTPC);
 139   fTrackCuts->SetMinRatioTPCclusters(0);
 140   fTrackCuts->SetTPCmodes(AliHFEextraCuts::kFound, AliHFEextraCuts::kFoundOverFindable);
 141   fTrackCuts->SetMinNClustersTPCPID(fNclustersTPCPID);
 142   fTrackCuts->SetMinNClustersITS(fNclustersITS);
 143   // Event cuts
 144   fTrackCuts->SetUseMixedVertex(kTRUE);
 145   fTrackCuts->SetVertexRange(10.);
 146   //printf("testa\n");
 147   fTrackCuts->Initialize();
 148   //printf("testb\n");
 149
 150   fExtraCuts = new AliHFEextraCuts("hfeExtraCuts","HFE Extra Cuts");
 151
 152   fHFEevent = new AliHFEreducedEvent;
 153   OpenFile(1);
 154   fHFEtree = new TTree("HFEtree", "HFE event tree");
 155   fHFEtree->Branch("HFEevent", "AliHFEreducedEvent", fHFEevent,128000,0);
 156   PostData(1, fHFEtree);
 157 }
 158
 159 void AliHFEreducedEventCreatorESD::UserExec(Option_t *){
 160   //
 161   // User Exec: Fill debug Tree
 162   //
 163
 164   // Get PID response
 165   AliPIDResponse *pid = NULL;
 166   AliInputEventHandler *handler = dynamic_cast<AliInputEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
 167   if(handler){
 168 //    printf("testb\n");
 169     pid = handler->GetPIDResponse();
 170   } else {
 171     AliError("No Handler");
 172   }
 173   if(!pid){
 174  //   printf("testc\n");
 175     AliError("No PID response");
 176     return;
 177   }
 178   if(!fInputEvent) {
 179     AliError("No Input event");
 180     return;
 181   }
 182
 183   if(fRemoveFirstEvent){
 184     if(fAnalysisUtils->IsFirstEventInChunk(fInputEvent)) return;
 185   }
 186
 187   AliDebug(1, "Event Selected");
 188
 189   AliESDtrack copyTrack;
 190
 191   // MC info
 192   Bool_t mcthere = dynamic_cast<AliMCEventHandler *>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()) != NULL;
 193   if(mcthere){
 194     fTrackCuts->SetMCEvent(fMCEvent);
 195     fSignalCuts->SetMCEvent(fMCEvent);
 196   }
 197
 198   fTrackCuts->SetRecEvent(fInputEvent);
 199
 200   if(!fTrackCuts->CheckEventCuts("fCutsEvRec", fInputEvent)){
 201     AliDebug(1, "Event rejected by the event cuts\n");
 202     return;
 203   }
 204
 205   // reject pile up in case of pp
 206   AliESDEvent *event = dynamic_cast<AliESDEvent *>(fInputEvent);
 207   if(event) {
 208     TString beamtype = event->GetBeamType();
 209     //printf("beamtype %s\n",(const char*)beamtype);
 210     if (strstr(beamtype,"p-p")) {
 211       //printf("Reject\n");
 212       if(fInputEvent->IsPileupFromSPD(3, 0.8, 3., 2., 5)){
 213               AliDebug(1, "Event flagged as pileup\n");
 214               return;
 215       }
 216     }
 217   }
 218
 219   if(!fExtraCuts){
 220     fExtraCuts = new AliHFEextraCuts("hfeExtraCuts","HFE Extra Cuts");
 221   }
 222   fExtraCuts->SetRecEventInfo(fInputEvent);
 223
 224
 225   // Make Reduced Event
 226   //AliHFEreducedEvent hfeevent;
 227   fHFEevent->~AliHFEreducedEvent();
 228   new(fHFEevent)AliHFEreducedEvent();
 229
 230   // Get run number
 231   fHFEevent->SetRunNumber(fInputEvent->GetRunNumber());
 232
 233   // Derive trigger
 234   UInt_t trigger = fInputHandler->IsEventSelected();
 235   if(trigger & AliVEvent::kMB) fHFEevent->SetMBTrigger();
 236   if((trigger & AliVEvent::kINT7)||(trigger & AliVEvent::kINT8)) fHFEevent->SetINTTrigger();
 237   if(trigger & AliVEvent::kCentral) fHFEevent->SetCentralTrigger();
 238   if(trigger & AliVEvent::kSemiCentral) fHFEevent->SetCentralTrigger();
 239   if(trigger & AliVEvent::kEMCEJE) fHFEevent->SetEMCALTrigger();
 240
 241   fTRDTriggerAnalysis->CalcTriggers(event);
 242   if(fTRDTriggerAnalysis->HasTriggeredConfirmed(AliTRDTriggerAnalysis::kHSE)) fHFEevent->SetTRDSETrigger();
 243   if(fTRDTriggerAnalysis->HasTriggeredConfirmed(AliTRDTriggerAnalysis::kHQU)) fHFEevent->SetTRDDQTrigger();
 244
 245   // Get Primary Vertex
 246   const AliVVertex *vertex = fInputEvent->GetPrimaryVertex();
 247   Double_t vtx[3];
 248   vertex->GetXYZ(vtx);
 249   fHFEevent->SetVX(vtx[0]);
 250   fHFEevent->SetVY(vtx[1]);
 251   fHFEevent->SetVZ(vtx[2]);
 252   Int_t ncontrib(vertex->GetNContributors());
 253   fHFEevent->SetNContribVertex(ncontrib);
 254   Double_t vcov[6];
 255   vertex->GetCovarianceMatrix(vcov);
 256   fHFEevent->SetVertexResolution(TMath::Sqrt(vcov[5]));
 257   // Get Primary Vertex from SPD
 258   const AliVVertex *vertexSPD = event->GetPrimaryVertexSPD();
 259   if(vertexSPD){
 260     memset(vtx, 0, sizeof(Double_t) *3);
 261     vertexSPD->GetXYZ(vtx);
 262     fHFEevent->SetVXSPD(vtx[0]);
 263     fHFEevent->SetVYSPD(vtx[1]);
 264     fHFEevent->SetVZSPD(vtx[2]);
 265     fHFEevent->SetNContribVertexSPD(vertexSPD->GetNContributors());
 266     memset(vcov, 0, sizeof(Double_t)*6);
 267     vertex->GetCovarianceMatrix(vcov);
 268     fHFEevent->SetVertexResolutionSPD(TMath::Sqrt(vcov[5]));
 269   }
 270
 271   // Get centrality
 272   AliCentrality *hicent = fInputEvent->GetCentrality();
 273   fHFEevent->SetCentrality(
 274     hicent->GetCentralityPercentile("V0M"),
 275     hicent->GetCentralityPercentile("V0A"),
 276     hicent->GetCentralityPercentile("V0C"),
 277     hicent->GetCentralityPercentile("TKL"),
 278     hicent->GetCentralityPercentile("TRK"),
 279     hicent->GetCentralityPercentile("ZNA"),
 280     hicent->GetCentralityPercentile("ZNC"),
 281     hicent->GetCentralityPercentile("CL0"),
 282     hicent->GetCentralityPercentile("CL1"),
 283     hicent->GetCentralityPercentile("CND")
 284   );
 285
 286   // Get VZERO Information
 287   AliVVZERO *vzeroinfo = fInputEvent->GetVZEROData();
 288   if(vzeroinfo) fHFEevent->SetV0Multiplicity(vzeroinfo->GetMTotV0A(), vzeroinfo->GetMTotV0C());
 289
 290   // Get ZDC Information
 291   AliVZDC *zdcinfo = fInputEvent->GetZDCData();
 292   if(zdcinfo) fHFEevent->SetZDCEnergy(zdcinfo->GetZNAEnergy(), zdcinfo->GetZNCEnergy(), zdcinfo->GetZPAEnergy(), zdcinfo->GetZPCEnergy());
 293
 294   // Set SPD multiplicity
 295   const AliMultiplicity *mult = event->GetMultiplicity();
 296   if(mult) fHFEevent->SetSPDMultiplicity(mult->GetNumberOfTracklets());
 297
 298   //
 299   // Loop on MC tracks only
 300   //
 301   AliMCParticle *mctrack(NULL);
 302   // Monte-Carlo info
 303   Int_t source(5);
 304   if(mcthere){
 305     for(Int_t itrack = 0; itrack < fMCEvent->GetNumberOfTracks(); itrack++) {
 306       mctrack = (AliMCParticle *)(fMCEvent->GetTrack(itrack));
 307       if(!mctrack) continue;
 308       AliHFEreducedMCParticle hfemcpart;
 309       if(!fTrackCuts->CheckParticleCuts(static_cast<UInt_t>(AliHFEcuts::kStepMCGenerated), mctrack))  continue;
 310       hfemcpart.SetSignal();
 311       // Kinematics
 312       hfemcpart.SetSignedPt(mctrack->Pt(), mctrack->Charge() > 0.);
 313       hfemcpart.SetP(mctrack->P());
 314       hfemcpart.SetEta(mctrack->Eta());
 315       hfemcpart.SetPhi(mctrack->Phi());
 316       hfemcpart.SetPdg(mctrack->PdgCode());
 317
 318       // Get Production Vertex in radial direction
 319       hfemcpart.SetProductionVertex(mctrack->Xv(),mctrack->Yv(),mctrack->Zv());
 320
 321       // Get Mother PDG code of the particle
 322       Int_t motherlabel = TMath::Abs(mctrack->GetMother());
 323       if(motherlabel >= 0 && motherlabel < fMCEvent->GetNumberOfTracks()){
 324         AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(motherlabel));
 325         if(mother) hfemcpart.SetMotherPdg(mother->PdgCode());
 326       }
 327
 328       // derive source
 329       source = 5;
 330       if(fSignalCuts->IsCharmElectron(mctrack)) source = 0;
 331       else if(fSignalCuts->IsBeautyElectron(mctrack)) source = 1;
 332       else if(fSignalCuts->IsGammaElectron(mctrack)) source = 2;
 333       else if(fSignalCuts->IsNonHFElectron(mctrack)) source = 3;
 334       else if(TMath::Abs(mctrack->PdgCode()) == 11) source = 4;
 335       else source = 5;
 336       hfemcpart.SetSource(source);
 337
 338       fHFEevent->AddMCParticle(&hfemcpart);
 339     }
 340   }
 341
 342   //
 343   // Loop on reconstructed tracks
 344   //
 345   TArrayI arraytrack(fInputEvent->GetNumberOfTracks());
 346   Int_t counterdc=0;
 347
 348   AliESDtrack *track = 0x0;
 349   for(Int_t itrack = 0; itrack < fInputEvent->GetNumberOfTracks(); itrack++){
 350     // Run track loop
 351     track = dynamic_cast<AliESDtrack *>(fInputEvent->GetTrack(itrack));
 352     if(!track) continue;
 353     // Cut track (Only basic track cuts)
 354     // printf("testv\n");
 355     if(!fTrackCuts->CheckParticleCuts(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepRecKineITSTPC, track)) continue;
 356     //
 357     //printf("testu\n");
 358
 359     // Kinematics
 360     AliHFEreducedTrack hfetrack;
 361     hfetrack.SetSignedPt(track->Pt(), track->Charge() > 0);
 362     hfetrack.SetP(track->P());
 363     hfetrack.SetEta(track->Eta());
 364     hfetrack.SetPhi(track->Phi());
 365     hfetrack.SetTPCmomentum(track->GetTPCmomentum());
 366
 367     // Track ID
 368     hfetrack.SetTrackID(track->GetID());
 369
 370     // status
 371     ULong_t status = track->GetStatus();
 372     if((status & AliVTrack::kITSrefit) == AliVTrack::kITSrefit) hfetrack.SetITSrefit();
 373     if((status & AliVTrack::kTPCrefit) == AliVTrack::kTPCrefit) hfetrack.SetTPCrefit();
 374     if((status & AliVTrack::kTOFpid) == AliVTrack::kTOFpid) hfetrack.SetTOFpid();
 375     //if((status & AliVTrack::kTOFmismatch) == AliVTrack::kTOFmismatch) hfetrack.SetTOFmismatch();
 376     if(IsTOFmismatch(track, pid)) hfetrack.SetTOFmismatch(); // New version suggested by Pietro Antonioli
 377     Bool_t isEMCAL(kFALSE);
 378     Int_t fClsId = track->GetEMCALcluster();
 379     if(fClsId >= 0) isEMCAL = kTRUE;
 380     AliDebug(2, Form("cluster ID: %d, EMCAL: %s", fClsId, isEMCAL ? "yes" : "no"));
 381     if(isEMCAL) hfetrack.SetEMCALpid();
 382     // no filter bits available for ESDs
 383
 384     // V0 information
 385     AliHFEreducedTrack::EV0PID_t v0pid = AliHFEreducedTrack::kV0undef;
 386     if(track->TestBit(BIT(14))) v0pid = AliHFEreducedTrack::kV0electron;
 387     else if(track->TestBit(BIT(15))) v0pid = AliHFEreducedTrack::kV0pion;
 388     else if(track->TestBit(BIT(16))) v0pid = AliHFEreducedTrack::kV0proton;
 389     hfetrack.SetV0PID(v0pid);
 390
 391     if(mcthere){
 392       // Fill Monte-Carlo Information
 393       Int_t label = TMath::Abs(track->GetLabel());
 394       if(label < fMCEvent->GetNumberOfTracks())
 395         mctrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(label));
 396       if(mctrack){
 397         AliDebug(2, "Associated MC particle found");
 398         if(fTrackCuts->CheckParticleCuts(static_cast<UInt_t>(AliHFEcuts::kStepMCGenerated), mctrack)) hfetrack.SetMCSignal();
 399         // Kinematics
 400         hfetrack.SetMCSignedPt(mctrack->Pt(),mctrack->Charge() > 0.);
 401         hfetrack.SetMCP(mctrack->P());
 402         hfetrack.SetMCEta(mctrack->Eta());
 403         hfetrack.SetMCPhi(mctrack->Phi());
 404         hfetrack.SetMCPDG(mctrack->PdgCode());
 405
 406         // Get Production Vertex in radial direction
 407         hfetrack.SetMCProdVtx(mctrack->Xv(),mctrack->Yv(),mctrack->Zv());
 408
 409         // Get Mother PDG code of the particle
 410         Int_t motherlabel = TMath::Abs(mctrack->GetMother());
 411         if(motherlabel >= 0 && motherlabel < fMCEvent->GetNumberOfTracks()){
 412           AliMCParticle *mother = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(motherlabel));
 413           if(mother) hfetrack.SetMCMotherPdg(mother->PdgCode());
 414         }
 415
 416         // derive source
 417         source = 5;
 418         if(fSignalCuts->IsCharmElectron(mctrack)) source = 0;
 419         else if(fSignalCuts->IsBeautyElectron(mctrack)) source = 1;
 420         else if(fSignalCuts->IsGammaElectron(mctrack)) source = 2;
 421         else if(fSignalCuts->IsNonHFElectron(mctrack)) source = 3;
 422         else if(TMath::Abs(mctrack->PdgCode()) == 11) source = 4;
 423         else source = 5;
 424         hfetrack.SetMCSource(source);
 425       } else {
 426         AliDebug(2, "Associated MC particle not found");
 427       }
 428     }
 429
 430     // HFE DCA
 431     Float_t dcaxy = -999.,
 432             dcaz = -999.;
 433     fExtraCuts->GetImpactParameters((AliVTrack *)track,dcaxy,dcaz);
 434     hfetrack.SetDCA(dcaxy, dcaz);
 435     Double_t hfeImpactParam(-999.), hfeImpactParamResol(-999.);
 436     fExtraCuts->GetHFEImpactParameters((AliVTrack *)track,hfeImpactParam,hfeImpactParamResol);
 437     hfetrack.SetHFEImpactParam(hfeImpactParam,hfeImpactParamResol);
 438
 439     // Different number of clusters definitions
 440     Int_t nclustersITS(track->GetITSclusters(NULL)),
 441           nclustersTPC(track->GetTPCNcls()),
 442           nclustersTPCall(track->GetTPCClusterMap().CountBits()),
 443           nclustersTPCshared(0);
 444     UChar_t nfindableTPC = track->GetTPCNclsF();
 445     const TBits &sharedTPC = track->GetTPCSharedMap();
 446     for(Int_t ibit = 0; ibit < 160; ibit++) if(sharedTPC.TestBitNumber(ibit)) nclustersTPCshared++;
 447     hfetrack.SetChi2PerTPCcluster(track->GetTPCchi2()/Double_t(nclustersTPC));
 448     hfetrack.SetITSnclusters(nclustersITS);
 449     hfetrack.SetTPCnclusters(nclustersTPC);
 450     hfetrack.SetTRDnclusters(track->GetTRDncls());
 451     hfetrack.SetTPCnclustersPID(track->GetTPCsignalN());
 452     hfetrack.SetTPCcrossedRows(track->GetTPCCrossedRows());
 453     hfetrack.SetTPCnclustersAll(nclustersTPCall);
 454     hfetrack.SetTPCsharedClusters(nclustersTPCshared);
 455     hfetrack.SetTPCclusterRatio(nfindableTPC ? static_cast<Float_t>(nclustersTPC)/static_cast<Float_t>(nfindableTPC) : 0);
 456     hfetrack.SetTPCclusterRatioAll(nfindableTPC ? static_cast<Float_t>(nclustersTPCall)/static_cast<Float_t>(nfindableTPC) : 0);
 457     UChar_t itsPixel = track->GetITSClusterMap();
 458     for(int ily = 0; ily < 6; ily++)
 459             if(TESTBIT(itsPixel, ily)) hfetrack.SetITScluster(ily);
 460
 461     // TRD related quantities (Yvonne)
 462     Int_t nslices = track->GetNumberOfTRDslices();
 463     hfetrack.SetTRDntrackletsPID(track->GetTRDntrackletsPID());
 464     hfetrack.SetTRDnslices(nslices);
 465     hfetrack.SetTRDchi2(track->GetTRDchi2());
 466     Int_t nslicetemp=0;
 467     for(Int_t iplane = 0; iplane < 6; iplane++){
 468             nslicetemp=0;
 469             for(Int_t isl = 0; isl < nslices; isl++){
 470               if(track->GetTRDntrackletsPID()>0){
 471                       if(track->GetTRDslice(iplane, isl)>0.001) nslicetemp++;
 472               }
 473             }
 474             if(nslicetemp > 0) hfetrack.SetTRDstatus(iplane);
 475     }
 476
 477
 478     //test for kink tracks
 479     if(fExtraCuts->IsKinkMother(track)) hfetrack.SetIsKinkMother();
 480     else if(fExtraCuts->IsKinkDaughter(track)) hfetrack.SetIsKinkDaughter();
 481
 482     // Double counted
 483     Int_t id(track->GetID());
 484     for(Int_t l=0; l < counterdc; l++){
 485       Int_t iTrack2 = arraytrack.At(l);
 486       if(iTrack2==id){
 487          hfetrack.SetDoubleCounted();
 488          break;
 489       }
 490     }
 491     // Add the id at this place
 492     arraytrack.AddAt(id,counterdc);
 493     counterdc++;
 494
 495     // PID
 496     hfetrack.SetTPCdEdx(track->GetTPCsignal());
 497     hfetrack.SetTPCsigmaEl(pid->NumberOfSigmasTPC(track, AliPID::kElectron));
 498     hfetrack.SetTOFsigmaEl(pid->NumberOfSigmasTOF(track, AliPID::kElectron));
 499     hfetrack.SetTOFmismatchProbability(pid->GetTOFMismatchProbability(track));
 500     hfetrack.SetITSsigmaEl(pid->NumberOfSigmasITS(track, AliPID::kElectron));
 501     // Eta correction
 502     copyTrack.~AliESDtrack();
 503     new(&copyTrack) AliESDtrack(*track);
 504     if(fTPCpid->HasCentralityCorrection()) fTPCpid->ApplyCentralityCorrection(&copyTrack, static_cast<Double_t>(ncontrib),AliHFEpidObject::kESDanalysis);
 505     if(fTPCpid->HasEtaCorrection()) fTPCpid->ApplyEtaCorrection(&copyTrack, AliHFEpidObject::kESDanalysis);
 506     hfetrack.SetTPCsigmaElCorrected(pid->NumberOfSigmasTPC(&copyTrack, AliPID::kElectron));
 507     hfetrack.SetTPCdEdxCorrected(copyTrack.GetTPCsignal());
 508     if(isEMCAL){
 509       AliDebug(2, "Adding EMCAL PID information");
 510       // EMCAL cluster
 511       Double_t emcalEnergyOverP = -1.,
 512                showershape[4] = {0.,0.,0.,0.};
 513       hfetrack.SetEMCALSigmaEl(pid->NumberOfSigmasEMCAL(track, AliPID::kElectron, emcalEnergyOverP, &showershape[0]));
 514       hfetrack.SetEMCALEoverP(emcalEnergyOverP);
 515       hfetrack.SetEMCALShowerShape(showershape);
 516     }
 517
 518     // Track finished, add NOW to the Event
 519     fHFEevent->AddTrack(&hfetrack);
 520     //printf("after\n");
 521   }
 522
 523   // Fill the debug tree
 524   //AliInfo(Form("Number of tracks: %d\n", fHFEevent->GetNumberOfTracks()));
 525   //AliInfo(Form("Number of MC particles: %d\n", fHFEevent->GetNumberOfMCParticles()));
 526   fHFEtree->Fill();
 527
 528   fEventNumber++;
 529   PostData(1, fHFEtree);
 530 }
 531
 532 void AliHFEreducedEventCreatorESD::Terminate(Option_t *){
 533   //
 534   // Terminate
 535   //
 536   AliInfo("terminating...\n");
 537
 538 }
 539
 540 Bool_t AliHFEreducedEventCreatorESD::IsTOFmismatch(const AliVTrack *const track, const AliPIDResponse *const pid) const {
 541   //
 542   // Is TOF mismatch
 543   //
 544   Double_t probs[AliPID::kSPECIESC];
 545   AliPIDResponse::EDetPidStatus status = pid->ComputeTOFProbability(track, AliPID::kSPECIESC, probs);
 546   return status == AliPIDResponse::kDetMismatch;
 547 }
 548