PWGHF/hfe/AliHFENonPhotonicElectron.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 /*************************************************************************************
  18  *                                                                                      *
  19  *      Class for the Selection of Non-Heavy-Flavour-Electrons trought          *
  20  *      the invariant mass method. The selection can be done from two                   *
  21  *      different algorithms, which can be choosed calling the function         *
  22  *              "SetAlgorithm(TString Algorithm)".                                      *
  23  *                                                                                      *
  24  *              Authors: R.Bailhache, C.A.Schmidt                                       *
  25  *                                                                                      *
  26  *************************************************************************************/
  27
  28 #include "TVector2.h"
  29 #include "THnSparse.h"
  30 #include "TMath.h"
  31 #include "TLorentzVector.h"
  32 #include "TParticle.h"
  33 #include "TList.h"
  34 #include "TDatabasePDG.h"
  35
  36 #include "AliVEvent.h"
  37 #include "AliMCEvent.h"
  38 #include "AliESDEvent.h"
  39 #include "AliMCParticle.h"
  40 #include "AliAODMCParticle.h"
  41 #include "AliAODEvent.h"
  42 #include "AliAODVertex.h"
  43 #include "AliAODTrack.h"
  44 #include "AliVTrack.h"
  45 #include "AliESDtrack.h"
  46 #include "AliESDtrackCuts.h"
  47 #include "AliPIDResponse.h"
  48
  49 #include "AliKFParticle.h"
  50 #include "AliKFVertex.h"
  51
  52 #include "AliHFEcuts.h"
  53 #include "AliHFEpid.h"
  54 #include "AliHFEpidQAmanager.h"
  55 #include "AliHFEtools.h"
  56
  57 #include "AliHFENonPhotonicElectron.h"
  58
  59 ClassImp(AliHFENonPhotonicElectron)
  60 //________________________________________________________________________
  61 AliHFENonPhotonicElectron::AliHFENonPhotonicElectron(const char *name, const Char_t *title)
  62   :TNamed(name, title)
  63   ,fMCEvent(NULL)
  64   ,fAODArrayMCInfo(NULL)
  65   ,fHFEBackgroundCuts(0)
  66   ,fPIDBackground(0x0)
  67   ,fPIDBackgroundQA(0)
  68   ,fAlgorithmMA(kTRUE)
  69   ,fUseFilterAOD(kTRUE)
  70   ,fFilter(-1)
  71   ,fChi2OverNDFCut(3.0)
  72   ,fMaxDCA(3.0)
  73   ,fMaxOpeningTheta(0.02)
  74   ,fMaxOpeningPhi(0.1)
  75   ,fMaxOpening3D(TMath::Pi())
  76   ,fMaxInvMass(0.6)
  77   ,fSetMassConstraint(kFALSE)
  78   ,fArraytrack(NULL)
  79   ,fCounterPoolBackground(0)
  80   ,fListOutput(NULL)
  81   ,fMCSourceee(NULL)
  82   ,fUSignee(NULL)
  83   ,fLSignee(NULL)
  84   ,fUSignAngle(NULL)
  85   ,fLSignAngle(NULL)
  86 {
  87   //
  88   // Constructor
  89   //
  90   fPIDBackground = new AliHFEpid("hfePidBackground");
  91   fPIDBackgroundQA = new AliHFEpidQAmanager;
  92 }
  93
  94 //________________________________________________________________________
  95 AliHFENonPhotonicElectron::AliHFENonPhotonicElectron()
  96   :TNamed()
  97   ,fMCEvent(0)
  98   ,fAODArrayMCInfo(NULL)
  99   ,fHFEBackgroundCuts(NULL)
 100   ,fPIDBackground(0x0)
 101   ,fPIDBackgroundQA(0)
 102   ,fAlgorithmMA(kTRUE)
 103   ,fUseFilterAOD(kTRUE)
 104   ,fFilter(-1)
 105   ,fChi2OverNDFCut(3.0)
 106   ,fMaxDCA(3.0)
 107   ,fMaxOpeningTheta(0.02)
 108   ,fMaxOpeningPhi(0.1)
 109   ,fMaxOpening3D(TMath::TwoPi())
 110   ,fMaxInvMass(0.6)
 111   ,fSetMassConstraint(kFALSE)
 112   ,fArraytrack(NULL)
 113   ,fCounterPoolBackground(0)
 114   ,fListOutput(NULL)
 115   ,fMCSourceee(NULL)
 116   ,fUSignee(NULL)
 117   ,fLSignee(NULL)
 118   ,fUSignAngle(NULL)
 119   ,fLSignAngle(NULL)
 120 {
 121   //
 122   // Constructor
 123   //
 124   fPIDBackground   = new AliHFEpid("hfePidBackground");
 125   fPIDBackgroundQA = new AliHFEpidQAmanager;
 126 }
 127
 128 //________________________________________________________________________
 129 AliHFENonPhotonicElectron::AliHFENonPhotonicElectron(const AliHFENonPhotonicElectron &ref)
 130   :TNamed(ref)
 131   ,fMCEvent(NULL)
 132   ,fAODArrayMCInfo(NULL)
 133   ,fHFEBackgroundCuts(ref.fHFEBackgroundCuts)
 134   ,fPIDBackground(ref.fPIDBackground)
 135   ,fPIDBackgroundQA(ref.fPIDBackgroundQA)
 136   ,fAlgorithmMA(ref.fAlgorithmMA)
 137   ,fUseFilterAOD(ref.fUseFilterAOD)
 138   ,fFilter(ref.fFilter)
 139   ,fChi2OverNDFCut(ref.fChi2OverNDFCut)
 140   ,fMaxDCA(ref.fMaxDCA)
 141   ,fMaxOpeningTheta(ref.fMaxOpeningTheta)
 142   ,fMaxOpeningPhi(ref.fMaxOpeningPhi)
 143   ,fMaxOpening3D(ref.fMaxOpening3D)
 144   ,fMaxInvMass(ref.fMaxInvMass)
 145   ,fSetMassConstraint(ref.fSetMassConstraint)
 146   ,fArraytrack(NULL)
 147   ,fCounterPoolBackground(0)
 148   ,fListOutput(ref.fListOutput)
 149   ,fMCSourceee(ref.fMCSourceee)
 150   ,fUSignee(ref.fUSignee)
 151   ,fLSignee(ref.fLSignee)
 152   ,fUSignAngle(ref.fUSignAngle)
 153   ,fLSignAngle(ref.fLSignAngle)
 154 {
 155   //
 156   // Copy Constructor
 157   //
 158   ref.Copy(*this);
 159 }
 160
 161 //____________________________________________________________
 162 AliHFENonPhotonicElectron &AliHFENonPhotonicElectron::operator=(const AliHFENonPhotonicElectron &ref){
 163   //
 164   // Assignment operator
 165   //
 166   if(this == &ref) ref.Copy(*this);
 167   return *this;
 168 }
 169
 170 //_________________________________________
 171 AliHFENonPhotonicElectron::~AliHFENonPhotonicElectron()
 172 {
 173   //
 174   // Destructor
 175   //
 176   if(fArraytrack) delete fArraytrack;
 177   if(fHFEBackgroundCuts) delete fHFEBackgroundCuts;
 178   if(fPIDBackground) delete fPIDBackground;
 179 }
 180
 181 //_____________________________________________________________________________________________
 182 void AliHFENonPhotonicElectron::Init()
 183 {
 184   //
 185   // Init
 186   //
 187
 188   if(!fListOutput) fListOutput = new TList;
 189   fListOutput->SetName("HFENonPhotonicElectron");
 190   fListOutput->SetOwner();
 191
 192   if(!fHFEBackgroundCuts) fHFEBackgroundCuts = new AliESDtrackCuts();
 193
 194   // Initialize PID
 195   if(!fPIDBackground) fPIDBackground = new AliHFEpid("default pid");
 196   if(fMCEvent || fAODArrayMCInfo) fPIDBackground->SetHasMCData(kTRUE);
 197   if(!fPIDBackground->GetNumberOfPIDdetectors()) fPIDBackground->AddDetector("TPC", 0);
 198   AliInfo("PID Background QA switched on");
 199   fPIDBackgroundQA->Initialize(fPIDBackground);
 200   fListOutput->Add(fPIDBackgroundQA->MakeList("HFENP_PID_Background"));
 201   fPIDBackground->SortDetectors();
 202
 203   Int_t nBinsPt = 24;
 204   Double_t binLimPt[25] = {0., 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1., 1.1, 1.2, 1.3, 1.4, 1.5, 1.75, 2., 2.25, 2.5, 3., 3.5, 4., 5., 6.};
 205
 206   Int_t nBinsC = 11;
 207   Double_t minC = 0.0;
 208   Double_t maxC = 11.0;
 209   Double_t binLimC[nBinsC+1];
 210   for(Int_t i=0; i<=nBinsC; i++) binLimC[i]=(Double_t)minC + (maxC-minC)/nBinsC*(Double_t)i ;
 211
 212   Int_t nBinsSource = 10;
 213   Double_t minSource = 0.;
 214   Double_t maxSource = 10.;
 215   Double_t binLimSource[nBinsSource+1];
 216   for(Int_t i=0; i<=nBinsSource; i++) binLimSource[i]=(Double_t)minSource + (maxSource-minSource)/nBinsSource*(Double_t)i ;
 217
 218   Int_t nBinsInvMass = 60;
 219   Double_t minInvMass = 0.;
 220   Double_t maxInvMass = 0.6;
 221   Double_t binLimInvMass[nBinsInvMass+1];
 222   for(Int_t i=0; i<=nBinsInvMass; i++) binLimInvMass[i]=(Double_t)minInvMass + (maxInvMass-minInvMass)/nBinsInvMass*(Double_t)i ;
 223
 224   Int_t nBinsPhi = 180;
 225   Double_t minPhi = 0.0;
 226   Double_t maxPhi = TMath::Pi();
 227   Double_t binLimPhi[nBinsPhi+1];
 228   for(Int_t i=0; i<=nBinsPhi; i++)
 229   {
 230     binLimPhi[i]=(Double_t)minPhi + (maxPhi-minPhi)/nBinsPhi*(Double_t)i ;
 231     AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
 232   }
 233
 234   Int_t nBinsAngle = 180;
 235   Double_t minAngle = 0.0;
 236   Double_t maxAngle = TMath::Pi();
 237   Double_t binLimAngle[nBinsAngle+1];
 238   for(Int_t i=0; i<=nBinsAngle; i++)
 239   {
 240     binLimAngle[i]=(Double_t)minAngle + (maxAngle-minAngle)/nBinsAngle*(Double_t)i ;
 241     AliDebug(2,Form("bin phi is %f for %d",binLimPhi[i],i));
 242   }
 243
 244   const Int_t nDimMCSource=3;
 245   Int_t nBinMCSource[nDimMCSource] = {nBinsC,nBinsPt,nBinsSource};
 246   fMCSourceee = new THnSparseF("fMCSourceee","fMCSourceee",nDimMCSource,nBinMCSource);
 247   fMCSourceee->SetBinEdges(0,binLimC);
 248   fMCSourceee->SetBinEdges(1,binLimPt);
 249   fMCSourceee->SetBinEdges(2,binLimSource);
 250 /*  fMCSourceee->GetAxis(0)->SetTitle("Multiplicity (arb.unit)");
 251   fMCSourceee->GetAxis(1)->SetTitle("#it{p}_{mc,T} (GeV/c)");
 252   fMCSourceee->GetAxis(2)->SetTitle("MC");*/
 253   fMCSourceee->Sumw2();
 254   AliDebug(2,"AliHFENonPhotonicElectron: fMCSourceee");
 255
 256   // ee invariant mass Unlike Sign
 257   const Int_t nDimUSign=6;
 258   Int_t nBinUSign[nDimUSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource,nBinsAngle};
 259   fUSignee = new THnSparseF("fUSignee","fUSignee",nDimUSign,nBinUSign);
 260   fUSignee->SetBinEdges(0,binLimPhi);
 261   fUSignee->SetBinEdges(1,binLimC);
 262   fUSignee->SetBinEdges(2,binLimPt);
 263   fUSignee->SetBinEdges(3,binLimInvMass);
 264   fUSignee->SetBinEdges(4,binLimSource);
 265   fUSignee->SetBinEdges(5,binLimAngle);
 266 /*  fUSignee->GetAxis(0)->SetTitle("#Delta #phi");
 267   fUSignee->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
 268   fUSignee->GetAxis(2)->SetTitle("#it{p}_{T} (GeV/c)");
 269   fUSignee->GetAxis(3)->SetTitle("m_{e^{#pm}e^{#mp}}} (GeV/#it{c}^{2}");
 270   fUSignee->GetAxis(4)->SetTitle("MC");
 271   fUSignee->GetAxis(5)->SetTitle("#alpha");*/
 272   fUSignee->Sumw2();
 273   AliDebug(2,"AliHFENonPhotonicElectron: fUSignee");
 274
 275   // ee invariant mass Like Sign
 276   const Int_t nDimLSign=6;
 277   Int_t nBinLSign[nDimLSign] = {nBinsPhi,nBinsC,nBinsPt,nBinsInvMass,nBinsSource,nBinsAngle};
 278   fLSignee = new THnSparseF("fLSignee","fLSignee",nDimLSign,nBinLSign);
 279   fLSignee->SetBinEdges(0,binLimPhi);
 280   fLSignee->SetBinEdges(1,binLimC);
 281   fLSignee->SetBinEdges(2,binLimPt);
 282   fLSignee->SetBinEdges(3,binLimInvMass);
 283   fLSignee->SetBinEdges(4,binLimSource);
 284   fLSignee->SetBinEdges(5,binLimAngle);
 285 /*  fLSignee->GetAxis(0)->SetTitle("#Delta #phi");
 286   fLSignee->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
 287   fLSignee->GetAxis(2)->SetTitle("#it{p}_{T} (GeV/c)");
 288   fLSignee->GetAxis(3)->SetTitle("m_{e^{#pm}e^{#pm}}} (GeV/#it{c}^{2}");
 289   fLSignee->GetAxis(4)->SetTitle("MC");
 290   fLSignee->GetAxis(5)->SetTitle("#alpha");*/
 291   fLSignee->Sumw2();
 292   AliDebug(2,"AliHFENonPhotonicElectron: fLSignee");
 293
 294   // ee angle Unlike Sign
 295   const Int_t nDimUSignAngle=3;
 296   Int_t nBinUSignAngle[nDimUSignAngle] = {nBinsAngle,nBinsC,nBinsSource};
 297   fUSignAngle = new THnSparseF("fUSignAngle","fUSignAngle",nDimUSignAngle,nBinUSignAngle);
 298   fUSignAngle->SetBinEdges(0,binLimAngle);
 299   fUSignAngle->SetBinEdges(1,binLimC);
 300   fUSignAngle->SetBinEdges(2,binLimSource);
 301 /*  fUSignAngle->GetAxis(0)->SetTitle("#alpha");
 302   fUSignAngle->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
 303   fUSignAngle->GetAxis(2)->SetTitle("MC");*/
 304   fUSignAngle->Sumw2();
 305   AliDebug(2,"AliHFENonPhotonicElectron: fUSignAngle");
 306
 307   // ee angle Like Sign
 308   const Int_t nDimLSignAngle=3;
 309   Int_t nBinLSignAngle[nDimLSignAngle] = {nBinsAngle,nBinsC,nBinsSource};
 310   fLSignAngle = new THnSparseF("fLSignAngle","fLSignAngle",nDimLSignAngle,nBinLSignAngle);
 311   fLSignAngle->SetBinEdges(0,binLimAngle);
 312   fLSignAngle->SetBinEdges(1,binLimC);
 313   fLSignAngle->SetBinEdges(2,binLimSource);
 314 /*  fLSignAngle->GetAxis(0)->SetTitle("#alpha");
 315   fLSignAngle->GetAxis(1)->SetTitle("Multiplicity (arb.unit)");
 316   fLSignAngle->GetAxis(2)->SetTitle("MC");*/
 317   fLSignAngle->Sumw2();
 318   AliDebug(2,"AliHFENonPhotonicElectron: fLSignAngle");
 319
 320
 321   fListOutput->Add(fMCSourceee);
 322   fListOutput->Add(fUSignee);
 323   fListOutput->Add(fLSignee);
 324   fListOutput->Add(fUSignAngle);
 325   fListOutput->Add(fLSignAngle);
 326
 327 }
 328
 329 //_____________________________________________________________________________________________
 330 void AliHFENonPhotonicElectron::InitRun(const AliVEvent *inputEvent,const AliPIDResponse *pidResponse)
 331 {
 332   //
 333   // Init run
 334   //
 335
 336   if(!pidResponse)
 337   {
 338     AliDebug(1, "Using default PID Response");
 339     Bool_t hasmc = kFALSE;
 340     if(fMCEvent || fAODArrayMCInfo) hasmc=kTRUE;
 341     pidResponse = AliHFEtools::GetDefaultPID(hasmc, inputEvent->IsA() == AliESDEvent::Class());
 342   }
 343
 344   if(!fPIDBackground) return;
 345   fPIDBackground->SetPIDResponse(pidResponse);
 346
 347   if(!fPIDBackground->IsInitialized())
 348   {
 349     // Initialize PID with the given run number
 350     fPIDBackground->InitializePID(inputEvent->GetRunNumber());
 351   }
 352
 353 }
 354
 355 //_____________________________________________________________________________________________
 356 Int_t AliHFENonPhotonicElectron::FillPoolAssociatedTracks(AliVEvent *inputEvent,Int_t binct)
 357 {
 358   //
 359   // Fill the pool of associated tracks
 360   // Return the number of associated tracks
 361   //
 362
 363   Int_t nbtracks = inputEvent->GetNumberOfTracks();
 364
 365   if( fArraytrack )
 366   {
 367     fArraytrack->~TArrayI();
 368     new(fArraytrack) TArrayI(nbtracks);
 369   }
 370   else
 371   {
 372     fArraytrack = new TArrayI(nbtracks);
 373   }
 374
 375   fCounterPoolBackground = 0;
 376
 377   for(Int_t k = 0; k < nbtracks; k++)
 378   {
 379     AliVTrack *track = (AliVTrack *) inputEvent->GetTrack(k);
 380     if(!track) continue;
 381
 382     // Track cuts
 383     Bool_t survivedbackground = kTRUE;
 384     AliAODEvent *aodeventu = dynamic_cast<AliAODEvent *>(inputEvent);
 385
 386     if(aodeventu)
 387     {
 388       // AOD
 389       AliAODTrack *aodtrack = dynamic_cast<AliAODTrack *>(track);
 390       if(aodtrack)
 391       {
 392         // filter
 393         if(fUseFilterAOD)
 394         {
 395           if(!(aodtrack->TestFilterBit(fFilter))) survivedbackground = kFALSE;
 396         }
 397
 398         // additional cuts
 399         if(survivedbackground)
 400         {
 401           AliESDtrack esdTrack(aodtrack);
 402
 403           // set the TPC cluster info
 404           esdTrack.SetTPCClusterMap(aodtrack->GetTPCClusterMap());
 405           esdTrack.SetTPCSharedMap(aodtrack->GetTPCSharedMap());
 406           esdTrack.SetTPCPointsF(aodtrack->GetTPCNclsF());
 407           AliAODVertex *vAOD = aodeventu->GetPrimaryVertex();
 408           Double_t bfield = aodeventu->GetMagneticField();
 409           Double_t pos[3],cov[6];
 410           vAOD->GetXYZ(pos);
 411           vAOD->GetCovarianceMatrix(cov);
 412           const AliESDVertex vESD(pos,cov,100.,100);
 413           esdTrack.RelateToVertex(&vESD,bfield,3.);
 414
 415           if(!fHFEBackgroundCuts->IsSelected(&esdTrack))
 416           {
 417             survivedbackground = kFALSE;
 418           }
 419         }
 420       }
 421     }
 422     else
 423     {
 424       // ESD
 425       AliESDtrack *esdtrack = dynamic_cast<AliESDtrack *>(track);
 426       if(esdtrack)
 427       {
 428         if(!fHFEBackgroundCuts->IsSelected(esdtrack)) survivedbackground = kFALSE;
 429       }
 430     }
 431
 432     // PID
 433     if(survivedbackground)
 434     {
 435       // PID track cuts
 436       AliHFEpidObject hfetrack2;
 437
 438       if(!aodeventu) hfetrack2.SetAnalysisType(AliHFEpidObject::kESDanalysis);
 439       else hfetrack2.SetAnalysisType(AliHFEpidObject::kAODanalysis);
 440
 441       hfetrack2.SetRecTrack(track);
 442       if(binct>-1)
 443       {
 444         hfetrack2.SetCentrality((Int_t)binct);
 445         AliDebug(2,Form("centrality %d and %d",binct,hfetrack2.GetCentrality()));
 446         hfetrack2.SetPbPb();
 447       }
 448
 449       if(fPIDBackground->IsSelected(&hfetrack2,0x0,"recTrackCont",fPIDBackgroundQA))
 450       {
 451         fArraytrack->AddAt(k,fCounterPoolBackground);
 452         fCounterPoolBackground++;
 453         AliDebug(2,Form("fCounterPoolBackground %d, track %d",fCounterPoolBackground,k));
 454       }
 455     }
 456   } // loop tracks
 457
 458   //printf(Form("Test pool: nbtrack %d, binct %d \n fCounterPoolBackground %d \n",nbtracks,binct,fCounterPoolBackground));
 459
 460   return fCounterPoolBackground;
 461 }
 462
 463 //_____________________________________________________________________________________________
 464 Int_t AliHFENonPhotonicElectron::LookAtNonHFE(Int_t iTrack1, AliVTrack *track1, AliVEvent *vEvent, Double_t weight, Int_t binct, Double_t deltaphi, Int_t source, Int_t indexmother)
 465 {
 466   //
 467   // Look At Non HFE
 468   //
 469
 470   /***********************************************************************************
 471    *                                                                                    *
 472    *    iTrack1:        index of the tagged electrons in AliVEvent                      *
 473    *    track1:         tagged electron                                                 *
 474    *    vEvent:         event                                                           *
 475    *    weight:         weight in pt if not realistic                                   *
 476    *    binct:          centrality bin                                                  *
 477    *    deltaphi:       phi-phi event plane for v2                                      *
 478    *    source:         MC sources                                                      *
 479    *    indexmother:    MC index mother                                                 *
 480    *                                                                                    *
 481    *                                                                                    *
 482    *    return -1  if  nothing                                                          *
 483    *    return  2  if  opposite         charge          within the mass range           *
 484    *    return  4  if      like         charge          within the mass range           *
 485    *    return  6  if  opposite & like charge           within the mass range           *
 486    *                                                                                    *
 487    ***********************************************************************************/
 488
 489   AliAODEvent *aodeventu = dynamic_cast<AliAODEvent*>(vEvent);
 490   Int_t taggedphotonic = -1;
 491
 492   AliDebug(2,Form("fCounterPoolBackground %d in LookAtNonHFE!!!",fCounterPoolBackground));
 493   if(!fArraytrack) return taggedphotonic;
 494   AliDebug(2,Form("process track %d",iTrack1));
 495
 496   //Set Fill-Arrays for THnSparse
 497   Double_t value_MCSource[3] = { binct, track1->Pt(), source};          //Centrality    Pt              Source
 498   Double_t valueee[6] = { deltaphi, binct, track1->Pt(), -1, source,-1};        //DeltaPhi      Centrality      Pt      InvariantMass   Source
 499   Double_t valueAngle[3] = {    -1,     binct,          source};                //Angle         Centrality      Source
 500
 501   Bool_t USignPhotonic = kFALSE;
 502   Bool_t LSignPhotonic = kFALSE;
 503   Bool_t hasdcaT1 = kFALSE;
 504   Bool_t hasdcaT2 = kFALSE;
 505
 506   Int_t pdg1 = CheckPdg(TMath::Abs(track1->GetLabel()));
 507   Int_t pdg2 = -100;
 508   Int_t numberfound = 0;
 509   Int_t iTrack2 = 0;
 510   Int_t source2 = 0;
 511   Int_t indexmother2 = -1;
 512   Int_t fPDGtrack1 = 0;
 513   Int_t fPDGtrack2 = 0;
 514
 515   Double_t eMass = TDatabasePDG::Instance()->GetParticle(11)->Mass(); //Electron mass in GeV
 516   Double_t bfield = vEvent->GetMagneticField();
 517   Double_t xt1 = 0; //radial position track 1 at the DCA point
 518   Double_t xt2 = 0; //radial position track 2 at the DCA point
 519   Double_t p1[3] = {0,0,0};
 520   Double_t p2[3] = {0,0,0};
 521   Double_t dca12 = 0;
 522
 523   Double_t angleESD = 0;
 524   Double_t invmassESD = 0;
 525
 526   Double_t chi2OverNDF = 0;
 527   Double_t width = 0;
 528   Double_t angleAOD = 0;
 529   Double_t invmassAOD = 0;
 530
 531   Float_t fCharge1 = 0;
 532   Float_t fCharge2 = 0;
 533
 534   TLorentzVector electron1;
 535   TLorentzVector electron2;
 536   TLorentzVector mother;
 537
 538   AliVTrack *track2;
 539   AliESDtrack *esdtrack1;
 540   AliESDtrack *esdtrack2;
 541   AliKFParticle *ktrack1;
 542   AliKFParticle *ktrack2;
 543   AliKFVertex primV(*(vEvent->GetPrimaryVertex()));
 544
 545   // FILL MCsource TODO: correct?
 546   fMCSourceee->Fill(&value_MCSource[0],weight);
 547
 548
 549   for(Int_t idex = 0; idex < fCounterPoolBackground; idex++)
 550   {
 551     iTrack2 = fArraytrack->At(idex);
 552     AliDebug(2,Form("track %d",iTrack2));
 553     track2 = (AliVTrack *)vEvent->GetTrack(iTrack2);
 554
 555     if(!track2)
 556     {
 557       printf("ERROR: Could not receive track %d", iTrack2);
 558       continue;
 559     }
 560
 561     if(iTrack2==iTrack1) continue;
 562     AliDebug(2,"Different");
 563
 564     fCharge1 = track1->Charge();                //Charge from track1
 565     fCharge2 = track2->Charge();                //Charge from track2
 566
 567     // Reset the MC info
 568     valueAngle[2] = source;
 569     valueee[4] = source;
 570
 571     // track cuts and PID already done
 572
 573     // if MC look
 574     if(fMCEvent || fAODArrayMCInfo)
 575     {
 576       source2    = FindMother(TMath::Abs(track2->GetLabel()), indexmother2);
 577       pdg2       = CheckPdg(TMath::Abs(track2->GetLabel()));
 578
 579       if(source2 >=0 )
 580       {
 581         if((indexmother2 == indexmother) && (source == source2) && ((pdg1*pdg2)<0.0))
 582         {
 583           if(source == kElectronfromconversion)
 584           {
 585             valueAngle[2] = kElectronfromconversionboth;
 586             valueee[4] = kElectronfromconversionboth;
 587             numberfound++;
 588           }
 589
 590           if(source == kElectronfrompi0)
 591           {
 592             valueAngle[2] = kElectronfrompi0both;
 593             valueee[4] = kElectronfrompi0both;
 594           }
 595
 596           if(source == kElectronfrometa)
 597           {
 598             valueAngle[2] = kElectronfrometaboth;
 599             valueee[4] = kElectronfrometaboth;
 600           }
 601         }
 602       }
 603     }
 604
 605     if(fAlgorithmMA && (!aodeventu))
 606     {
 607       /**                               *
 608        *        ESD-Analysis            *
 609        **                               */
 610
 611       esdtrack1 = dynamic_cast<AliESDtrack *>(track1);                  //ESD-track1
 612       esdtrack2 = dynamic_cast<AliESDtrack *>(track2);                  //ESD-track2
 613       if((!esdtrack1) || (!esdtrack2)) continue;
 614
 615       dca12 = esdtrack2->GetDCA(esdtrack1,bfield,xt2,xt1);              //DCA track1-track2
 616
 617       if(dca12 > fMaxDCA) continue;                                     //! Cut on DCA
 618
 619       //Momento of the track extrapolated to DCA track-track
 620       hasdcaT1 = esdtrack1->GetPxPyPzAt(xt1,bfield,p1);         //Track1
 621       hasdcaT2 = esdtrack2->GetPxPyPzAt(xt2,bfield,p2);         //Track2
 622       if(!hasdcaT1 || !hasdcaT2) AliWarning("It could be a problem in the extrapolation");
 623
 624       electron1.SetXYZM(esdtrack1->Px(), esdtrack1->Py(), esdtrack1->Pz(), eMass);
 625       electron2.SetXYZM(esdtrack2->Px(), esdtrack2->Py(), esdtrack2->Pz(), eMass);
 626
 627       mother      = electron1 + electron2;
 628       invmassESD  = mother.M();
 629       angleESD    = TVector2::Phi_0_2pi(electron1.Angle(electron2.Vect()));
 630
 631       valueAngle[0] = angleESD;
 632       valueee[5] = angleESD;
 633
 634       if((fCharge1*fCharge2)>0.0)       fLSignAngle->Fill(&valueAngle[0],weight);
 635       else                              fUSignAngle->Fill(&valueAngle[0],weight);
 636
 637       if(angleESD > fMaxOpening3D) continue;                             //! Cut on Opening Angle
 638       valueee[3] = invmassESD;
 639
 640       if((fCharge1*fCharge2)>0.0)       fLSignee->Fill(&valueee[0],weight);
 641       else                              fUSignee->Fill(&valueee[0],weight);
 642
 643       if(invmassESD > fMaxInvMass) continue;                            //! Cut on Invariant Mass
 644
 645       if((fCharge1*fCharge2)>0.0)       LSignPhotonic=kTRUE;
 646       else                              USignPhotonic=kTRUE;
 647     }
 648     else
 649     {
 650       /**                               *
 651        *        AliKF-Analysis          *
 652        **                               */
 653
 654       fPDGtrack1 = 11;
 655       fPDGtrack2 = 11;
 656
 657       if(fCharge1>0) fPDGtrack1 = -11;
 658       if(fCharge2>0) fPDGtrack2 = -11;
 659
 660       ktrack1 = new AliKFParticle(*track1, fPDGtrack1);
 661       ktrack2 = new AliKFParticle(*track2, fPDGtrack2);
 662       AliKFParticle recoGamma(*ktrack1,*ktrack2);
 663
 664       if(recoGamma.GetNDF()<1) continue;                                //! Cut on Reconstruction
 665
 666       chi2OverNDF = recoGamma.GetChi2()/recoGamma.GetNDF();
 667       if(TMath::Sqrt(TMath::Abs(chi2OverNDF))>fChi2OverNDFCut) continue;
 668
 669       // DCA
 670       //Double_t dca12 = ktrack1.GetDistanceFromParticle(ktrack2);
 671       //if(dca12 > fMaxDCA) continue;
 672
 673       // if set mass constraint
 674       if(fSetMassConstraint) //&& pVtx)
 675       {
 676         primV += recoGamma;
 677         primV -= *ktrack1;
 678         primV -= *ktrack2;
 679         recoGamma.SetProductionVertex(primV);
 680         recoGamma.SetMassConstraint(0,0.0001);
 681       }
 682
 683       recoGamma.GetMass(invmassAOD,width);
 684       angleAOD = ktrack1->GetAngle(*ktrack2);
 685
 686       valueAngle[0] = angleAOD;
 687       valueee[5] = angleAOD;
 688
 689       if((fCharge1*fCharge2)>0.0)       fLSignAngle->Fill(&valueAngle[0],weight);
 690       else                              fUSignAngle->Fill(&valueAngle[0],weight);
 691
 692       if(angleAOD > fMaxOpening3D) continue;                            //! Cut on Opening Angle
 693
 694       valueee[3] = invmassAOD;
 695
 696       if((fCharge1*fCharge2)>0.0)       fLSignee->Fill(&valueee[0],weight);
 697       else                              fUSignee->Fill(&valueee[0],weight);
 698
 699       if(invmassAOD > fMaxInvMass) continue;                            //! Cut on Invariant Mass
 700
 701       if((fCharge1*fCharge2)>0.0)       LSignPhotonic=kTRUE;
 702       else                              USignPhotonic=kTRUE;
 703     }
 704   }
 705
 706   if( USignPhotonic &&  LSignPhotonic) taggedphotonic = 6;
 707   if(!USignPhotonic &&  LSignPhotonic) taggedphotonic = 4;
 708   if( USignPhotonic && !LSignPhotonic) taggedphotonic = 2;
 709
 710   return taggedphotonic;
 711 }
 712
 713 //_________________________________________________________________________
 714 Int_t AliHFENonPhotonicElectron::FindMother(Int_t tr, Int_t &indexmother){
 715   //
 716   // Find the mother if MC
 717   //
 718
 719   if(!fMCEvent && !fAODArrayMCInfo) return 0;
 720
 721   Int_t pdg = CheckPdg(tr);
 722   if(TMath::Abs(pdg)!= 11)
 723   {
 724     indexmother = -1;
 725     return kNoElectron;
 726   }
 727
 728   indexmother = IsMotherGamma(tr);
 729   if(indexmother > 0) return kElectronfromconversion;
 730   indexmother = IsMotherPi0(tr);
 731   if(indexmother > 0) return kElectronfrompi0;
 732   indexmother = IsMotherC(tr);
 733   if(indexmother > 0) return kElectronfromC;
 734   indexmother = IsMotherB(tr);
 735   if(indexmother > 0) return kElectronfromB;
 736   indexmother = IsMotherEta(tr);
 737   if(indexmother > 0) return kElectronfrometa;
 738
 739   return kElectronfromother;
 740 }
 741
 742 //________________________________________________________________________________________________
 743 Int_t AliHFENonPhotonicElectron::CheckPdg(Int_t tr) {
 744
 745   //
 746   // Return the pdg of the particle
 747   //
 748
 749   Int_t pdgcode = -1;
 750   if(tr < 0) return pdgcode;
 751
 752   if(fMCEvent)
 753   {
 754     AliVParticle *mctrack = fMCEvent->GetTrack(tr);
 755     if(!mctrack) return -1;
 756     AliMCParticle *mctrackesd = NULL;
 757     if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return pdgcode;
 758     pdgcode = mctrackesd->PdgCode();
 759   }
 760
 761   if(fAODArrayMCInfo)
 762   {
 763     if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
 764     AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
 765     if(!mctrackaod) return pdgcode;
 766     pdgcode = mctrackaod->GetPdgCode();
 767   }
 768
 769   return pdgcode;
 770 }
 771
 772 //_______________________________________________________________________________________________
 773 Int_t AliHFENonPhotonicElectron::IsMotherGamma(Int_t tr) {
 774
 775   //
 776   // Return the lab of gamma mother or -1 if not gamma
 777   //
 778
 779   if(tr < 0) return -1;
 780
 781   if(fMCEvent)
 782   {
 783     AliVParticle *mctrack = fMCEvent->GetTrack(tr);
 784     if(!mctrack) return -1;
 785     AliMCParticle *mctrackesd = NULL;
 786     if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
 787     TParticle *particle = 0x0;
 788     particle = mctrackesd->Particle();
 789
 790     // Take mother
 791     if(!particle) return -1;
 792     Int_t imother   = particle->GetFirstMother();
 793     if(imother < 0) return -1;
 794     AliMCParticle *mothertrack = NULL;
 795     if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
 796     TParticle * mother = mothertrack->Particle();
 797     if(!mother) return -1;
 798
 799     // Check gamma
 800     Int_t pdg = mother->GetPdgCode();
 801     if(TMath::Abs(pdg) == 22) return imother;
 802     if(TMath::Abs(pdg) == 11)
 803     {
 804       return IsMotherGamma(imother);
 805     }
 806
 807     return -1;
 808   }
 809
 810   if(fAODArrayMCInfo)
 811   {
 812     if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
 813     AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
 814     if(!mctrackaod) return -1;
 815
 816     // Take mother
 817     Int_t imother = mctrackaod->GetMother();
 818     if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
 819     AliAODMCParticle *mothertrack = NULL;
 820     if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
 821
 822     // Check gamma
 823     Int_t pdg = mothertrack->GetPdgCode();
 824     if(TMath::Abs(pdg) == 22) return imother;
 825     if(TMath::Abs(pdg) == 11)
 826     {
 827       return IsMotherGamma(imother);
 828     }
 829
 830     return -1;
 831   }
 832
 833   return -1;
 834 }
 835
 836 //________________________________________________________________________________________________
 837 Int_t AliHFENonPhotonicElectron::IsMotherPi0(Int_t tr) {
 838
 839   //
 840   // Return the lab of pi0 mother or -1 if not pi0
 841   //
 842
 843   if(tr < 0) return -1;
 844
 845   if(fMCEvent)
 846   {
 847     AliVParticle *mctrack = fMCEvent->GetTrack(tr);
 848     if(!mctrack) return -1;
 849     AliMCParticle *mctrackesd = NULL;
 850     if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
 851     TParticle *particle = 0x0;
 852     particle = mctrackesd->Particle();
 853     // Take mother
 854     if(!particle) return -1;
 855     Int_t imother   = particle->GetFirstMother();
 856     if(imother < 0) return -1;
 857     AliMCParticle *mothertrack = NULL;
 858     if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
 859     TParticle * mother = mothertrack->Particle();
 860     if(!mother) return -1;
 861     // Check pi0
 862     Int_t pdg = mother->GetPdgCode();
 863     if(TMath::Abs(pdg) == 111) return imother;
 864     if(TMath::Abs(pdg) == 11)
 865     {
 866       return IsMotherPi0(imother);
 867     }
 868
 869     return -1;
 870   }
 871
 872   if(fAODArrayMCInfo)  {
 873
 874     if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
 875     AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
 876     if(!mctrackaod) return -1;
 877
 878     // Take mother
 879     Int_t imother = mctrackaod->GetMother();
 880     if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
 881     AliAODMCParticle *mothertrack = NULL;
 882     if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
 883     // Check pi0
 884     Int_t pdg = mothertrack->GetPdgCode();
 885     if(TMath::Abs(pdg) == 111) return imother;
 886     if(TMath::Abs(pdg) == 11)
 887     {
 888       return IsMotherPi0(imother);
 889     }
 890
 891     return -1;
 892   }
 893
 894   return -1;
 895 }
 896 //________________________________________________________________________________________________
 897 Int_t AliHFENonPhotonicElectron::IsMotherC(Int_t tr) {
 898
 899   //
 900   // Return the lab of signal mother or -1 if not from C
 901   //
 902
 903   if(tr < 0) return -1;
 904
 905   if(fMCEvent)
 906   {
 907     AliVParticle *mctrack = fMCEvent->GetTrack(tr);
 908     if(!mctrack) return -1;
 909     AliMCParticle *mctrackesd = NULL;
 910     if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
 911     TParticle *particle = 0x0;
 912     particle = mctrackesd->Particle();
 913
 914     // Take mother
 915     if(!particle) return -1;
 916     Int_t imother   = particle->GetFirstMother();
 917     if(imother < 0) return -1;
 918     AliMCParticle *mothertrack = NULL;
 919     if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
 920     TParticle * mother = mothertrack->Particle();
 921     if(!mother) return -1;
 922
 923     // Check C
 924     Int_t pdg = mother->GetPdgCode();
 925     if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;
 926     if(TMath::Abs(pdg) == 11)
 927     {
 928       return IsMotherC(imother);
 929     }
 930
 931     return -1;
 932   }
 933
 934   if(fAODArrayMCInfo)
 935   {
 936     if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
 937     AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
 938     if(!mctrackaod) return -1;
 939
 940     // Take mother
 941     Int_t imother = mctrackaod->GetMother();
 942     if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
 943     AliAODMCParticle *mothertrack = NULL;
 944     if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
 945
 946     // Check C
 947     Int_t pdg = mothertrack->GetPdgCode();
 948     if((TMath::Abs(pdg)==411) || (TMath::Abs(pdg)==421) || (TMath::Abs(pdg)==431) || (TMath::Abs(pdg)==4122) || (TMath::Abs(pdg)==4132) || (TMath::Abs(pdg)==4232) || (TMath::Abs(pdg)==43320)) return imother;
 949     if(TMath::Abs(pdg) == 11)
 950     {
 951       return IsMotherC(imother);
 952     }
 953
 954     return -1;
 955   }
 956
 957   return -1;
 958 }
 959 //_______________________________________________________________________________________________
 960 Int_t AliHFENonPhotonicElectron::IsMotherB(Int_t tr) {
 961
 962   //
 963   // Return the lab of signal mother or -1 if not B
 964   //
 965
 966   if(tr < 0) return -1;
 967
 968   if(fMCEvent)
 969   {
 970     AliVParticle *mctrack = fMCEvent->GetTrack(tr);
 971     if(!mctrack) return -1;
 972     AliMCParticle *mctrackesd = NULL;
 973     if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
 974     TParticle *particle = 0x0;
 975     particle = mctrackesd->Particle();
 976
 977     // Take mother
 978     if(!particle) return -1;
 979     Int_t imother   = particle->GetFirstMother();
 980     if(imother < 0) return -1;
 981     AliMCParticle *mothertrack = NULL;
 982     if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
 983     TParticle * mother = mothertrack->Particle();
 984     if(!mother) return -1;
 985
 986     // Check B
 987     Int_t pdg = mother->GetPdgCode();
 988     if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother;
 989     if(TMath::Abs(pdg) == 11)
 990     {
 991       return IsMotherB(imother);
 992     }
 993
 994     return -1;
 995   }
 996
 997   if(fAODArrayMCInfo)
 998   {
 999     if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
1000     AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
1001     if(!mctrackaod) return -1;
1002
1003     // Take mother
1004     Int_t imother = mctrackaod->GetMother();
1005     if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
1006     AliAODMCParticle *mothertrack = NULL;
1007     if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
1008     // Check B
1009     Int_t pdg = mothertrack->GetPdgCode();
1010     if((TMath::Abs(pdg)==511) || (TMath::Abs(pdg)==521) || (TMath::Abs(pdg)==531) || (TMath::Abs(pdg)==5122) || (TMath::Abs(pdg)==5132) || (TMath::Abs(pdg)==5232) || (TMath::Abs(pdg)==53320)) return imother;
1011     if(TMath::Abs(pdg) == 11)
1012     {
1013       return IsMotherB(imother);
1014     }
1015
1016     return -1;
1017   }
1018
1019   return -1;
1020 }
1021
1022 //_______________________________________________________________________________________________
1023 Int_t AliHFENonPhotonicElectron::IsMotherEta(Int_t tr) {
1024
1025   //
1026   // Return the lab of eta mother or -1 if not eta
1027   //
1028
1029   if(tr < 0) return -1;
1030
1031   if(fMCEvent)
1032   {
1033     AliVParticle *mctrack = fMCEvent->GetTrack(tr);
1034     if(!mctrack) return -1;
1035     AliMCParticle *mctrackesd = NULL;
1036     if(!(mctrackesd = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(tr))))) return -1;
1037     TParticle *particle = 0x0;
1038     particle = mctrackesd->Particle();
1039
1040     // Take mother
1041     if(!particle) return -1;
1042     Int_t imother   = particle->GetFirstMother();
1043     if(imother < 0) return -1;
1044     AliMCParticle *mothertrack = NULL;
1045     if(!(mothertrack = dynamic_cast<AliMCParticle *>(fMCEvent->GetTrack(TMath::Abs(imother))))) return -1;
1046     TParticle * mother = mothertrack->Particle();
1047     if(!mother) return -1;
1048
1049     // Check eta
1050     Int_t pdg = mother->GetPdgCode();
1051     if(TMath::Abs(pdg) == 221) return imother;
1052     if(TMath::Abs(pdg) == 11)
1053     {
1054       return IsMotherEta(imother);
1055     }
1056
1057     return -1;
1058   }
1059
1060   if(fAODArrayMCInfo)
1061   {
1062     if((tr+1)>fAODArrayMCInfo->GetEntriesFast()) return -1;
1063     AliAODMCParticle *mctrackaod = (AliAODMCParticle *) fAODArrayMCInfo->At(tr);
1064     if(!mctrackaod) return -1;
1065
1066     // Take mother
1067     Int_t imother = mctrackaod->GetMother();
1068     if(imother < 0 || ((imother+1)>fAODArrayMCInfo->GetEntriesFast())) return -1;
1069     AliAODMCParticle *mothertrack = NULL;
1070     if(!(mothertrack = dynamic_cast<AliAODMCParticle *>(fAODArrayMCInfo->At(TMath::Abs(imother))))) return -1;
1071
1072     // Check eta
1073     Int_t pdg = mothertrack->GetPdgCode();
1074     if(TMath::Abs(pdg) == 221) return imother;
1075     if(TMath::Abs(pdg) == 11)
1076     {
1077       return IsMotherEta(imother);
1078     }
1079
1080     return -1;
1081   }
1082
1083   return -1;
1084 }