QA of TRD tracking performance
authorcblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 5 Aug 2008 16:41:09 +0000 (16:41 +0000)
committercblume <cblume@f7af4fe6-9843-0410-8265-dc069ae4e863>
Tue, 5 Aug 2008 16:41:09 +0000 (16:41 +0000)
12 files changed:
TRD/TRDqaRecLinkDef.h [new file with mode: 0644]
TRD/libTRDqaRec.pkg [new file with mode: 0644]
TRD/qaRec/AliTRDtrackInfo/AliTRDtrackInfo.cxx [new file with mode: 0644]
TRD/qaRec/AliTRDtrackInfo/AliTRDtrackInfo.h [new file with mode: 0644]
TRD/qaRec/AliTRDtrackInfoGen.cxx [new file with mode: 0644]
TRD/qaRec/AliTRDtrackInfoGen.h [new file with mode: 0644]
TRD/qaRec/AliTRDtrackingEfficiency.cxx [new file with mode: 0644]
TRD/qaRec/AliTRDtrackingEfficiency.h [new file with mode: 0644]
TRD/qaRec/AliTRDtrackingEfficiencyCombined.cxx [new file with mode: 0644]
TRD/qaRec/AliTRDtrackingEfficiencyCombined.h [new file with mode: 0644]
TRD/qaRec/CreateESDChain.C [new file with mode: 0644]
TRD/qaRec/run.C [new file with mode: 0644]

diff --git a/TRD/TRDqaRecLinkDef.h b/TRD/TRDqaRecLinkDef.h
new file mode 100644 (file)
index 0000000..2d63f4b
--- /dev/null
@@ -0,0 +1,14 @@
+#ifdef __CINT__
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+
+#pragma link C++ class  AliTRDtrackInfoGen+;
+#pragma link C++ class  AliTRDtrackingEfficiency+;
+#pragma link C++ class  AliTRDtrackingEfficiencyCombined+;
+#pragma link C++ class  AliTRDtrackInfo+;
+
+#endif
diff --git a/TRD/libTRDqaRec.pkg b/TRD/libTRDqaRec.pkg
new file mode 100644 (file)
index 0000000..b59f988
--- /dev/null
@@ -0,0 +1,17 @@
+#-*- Mode: Makefile -*-
+
+SRCS= qaRec/AliTRDtrackInfoGen.cxx \
+      qaRec/AliTRDtrackingEfficiency.cxx \
+      qaRec/AliTRDtrackingEfficiencyCombined.cxx \
+      qaRec/AliTRDtrackInfo/AliTRDtrackInfo.cxx
+
+HDRS= $(SRCS:.cxx=.h)
+
+DHDR= TRDqaRecLinkDef.h
+
+
+ifeq (win32gcc,$(ALICE_TARGET))
+PACKSOFLAGS:= $(SOFLAGS) -L$(ALICE_ROOT)/lib/tgt_$(ALICE_TARGET) \
+                         -lANALYSIS -lANALYSISalice -lESD \
+#                         -L$(shell root-config --libdir) -lGui
+endif
diff --git a/TRD/qaRec/AliTRDtrackInfo/AliTRDtrackInfo.cxx b/TRD/qaRec/AliTRDtrackInfo/AliTRDtrackInfo.cxx
new file mode 100644 (file)
index 0000000..f02e32d
--- /dev/null
@@ -0,0 +1,273 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id: AliTRDtrackInfo.cxx 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+//  Authors:                                                              //
+//    Markus Fasel <M.Fasel@gsi.de>                                       //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include "AliTrackReference.h"
+#include "AliExternalTrackParam.h"
+#include "AliTRDseedV1.h"
+#include "AliTRDtrackV1.h"
+
+#include <cstdio>
+#include <cstdlib>
+
+#include "AliTRDtrackInfo.h"
+
+ClassImp(AliTRDtrackInfo)
+
+//___________________________________________________
+AliTRDtrackInfo::AliTRDtrackInfo():
+  TObject()
+  ,fPDG(0)
+  ,fStatus(0)
+  ,fId(-1)
+  ,fLabel(0)
+  ,fNClusters(0)
+  ,fNTrackRefs(0)
+  ,fTRDtrack(0x0)
+  ,fOP(0x0)
+{
+  //
+  // Default constructor
+  //
+
+  // Set 0-Pointers
+  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
+}
+
+//___________________________________________________
+AliTRDtrackInfo::AliTRDtrackInfo(Int_t pdg):
+  TObject()
+  ,fPDG(pdg)
+  ,fStatus(0)
+  ,fId(-1)
+  ,fLabel(0)
+  ,fNClusters(0)
+  ,fNTrackRefs(0)
+  ,fTRDtrack(0x0)
+  ,fOP(0x0)
+{
+  //
+  // PDG constructor
+  //
+
+  // Set 0-Pointers
+  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
+}
+
+//___________________________________________________
+AliTRDtrackInfo::AliTRDtrackInfo(const AliTRDtrackInfo &trdInfo):
+  TObject((const TObject&)trdInfo)  
+  ,fPDG(trdInfo.fPDG)
+  ,fStatus(trdInfo.fStatus)
+  ,fId(trdInfo.fId)
+  ,fLabel(trdInfo.fLabel)
+  ,fNClusters(trdInfo.fNClusters)
+  ,fNTrackRefs(trdInfo.fNTrackRefs)
+  ,fTRDtrack(0x0)
+  ,fOP(0x0)
+{
+  //
+  // copy Entries
+  //
+
+  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
+  for(Int_t ien = 0; ien < 12; ien++){
+       if(trdInfo.fTrackRefs[ien])
+      fTrackRefs[ien] = new AliTrackReference(*(trdInfo.fTrackRefs[ien]));
+  }
+  if(trdInfo.fOP) fOP = new AliExternalTrackParam(*trdInfo.fOP);
+  if(trdInfo.fTRDtrack) fTRDtrack = new AliTRDtrackV1(*trdInfo.fTRDtrack);
+}
+
+//___________________________________________________
+AliTRDtrackInfo::~AliTRDtrackInfo()
+{
+  //
+  // Destructor
+  //
+
+  if(fOP) delete fOP;
+  for(Int_t ien = 0; ien < 12; ien++){
+    if(fTrackRefs[ien]) delete fTrackRefs[ien];
+  }
+  if(fTRDtrack) delete fTRDtrack;
+}
+
+
+//___________________________________________________
+AliTRDtrackInfo& AliTRDtrackInfo::operator=(const AliTRDtrackInfo &trdInfo)
+{
+  //
+  // = Operator
+  //
+
+  fPDG    = trdInfo.fPDG;
+  fStatus = trdInfo.fStatus;
+  fId     = trdInfo.fId;
+  fLabel  = trdInfo.fLabel;
+  fNClusters  = trdInfo.fNClusters;
+  fNTrackRefs = trdInfo.fNTrackRefs;
+
+  // copy Entries
+  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
+  for(Int_t ien = 0; ien < 12; ien++){
+    if(trdInfo.fTrackRefs[ien])
+       if(!fTrackRefs[ien])
+               fTrackRefs[ien] = new AliTrackReference(*(trdInfo.fTrackRefs[ien]));
+       else
+               new(&fTrackRefs[ien]) AliTrackReference(*(trdInfo.fTrackRefs[ien]));
+  }
+  if(trdInfo.fOP){
+       if(!fOP)
+               fOP = new AliExternalTrackParam(*trdInfo.fOP);
+       else
+               new(fOP) AliExternalTrackParam(*trdInfo.fOP);
+  }
+  if(trdInfo.fTRDtrack){
+       if(!fTRDtrack)
+               fTRDtrack = new AliTRDtrackV1(*trdInfo.fTRDtrack);
+       else
+               new(fTRDtrack) AliTRDtrackV1(*trdInfo.fTRDtrack);
+  }
+
+  return *this;
+}
+
+//___________________________________________________
+void AliTRDtrackInfo::Delete(const Option_t *)
+{
+  //
+  // Delete
+  //
+
+  fPDG    = 0;
+  fStatus = 0;
+  fId     = -1;
+  fLabel  = 0;
+  fNClusters  = 0;
+  fNTrackRefs = 0;
+  if(fOP) delete fOP; fOP = 0x0;
+  if(fTRDtrack) delete fTRDtrack; fTRDtrack = 0x0;
+  for(Int_t ien = 0; ien < 12; ien++){
+    if(fTrackRefs[ien])
+      delete fTrackRefs[ien];
+  }
+  memset(fTrackRefs, 0, sizeof(AliTrackReference *) * 12);
+}
+
+
+//___________________________________________________
+void AliTRDtrackInfo::SetTRDtrack(const AliTRDtrackV1 *track)
+{
+  //
+  // Set the TRD track
+  //
+
+  if(!fTRDtrack)
+       fTRDtrack = new AliTRDtrackV1(*track);  
+  else
+       new(fTRDtrack)AliTRDtrackV1(*track);
+       fTRDtrack->SetOwner();
+  // Make a copy for the object in order to avoid ownership problems
+}
+
+//___________________________________________________
+void AliTRDtrackInfo::AddTrackRef(const AliTrackReference *tref)
+{
+  //
+  // Add track reference
+  //
+
+  if(fNTrackRefs >= 12){ 
+    SetCurved();
+    return;
+  }
+  // Make a copy for the object in order to avoid ownership problems
+  fTrackRefs[fNTrackRefs++] = new AliTrackReference(*tref);
+}
+
+//___________________________________________________
+AliTRDseedV1* AliTRDtrackInfo::GetTracklet(Int_t idx) const 
+{
+  //
+  // Returns a tracklet
+  //
+
+       if(!fTRDtrack) return 0x0;
+  return idx < 6 ? const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(idx)) : 0x0;
+}
+
+//___________________________________________________
+AliTrackReference * AliTRDtrackInfo::GetTrackRef(Int_t idx) const
+{
+  //
+  // Returns a track reference
+  //
+
+  return idx < 12 ? fTrackRefs[idx] : 0x0;
+}
+
+//___________________________________________________
+Int_t AliTRDtrackInfo::GetNumberOfClusters() const
+{
+  //
+  // Returns the number of clusters
+  //
+
+  Int_t n = 0;
+       if(!fTRDtrack) return 0;
+  if(fTRDtrack->GetNumberOfTracklets() == 0) return n;
+  AliTRDseedV1 *tracklet = 0x0;
+  for(Int_t ip=0; ip<6; ip++){
+       if(!(tracklet = const_cast<AliTRDseedV1 *>(fTRDtrack->GetTracklet(ip)))) continue;
+    n+=tracklet->GetN();
+  }
+  return n;
+}
+
+
+//___________________________________________________
+void  AliTRDtrackInfo::SetOuterParam(const AliExternalTrackParam *op)
+{
+  //
+  // Set outer track parameters
+  //
+
+  if(!op) return;
+  if(fOP) new(fOP) AliExternalTrackParam(*op);
+  else
+       fOP = new AliExternalTrackParam(*op);
+}
+
+//___________________________________________________
+Int_t AliTRDtrackInfo::GetNTracklets() const
+{
+  //
+  // Return the number of tracklets
+  //
+
+       if(!fTRDtrack) return 0x0;
+       return fTRDtrack->GetNumberOfTracklets();
+}
diff --git a/TRD/qaRec/AliTRDtrackInfo/AliTRDtrackInfo.h b/TRD/qaRec/AliTRDtrackInfo/AliTRDtrackInfo.h
new file mode 100644 (file)
index 0000000..abffa4b
--- /dev/null
@@ -0,0 +1,80 @@
+#ifndef ALITRDTRACKINFO_H
+#define ALITRDTRACKINFO_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id: AliTRDtrackInfo.h 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#ifndef Root_TObject
+#include "TObject.h"
+#endif
+
+class AliTRDseedV1;
+class AliTRDtrackV1;
+class AliTrackReference;
+class AliExternalTrackParam;
+class AliTRDtrackInfo : public TObject{
+public:
+  AliTRDtrackInfo();
+  AliTRDtrackInfo(Int_t pdg);
+  AliTRDtrackInfo(const AliTRDtrackInfo &);
+  ~AliTRDtrackInfo();
+  
+//  void               Clear(const Option_t *){}
+  void               Delete(const Option_t *);
+  
+  AliTRDtrackInfo&   operator=(const AliTRDtrackInfo &);
+  
+  void               AddTrackRef(const AliTrackReference *trackRef);
+  
+  Int_t              GetTrackId() { return fId;}
+  Int_t              GetNumberOfClusters() const;
+  Int_t              GetNumberOfClustersRefit() const {return fNClusters;}
+  Int_t              GetNTracklets() const;
+  Int_t              GetNTrackRefs() const {return fNTrackRefs;} 
+  Int_t              GetLabel() const { return fLabel; }
+  Int_t              GetPDG() const { return fPDG; }
+  ULong_t            GetStatus() const {return fStatus;}
+  AliTRDseedV1*      GetTracklet(Int_t entry) const;
+  AliTRDtrackV1 *               GetTRDtrack() const { return fTRDtrack; }
+  AliTrackReference* GetTrackRef(Int_t entry) const;
+  AliExternalTrackParam* GetOuterParam() const {return fOP;}
+
+  Bool_t             IsCurved() const {return TestBit(kCurv);}
+  Bool_t                                                IsPrimary() const {return TestBit(kPrim);}
+       Bool_t                                           HasESDtrack() const{return ((fTRDtrack != 0x0) ||(fOP != 0));}
+
+  void               SetCurved(Bool_t curv = kTRUE) {SetBit(kCurv, curv);}
+  void               SetLabel(Int_t lab) { fLabel = lab; }
+  void               SetNumberOfClustersRefit(Int_t n) {fNClusters = n;}
+  void               SetPDG(Int_t pdg) { fPDG = pdg; }
+  void                                                  SetPrimary(Bool_t prim = kTRUE) {SetBit(kPrim, prim);}
+  void               SetOuterParam(const AliExternalTrackParam *op);
+  void               SetStatus(ULong_t stat) {fStatus = stat;}
+  void               SetTrackId(Int_t id) {fId = id;}
+  void               SetTRDtrack(const AliTRDtrackV1 *track);
+  
+private:
+       enum{
+               kCurv = 1,
+               kPrim = 2
+       };
+  Int_t              fPDG;             // particle code
+  ULong_t            fStatus;          // ESD track status
+  Int_t              fId;              // ESD track id
+  Int_t              fLabel;           // MC label  
+  Int_t              fNClusters;       // Numer of clusters from refit
+  Int_t              fNTrackRefs;      // number of track refs
+  AliTrackReference  *fTrackRefs[12];  // no of track refs
+  AliTRDtrackV1      *fTRDtrack;       // no of tracklets
+  AliExternalTrackParam *fOP;          // outer param if no tracklets
+  ClassDef(AliTRDtrackInfo, 1)          // TRD track info
+};
+#endif
diff --git a/TRD/qaRec/AliTRDtrackInfoGen.cxx b/TRD/qaRec/AliTRDtrackInfoGen.cxx
new file mode 100644 (file)
index 0000000..6adccdd
--- /dev/null
@@ -0,0 +1,361 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id: AliTRDtrackInfoGen.cxx 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+//  Authors:                                                              //
+//    Markus Fasel <M.Fasel@gsi.de>                                       //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include <TClonesArray.h>
+#include <TObjArray.h>
+#include <TObject.h>
+#include <TH1F.h>
+#include <TFile.h>
+#include <TTree.h>
+#include <TROOT.h>
+#include <TChain.h>
+#include <TParticle.h>
+
+#include "AliLog.h"
+#include "AliAnalysisManager.h"
+#include "AliESDEvent.h"
+#include "AliMCEvent.h"
+#include "AliESDInputHandler.h"
+#include "AliMCEventHandler.h"
+
+#include "AliESDfriend.h"
+#include "AliESDfriendTrack.h"
+#include "AliESDtrack.h"
+#include "AliMCParticle.h"
+#include "AliStack.h"
+#include "AliTRDtrackV1.h"
+#include "AliTrackReference.h"
+#include "AliTRDgeometry.h"
+#include "TTreeStream.h"
+
+#include <cstdio>
+#include <cstring>
+
+#include "AliTRDtrackInfoGen.h"
+#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
+
+ClassImp(AliTRDtrackInfoGen)
+
+
+//____________________________________________________________________
+AliTRDtrackInfoGen::AliTRDtrackInfoGen(const Char_t *name):
+  AliAnalysisTask(name, "")
+  ,fESD(0x0)
+  ,fMC(0x0)
+  ,fESDfriend(0x0)
+  ,fTrackInfo(0x0)
+  ,fObjectContainer(0x0)
+  ,fDebugLevel(0)
+  ,fDebugStream(0x0)
+  //,fTree(0x0) // temporary
+{
+  //
+  // Default constructor
+  //
+
+  DefineInput(0, TChain::Class());
+  DefineOutput(0, TObjArray::Class());
+  //DefineOutput(1, TTree::Class());
+}
+
+//____________________________________________________________________
+void AliTRDtrackInfoGen::ConnectInputData(Option_t *)
+{
+  //
+  // Link the Input Data
+  //
+       TTree *tree = dynamic_cast<TChain*>(GetInputData(0));
+       if(!tree){
+    printf("ERROR - ESD event not found");
+       } else {
+               tree->SetBranchStatus("Tracks", 1);
+               tree->SetBranchStatus("ESDfriend*",1);
+  }
+
+  AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*>(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+  if(!esdH){
+               printf("ERROR - ESD input handler not found");
+  } else {
+    fESD = esdH->GetEvent();
+    if(!fESD){
+                 printf("ERROR - ESD event not found");
+    } else {
+      esdH->SetActiveBranches("ESDfriend*");
+      fESDfriend = (AliESDfriend *)fESD->FindListObject("AliESDfriend");
+      printf("fESDfriend = %p\n", (void*)fESDfriend);
+    }
+  }
+  AliMCEventHandler *mcH = dynamic_cast<AliMCEventHandler*>(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler());
+  if(!mcH){ 
+    AliError("MC input handler not found");
+  } else {
+    fMC = mcH->MCEvent();
+  }
+}
+
+//____________________________________________________________________
+void AliTRDtrackInfoGen::CreateOutputObjects()
+{      
+  //
+  // Create Output Containers (TObjectArray containing 1D histograms)
+  //
+  fTrackInfo = new AliTRDtrackInfo();
+  fObjectContainer = new TObjArray(1000);
+
+/*  OpenFile(1, "RECREATE");
+  fTree = new TTree("trd", "extract of the TRD detector");
+  fTree->Branch("info",  &fTrackInfo);
+  printf("output tree build in %s\n", fTree->GetDirectory()->GetName());*/
+}
+
+//____________________________________________________________________
+void AliTRDtrackInfoGen::Exec(Option_t *){
+       //
+       // Run the Analysis
+       //
+       if(!fESD){
+               puts("Error: ESD not found");
+               return;
+       }
+       if(!fESDfriend){
+               puts("Error: ESD friend not found");
+               return;
+       }
+       if(!fMC){
+               puts("Error: Monte Carlo Event not available");
+               return;
+       }
+       fObjectContainer->Delete();
+       fESD->SetESDfriend(fESDfriend);
+
+       AliStack *mStack = fMC->Stack();
+       if(!mStack){
+               puts("Error: Cannot get the Monte Carlo Stack");
+               return;
+       }
+       Bool_t *trackMap = new Bool_t[fMC->GetNumberOfTracks()];
+       memset(trackMap, 0, sizeof(Bool_t) * fMC->GetNumberOfTracks());
+
+       Int_t nTRD = 0, nTPC = 0, nclsTrklt;
+       Int_t nTracks = fESD->GetNumberOfTracks();
+       if(fDebugLevel>=1) printf("%3d Tracks: ESD[%d] MC[%d]\n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), nTracks, mStack->GetNtrack());
+       AliESDtrack *esdTrack = 0x0;
+       AliESDfriendTrack *esdFriendTrack = 0x0;
+       TObject *calObject = 0x0;
+       AliTRDtrackV1 *track = 0x0;
+       for(Int_t itrk = 0; itrk < nTracks; itrk++){
+               esdTrack = fESD->GetTrack(itrk);
+               if(fDebugLevel>=2) printf("\n%3d ITS[%d] TPC[%d] TRD[%d]\n", itrk, esdTrack->GetNcls(0), esdTrack->GetNcls(1), esdTrack->GetNcls(2));
+               if(esdTrack->GetNcls(1)) nTPC++;
+               if(esdTrack->GetNcls(2)) nTRD++;
+
+               // look at esternal track param
+               const AliExternalTrackParam *op = esdTrack->GetOuterParam();
+               Double_t xyz[3];
+               if(op){
+                               op->GetXYZ(xyz);
+                               op->Global2LocalPosition(xyz, op->GetAlpha());
+                               if(fDebugLevel>=2) printf("op @ X[%7.3f]\n", xyz[0]);
+               }
+/*             if(xyz[0] < 270.){ 
+                               printf("TPC track missing\n");
+               } else { */
+//                     op->GetXYZAt(ref->LocalX(), AliTracker::GetBz(), xyz);
+//                     op->Global2LocalPosition(xyz, op->GetAlpha());
+//                     dy = ref->LocalY()- xyz[1];
+//                     dz = ref->Z() - xyz[2];
+//     }
+
+               // read MC info
+               Int_t label = esdTrack->GetLabel();
+               if(label < fMC->GetNumberOfTracks()) trackMap[TMath::Abs(label)] = kTRUE; // register the track
+               //if (TMath::Abs(label) > mStack->GetNtrack()) continue; 
+               AliMCParticle *mcParticle = fMC->GetTrack(TMath::Abs(label));
+               Int_t fPdg = mcParticle->Particle()->GetPdgCode();
+               Int_t nRefs = mcParticle->GetNumberOfTrackReferences();
+               Int_t iref = 0; AliTrackReference *ref = 0x0; 
+               while(iref<nRefs){
+                       ref = mcParticle->GetTrackReference(iref);
+                       if(ref->LocalX() > 250.) break;
+                       //printf("\ttrackRef[%2d] @ %7.3f\n", iref, ref->LocalX());
+                       iref++;
+               }
+               if(iref == nRefs){
+//                     if(!esdTrack->GetNcls(2)) continue;
+/*                     printf("No TRD Track References in the Track [%d] I\n", itrk);
+                       printf("Label = %d ITS[%d] TPC[%d] TRD[%d]\n", label, esdTrack->GetITSLabel(), esdTrack->GetTPCLabel(), esdTrack->GetTRDLabel());
+                       Int_t kref = 0;
+                       while(kref<nRefs){
+                               ref = mcParticle->GetTrackReference(kref);
+                               printf("\ttrackRef[%2d] @ %7.3f\n", kref, ref->LocalX());
+                               kref++;
+                       }*/
+               }
+
+               new(fTrackInfo) AliTRDtrackInfo(fPdg);
+               fTrackInfo->SetPrimary(mcParticle->Particle()->IsPrimary());
+               Int_t jref = iref;//, kref = 0;
+               while(jref<nRefs){
+                       ref = mcParticle->GetTrackReference(jref);
+                       if(ref->LocalX() > 370.) break;
+                       if(fDebugLevel>=3) printf("\ttrackRef[%2d (%2d)] @ %7.3f OK\n", jref-iref, jref, ref->LocalX());
+                       fTrackInfo->AddTrackRef(ref);
+                       jref++;
+               }
+               if(!fTrackInfo->GetNTrackRefs()){ 
+                               //if(!esdTrack->GetNcls(2)) continue;
+/*             printf("No TRD Track References in the Track [%d] II\n", itrk);
+               printf("Label = %d ITS[%d] TPC[%d] TRD[%d]\n", label, esdTrack->GetITSLabel(), esdTrack->GetTPCLabel(), esdTrack->GetTRDLabel());
+               Int_t kref = 0;
+               while(kref<nRefs){
+                       ref = mcParticle->GetTrackReference(kref);
+                       printf("\ttrackRef[%2d] @ %7.3f\n", kref, ref->LocalX());
+                       kref++;
+               }*/
+               }
+               if(fDebugLevel>=2) printf("NtrackRefs[%d(%d)]\n", fTrackInfo->GetNTrackRefs(), nRefs);
+
+               // copy some relevant info to TRD track info
+               fTrackInfo->SetStatus(esdTrack->GetStatus());
+               fTrackInfo->SetTrackId(esdTrack->GetID());
+               fTrackInfo->SetLabel(label);
+               fTrackInfo->SetNumberOfClustersRefit(esdTrack->GetNcls(2));
+               nclsTrklt = 0;
+       
+
+               // read REC info
+               esdFriendTrack = fESDfriend->GetTrack(itrk);
+               if(esdFriendTrack){
+                       Int_t icalib = 0;
+                       while((calObject = esdFriendTrack->GetCalibObject(icalib++))){
+                               if(strcmp(calObject->IsA()->GetName(),"AliTRDtrackV1") != 0) continue; // Look for the TRDtrack
+                               track = dynamic_cast<AliTRDtrackV1*>(calObject);
+                               //nTRD++;
+                               if(!esdTrack->GetNcls(2)){
+                                       printf("No TRD clusters but track\n");
+                                       nTRD++;  
+                               }  
+                               if(!track) continue;
+                               if(fDebugLevel>=3) printf("TRD track OK\n");
+                               fTrackInfo->SetTRDtrack(track);
+                               break;
+                       }
+                       if(fDebugLevel>=2) printf("Ntracklets[%d]\n", fTrackInfo->GetNTracklets());
+               } else if(fDebugLevel>=2) printf("No ESD friends\n");
+               if(op) fTrackInfo->SetOuterParam(op);
+
+               if(fDebugLevel >= 1){
+                       Int_t ncls = esdTrack->GetNcls(2);
+                       nclsTrklt = fTrackInfo->GetNumberOfClusters();
+                       (*fDebugStream) << "trackInfo"
+                       << "ncls="       << ncls
+                       << "Ncls="       << nclsTrklt
+                       << "TrackInfo.=" << fTrackInfo
+                       << "\n";
+               }
+       
+               fObjectContainer->Add(new AliTRDtrackInfo(*fTrackInfo));
+               fTrackInfo->Delete("");
+       }
+       
+       // Insert also MC tracks which are passing TRD where the track is not reconstructed
+       if(fDebugLevel>=1) printf("%3d Tracks: TPC[%d] TRD[%d]\n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), nTPC, nTRD);
+       if(fDebugLevel > 10){
+               printf("Output of the MC track map:\n");
+               for(Int_t itk = 0; itk < fMC->GetNumberOfTracks();  itk++)
+                       printf("trackMap[%d] = %s\n", itk, trackMap[itk] == kTRUE ? "TRUE" : "kFALSE");
+  }
+       
+       for(Int_t itk = 0; itk < fMC->GetNumberOfTracks(); itk++){
+               if(fDebugLevel >=2 ) printf("Number of MC tracks: %d\n", fMC->GetNumberOfTracks());
+               if(trackMap[itk]) continue;
+               AliMCParticle *mcParticle = fMC->GetTrack(TMath::Abs(itk));
+               Int_t fPdg = mcParticle->Particle()->GetPdgCode();
+               Int_t nRefs = mcParticle->GetNumberOfTrackReferences();
+               Int_t iref = 0; AliTrackReference *ref = 0x0; 
+               Int_t nRefsTRD = 0;
+               new(fTrackInfo) AliTRDtrackInfo(fPdg);
+               while(iref<nRefs){
+                       ref = mcParticle->GetTrackReference(iref);
+                       if(fDebugLevel > 3) printf("\ttrackRef[%2d] @ %7.3f", iref, ref->LocalX());
+                       if(ref->LocalX() > 250. && ref->LocalX() < 370.){
+                               if(fDebugLevel > 3) printf(" OK\n");
+                               fTrackInfo->AddTrackRef(ref);
+                               nRefsTRD++;
+                       }
+                       else
+                               if(fDebugLevel > 3) printf("\n");
+                       iref++;
+               }
+               if(!nRefsTRD){
+                       // In this stage we at least require 1 hit inside TRD. What will be done with this tracks is a task for the 
+                       // analysis job
+                       fTrackInfo->Delete("");
+                       continue;
+               }
+               fTrackInfo->SetPrimary(mcParticle->Particle()->IsPrimary());
+               fTrackInfo->SetLabel(itk);
+               if(fDebugLevel >= 1){
+                       Int_t ncls = esdTrack->GetNcls(2);
+                       (*fDebugStream) << "trackInfo"
+                       << "ntrackRefs=" << ncls
+                       << "TrackInfo.=" << fTrackInfo
+                       << "\n";
+               }
+               if(fDebugLevel > 2)printf("Registering rejected MC track with label %d\n", itk);
+               fObjectContainer->Add(new AliTRDtrackInfo(*fTrackInfo));
+               fTrackInfo->Delete("");
+       }
+       delete[] trackMap;
+       PostData(0, fObjectContainer);
+}
+
+//____________________________________________________________________
+void  AliTRDtrackInfoGen::SetDebugLevel(Int_t level) 
+{
+  //
+  // Set the debug level
+  //
+
+  fDebugLevel = level;
+  if(fDebugLevel<=0) return;
+  if(fDebugStream) return;  
+  fDebugStream = new TTreeSRedirector("TRD.TrackInfoDebug.root");
+}
+
+//____________________________________________________________________
+void AliTRDtrackInfoGen::Terminate(Option_t *)
+{
+  //
+  // Stays empty because we are only interested in the tree
+  //
+
+  if(fDebugLevel>=1) printf("Terminate:\n");
+  //TFile *f =((TFile*)gROOT->FindObject("TRD.TrackInfo.root"));
+  //f->cd(); f->Write(); f->Close();
+
+  if(fDebugStream) delete fDebugStream;
+}
diff --git a/TRD/qaRec/AliTRDtrackInfoGen.h b/TRD/qaRec/AliTRDtrackInfoGen.h
new file mode 100644 (file)
index 0000000..554110b
--- /dev/null
@@ -0,0 +1,54 @@
+#ifndef ALITRDTRACKINFOGEN_H
+#define ALITRDTRACKINFOGEN_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id: AliTRDtrackInfoGen.h 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include "AliAnalysisTask.h"
+
+class AliESDEvent;
+class AliMCEvent;
+class AliESDfriend;
+class AliTRDtrackInfo;
+class TObjArray;
+class TTreeSRedirector;
+
+class AliTRDtrackInfoGen : public AliAnalysisTask{
+public:
+
+  AliTRDtrackInfoGen(const Char_t *name = "TRD Track Info");
+  ~AliTRDtrackInfoGen(){};
+  
+  void  ConnectInputData(Option_t *);
+  void  CreateOutputObjects();
+  Int_t GetDebugLevel() const {return fDebugLevel;} 
+  void  Exec(Option_t *);
+  void  SetDebugLevel(Int_t level);
+  void  Terminate(Option_t *);
+
+private:
+
+  AliTRDtrackInfoGen(const AliTRDtrackInfoGen&);
+  AliTRDtrackInfoGen& operator=(const AliTRDtrackInfoGen&);
+
+private:
+
+  AliESDEvent      *fESD;                  // ESD event
+  AliMCEvent       *fMC;                   // MC event
+  AliESDfriend     *fESDfriend;            // ESD friends
+  AliTRDtrackInfo  *fTrackInfo;            // Track info
+  TObjArray        *fObjectContainer;      // Object container
+
+  Int_t            fDebugLevel;            // Debug level
+  TTreeSRedirector *fDebugStream;          // Debug stream
+
+  ClassDef(AliTRDtrackInfoGen, 1)          // entry to TRD analysis
+};
+#endif
diff --git a/TRD/qaRec/AliTRDtrackingEfficiency.cxx b/TRD/qaRec/AliTRDtrackingEfficiency.cxx
new file mode 100644 (file)
index 0000000..92c3842
--- /dev/null
@@ -0,0 +1,272 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id: AliTRDtrackingEfficiency.cxx 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+//  Authors:                                                              //
+//    Markus Fasel <M.Fasel@gsi.de>                                       //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include <TClonesArray.h>
+#include <TObjArray.h>
+#include <TList.h>
+#include <TProfile.h>
+
+#include "AliESDtrack.h"
+#include "AliTrackReference.h"
+#include "AliExternalTrackParam.h"
+#include "AliTracker.h"
+#include "AliMagFMaps.h"
+#include "AliAnalysisManager.h"
+
+#include "AliTRDtrackingEfficiency.h"
+#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
+
+ClassImp(AliTRDtrackingEfficiency)
+
+//____________________________________________________________________
+AliTRDtrackingEfficiency::AliTRDtrackingEfficiency(const Char_t *name):
+  AliAnalysisTask(name, "")
+  ,fObjectContainer(0x0)
+  ,fTracks(0x0)
+  ,fMissed(0x0)
+  ,fDebugLevel(1)
+  ,fDebugStream(0x0)
+{
+  //
+  // Default constructor
+  //
+
+  DefineInput(0, TObjArray::Class());
+  DefineOutput(0, TList::Class());
+}
+
+//____________________________________________________________________
+void  AliTRDtrackingEfficiency::ConnectInputData(Option_t *)
+{
+  //
+  // Connect input data
+  //
+
+       fTracks = dynamic_cast<TObjArray*>(GetInputData(0));
+}
+
+//____________________________________________________________________
+void  AliTRDtrackingEfficiency::CreateOutputObjects()
+{
+  //
+  // Create output objects
+  //
+
+  OpenFile(0, "RECREATE");
+  const Int_t nbins = 11;
+  Float_t xbins[nbins+1] = {.5, .7, .9, 1.3, 1.7, 2.4, 3.5, 4.5, 5.5, 7., 9., 11.};
+  fObjectContainer = new TList();
+  fObjectContainer->Add(new TProfile("h", "", nbins, xbins));
+} 
+
+//____________________________________________________________________
+void AliTRDtrackingEfficiency::Exec(Option_t *)
+{
+  //
+  // Do it
+  //
+
+  if(!AliTracker::GetFieldMap()){
+    AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
+    AliTracker::SetFieldMap(field, kTRUE);
+  }
+  TProfile *h = (TProfile*)fObjectContainer->At(0);    
+       Int_t labelsacc[10000]; memset(labelsacc, 0, sizeof(Int_t) * 10000);
+       
+  if(!fMissed){ 
+    fMissed = new TClonesArray("AliTRDtrackInfo", 10);
+    fMissed->SetOwner();
+  }
+
+  Float_t mom;
+  Int_t selection[10000], nselect = 0;
+  ULong_t status;
+  Int_t nTRD = 0, nTPC = 0, nMiss = 0;
+  AliTRDtrackInfo     *track = 0x0;
+  AliTrackReference     *ref = 0x0;
+  AliExternalTrackParam *esd = 0x0;
+  for(Int_t itrk=0; itrk<fTracks->GetEntriesFast(); itrk++){
+    track = (AliTRDtrackInfo*)fTracks->UncheckedAt(itrk);
+               if(!track->HasESDtrack()) continue;
+    status = track->GetStatus();
+    if(!(status&AliESDtrack::kTPCout)){ 
+      // missing TPC propagation - interesting for SA
+      continue;
+    }
+
+    nTPC++;
+    selection[nselect++]=itrk;
+
+    //Int_t n = track->GetNumberOfClusters(); 
+    // where are this tracklets ???
+    //if(ncls0 > ncls1) printf("%3d ESD[%3d] TRD[%3d|%3d]\n", itrk, ncls0, ncls1, n);
+    if(track->GetNumberOfClustersRefit()){ 
+      ref = track->GetTrackRef(0);
+      esd = track->GetOuterParam();
+      mom = ref ? ref->P(): esd->P();
+      //printf("FOUND Id[%d] mom[%f]\n", itrk, mom);
+      h->Fill(mom, 1.);
+                       labelsacc[nTRD] = track->GetLabel();
+      nTRD++;
+      continue;
+    }
+
+
+    Int_t nrefs = track->GetNTrackRefs();
+    if(nrefs<=1){ 
+      nTPC--;
+      // we don't have MC info.
+      // we should discard this track 
+      continue;
+    }
+
+    Float_t xmed, xleng;
+    Int_t iref = 1; Bool_t found = kFALSE;
+    while((ref = track->GetTrackRef(iref))){
+      xmed = .5*(ref->LocalX() + track->GetTrackRef(iref-1)->LocalX());
+      xleng= (ref->LocalX() - track->GetTrackRef(iref-1)->LocalX());
+      if(TMath::Abs(xmed - 298.5) < .5 &&
+        TMath::Abs(xleng - 3.7) < .1){ 
+        found = kTRUE;
+        break;
+      }
+      iref++;
+    }
+    if(!found){ 
+      nTPC--;
+      // track missing first layer. Maybe interesting for SA.
+      continue;
+    }
+    nselect--;
+    new ((*fMissed)[nMiss]) AliTRDtrackInfo(*track);
+    nMiss++;
+  }
+  if(fDebugLevel>=1) printf("%3d Tracks: ESD[%3d] TPC[%3d] TRD[%3d | %5.2f%%] Off[%d]\n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), fTracks->GetEntriesFast(), nTPC, nTRD, nTPC ? 1.E2*nTRD/float(nTPC) : 0., fMissed->GetEntriesFast());
+
+
+  // Find double tracks
+  Float_t threshold = 10.;
+  AliTrackReference *refMiss = 0x0;
+  AliExternalTrackParam *op = 0x0;
+  AliTRDtrackInfo       *tt = 0x0;
+  for(Int_t imiss=0; imiss<nMiss; imiss++){
+    //printf("Searching missing %d ...\n", imiss);
+
+    // get outer param of missed
+    tt = (AliTRDtrackInfo*)fMissed->UncheckedAt(imiss);
+    op = tt->GetOuterParam();
+    Double_t alpha = op->GetAlpha(), cosa = TMath::Cos(alpha), sina = TMath::Sin(alpha);
+
+    Double_t xyz[3], x0, y0, z0, x, y, z, dx, dy, dz, d;
+
+    Bool_t FOUND = kFALSE;
+    for(Int_t iselect=0; iselect<nselect; iselect++){
+      track = (AliTRDtrackInfo*)fTracks->UncheckedAt(selection[iselect]);
+
+      // check first MC ... if available
+      d = 0;
+      for(Int_t iref=0; iref<track->GetNTrackRefs(); iref++){
+        if(!(ref = track->GetTrackRef(iref))) continue;
+        if((refMiss = tt->GetTrackRef(iref))){
+          dy = ref->LocalY() - refMiss->LocalY();
+          dz = ref->Z() - refMiss->Z();
+        } else {
+          // compare missOP with refTrackRef in LTC
+          x0 = ref->LocalX();
+          op->GetYAt(x0, AliTracker::GetBz(), y0);
+          op->GetZAt(x0, AliTracker::GetBz(), z0);
+          dy = y0 - ref->LocalY();
+          dz = z0 - ref->Z();
+        }
+        d += (dy*dy + dz*dz);
+      }
+      //printf("\td[%d] = %f N[%d]\n", selection[iselect], d, track->GetNTrackRefs());
+      if((track->GetNTrackRefs())){ 
+        d /= track->GetNTrackRefs();
+        if(d < threshold){
+          //printf("\t\tFound %2d in ref[%3d] : d[%f]\n", imiss, selection[iselect], d/track->GetNTrackRefs());
+          FOUND = kTRUE; break;
+        }
+      }
+
+      // process outer param ... always available
+      // compare missOP with OP in GTC
+      esd = track->GetOuterParam();
+      esd->GetXYZ(xyz);
+      x0 = esd->GetX();
+      op->GetYAt(x0, AliTracker::GetBz(), y0);
+      op->GetZAt(x0, AliTracker::GetBz(), z0);
+      x = x0*cosa - y0*sina;
+      y = x0*sina + y0*cosa;
+      z = z0;
+      dx=xyz[0]-x;
+      dy=xyz[1]-y;
+      dz=xyz[2]-z;
+      d = dx*dx+dy*dy+dz*dz;
+      //printf("\td[%d] = %f op\n", selection[iselect], d);
+      if(d < threshold){
+        //printf("\t\tFound %2d in op[%3d]  : d[%f] dx[%5.2f] dy[%5.2f] dz[%5.2f]\n", imiss, selection[iselect], d, dx, dy, dz);
+        FOUND = kTRUE; break;
+      }
+    }
+    if(FOUND) nTPC--;
+    else{ 
+      ref = tt->GetTrackRef(0);
+      mom = ref ? ref->P(): op->P();
+      h->Fill(mom, 0.);
+      printf("\tNOT FOUND Id[%d] Mom[%f]\n", tt->GetTrackId(), mom);
+    }
+  }
+
+  //if(fDebugLevel>=1)
+  printf("%3d Tracks: ESD[%3d] TPC[%3d] TRD[%3d | %5.2f%%] Off[%d]\n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), fTracks->GetEntriesFast(), nTPC, nTRD, 1.E2*nTRD/float(nTPC), fMissed->GetEntriesFast());
+
+  //fMissed->Delete();
+       // check for double countings
+       Int_t indices[10000]; memset(indices, 0, sizeof(Int_t) * 10000);
+       TMath::Sort(nTRD, labelsacc, indices);
+       if(fDebugLevel > 2){
+       for(Int_t itk = 0; itk < nTRD - 1; itk++)
+               if(labelsacc[indices[itk]] ==labelsacc[indices[itk + 1]]) printf("Double counted MC track: %d\n", labelsacc[indices[itk]]);
+       }
+  PostData(0, fObjectContainer);
+}
+
+//____________________________________________________________________
+void AliTRDtrackingEfficiency::Terminate(Option_t *)
+{
+  //
+  // Terminate
+  //
+
+  fObjectContainer = dynamic_cast<TList*>(GetOutputData(0));
+  if (!fObjectContainer) {
+    Printf("ERROR: list not available");
+    return;
+  }
+
+}
+
diff --git a/TRD/qaRec/AliTRDtrackingEfficiency.h b/TRD/qaRec/AliTRDtrackingEfficiency.h
new file mode 100644 (file)
index 0000000..59e7d88
--- /dev/null
@@ -0,0 +1,48 @@
+#ifndef ALITRDTRACKINGEFFICIENCY_H
+#define ALITRDTRACKINGEFFICIENCY_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id: AliTRDtrackingEfficiency.h 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include "AliAnalysisTask.h"
+
+class TObjArray;
+class TList;
+class TClonesArray;
+class TTreeSRedirector;
+class AliTRDtrackingEfficiency : public AliAnalysisTask
+{
+public:
+  AliTRDtrackingEfficiency(const Char_t *name = "TRD Tracking efficiency");
+  ~AliTRDtrackingEfficiency(){};
+  void  ConnectInputData(Option_t *);
+  void  CreateOutputObjects();
+  Int_t GetDebugLevel() const {return fDebugLevel;} 
+  void  Exec(Option_t *);
+  void  SetDebugLevel(Int_t debug){fDebugLevel = debug;}
+  void  Terminate(Option_t *);
+
+private:
+  AliTRDtrackingEfficiency(const AliTRDtrackingEfficiency&);
+  AliTRDtrackingEfficiency& operator=(const AliTRDtrackingEfficiency&);
+
+private:
+  TList        *fObjectContainer;       // Container
+  TObjArray        *fTracks;            // Array of tracks
+  TClonesArray     *fMissed;            // Missed ?
+
+  Int_t            fDebugLevel;         // Debug level
+  TTreeSRedirector *fDebugStream;       // Debug stream
+
+  ClassDef(AliTRDtrackingEfficiency, 1) // TRD tracking efficiency
+};
+
+#endif
+
diff --git a/TRD/qaRec/AliTRDtrackingEfficiencyCombined.cxx b/TRD/qaRec/AliTRDtrackingEfficiencyCombined.cxx
new file mode 100644 (file)
index 0000000..50feb9f
--- /dev/null
@@ -0,0 +1,266 @@
+/**************************************************************************
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ *                                                                        *
+ * Author: The ALICE Off-line Project.                                    *
+ * Contributors are mentioned in the code where appropriate.              *
+ *                                                                        *
+ * Permission to use, copy, modify and distribute this software and its   *
+ * documentation strictly for non-commercial purposes is hereby granted   *
+ * without fee, provided that the above copyright notice appears in all   *
+ * copies and that both the copyright notice and this permission notice   *
+ * appear in the supporting documentation. The authors make no claims     *
+ * about the suitability of this software for any purpose. It is          *
+ * provided "as is" without express or implied warranty.                  *
+ **************************************************************************/
+
+/* $Id: AliTRDtrackingEfficiencyCombined.cxx 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+//  Authors:                                                              //
+//    Markus Fasel <M.Fasel@gsi.de>                                       //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include <TObjArray.h>
+#include <TProfile.h>
+#include <TMath.h>
+#include <TCanvas.h>
+
+#include "AliMagFMaps.h"
+#include "AliTracker.h"
+#include "AliTRDseedV1.h"
+#include "AliTRDtrackV1.h"
+#include "AliTrackReference.h"
+#include "TTreeStream.h"
+
+#include "AliAnalysisManager.h"
+#include "AliTRDtrackInfo/AliTRDtrackInfo.h"
+#include "AliTRDtrackingEfficiencyCombined.h"
+
+ClassImp(AliTRDtrackingEfficiencyCombined)
+
+//_____________________________________________________________________________
+AliTRDtrackingEfficiencyCombined::AliTRDtrackingEfficiencyCombined(const char *name):
+  AliAnalysisTask(name, " "),
+  fObjectContainer(0x0),
+  fTrackInfos(0x0),
+  fDebugLevel(0),
+  fDebugStream(0x0)
+{
+  //
+  // Default constructor
+  //
+
+  DefineInput(0, TObjArray::Class());
+  DefineOutput(0, TObjArray::Class());
+}
+
+//_____________________________________________________________________________
+void AliTRDtrackingEfficiencyCombined::ConnectInputData(Option_t *){
+  //
+  // Connect input data
+  //
+
+  fTrackInfos = dynamic_cast<TObjArray *>(GetInputData(0));
+}
+
+//_____________________________________________________________________________
+void AliTRDtrackingEfficiencyCombined::CreateOutputObjects(){
+  //
+  // Create output objects
+  //
+
+  OpenFile(0, "RECREATE");
+
+  const Int_t nbins = 11;
+  Float_t xbins[nbins+1] = {.5, .7, .9, 1.3, 1.7, 2.4, 3.5, 4.5, 5.5, 7., 9., 11.};
+  
+  fObjectContainer = new TObjArray();
+  fObjectContainer->Add(new TProfile("trEffComb", "Combined Tracking Efficiency", nbins, xbins));
+  fObjectContainer->Add(new TProfile("trContComb", "Combined Tracking Contamination", nbins, xbins));
+  if(fDebugLevel) fDebugStream = new TTreeSRedirector("TRDeffDebug.root");
+}
+
+//_____________________________________________________________________________
+void AliTRDtrackingEfficiencyCombined::Exec(Option_t *){
+  //
+  // Do it
+  //
+
+       const Float_t kAlpha = 0.349065850;
+       Int_t naccepted = 0, nrejected = 0, ndoublecounted = 0;
+       Int_t labelsacc[10000];
+       Int_t labelsrej[10000];
+       Float_t momacc[10000];
+       Float_t momrej[10000];
+       if(!AliTracker::GetFieldMap()){
+               AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
+               AliTracker::SetFieldMap(field, kTRUE);
+       }
+       TProfile *efficiency = (TProfile *)fObjectContainer->UncheckedAt(0);
+       TProfile *contamination = (TProfile *)fObjectContainer->UncheckedAt(1);
+       
+       Int_t nTrackInfos = fTrackInfos->GetEntriesFast();
+       Double_t mom = 0;
+       AliTRDtrackV1 *TRDtrack = 0x0;
+       AliTRDtrackInfo *trkInf = 0x0;
+       AliTrackReference *trackRef = 0x0;
+       for(Int_t itinf = 0; itinf < nTrackInfos; itinf++){
+               mom = 0.;
+               trkInf = dynamic_cast<AliTRDtrackInfo *>(fTrackInfos->UncheckedAt(itinf));
+               if(!trkInf) continue;
+               if((TRDtrack = trkInf->GetTRDtrack()) || trkInf->GetNumberOfClustersRefit()){
+                       // check if allready found by the tracker
+                       Bool_t found = kFALSE;
+                       for(Int_t il = 0; il < naccepted; il++){
+                               if(labelsacc[il] == trkInf->GetLabel()) found = kTRUE;
+                       }
+                       if(found){
+                               mom =  trackRef ? trackRef->P() : TRDtrack->P();
+                               contamination->Fill(mom, 1);
+                               ndoublecounted++;
+                               continue;
+                       }
+                       if(trkInf->GetNTrackRefs()){
+                               Int_t iref = 0;
+                               while(!(trackRef = trkInf->GetTrackRef(iref++)));
+                       }
+                       if(!trackRef) printf("Error: Track Reference missing for Track %d\n", TRDtrack->GetLabel());
+                       mom =  trackRef ? trackRef->P() : TRDtrack->P();
+//           Accept track
+                       if(fDebugLevel > 3)printf("Accept track\n");
+                       momacc[naccepted] = mom; 
+                       labelsacc[naccepted++] = trkInf->GetLabel();
+/*                     printf("Reconstructed: event %3d Tracks: MC[%d] ESD[%d] NRefs[%d]\n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), trkInf->GetLabel(),  trkInf->GetTrackId(), trkInf->GetNTrackRefs());*/
+    } else{
+                       if(fDebugLevel>10) printf("Analysing Track References\n");
+      // Check if track is findable
+                       Float_t xmin = 10000.0, xmax = 0.0; 
+                       Float_t ymin = 0.0, ymax = 0.0;
+                       Float_t zmin = 0.0, zmax = 0.0;
+                       Float_t lastx = 0.0, x = 0.0;
+                       Int_t nLayers = 0;
+/*                     trackRef = trkInf->GetTrackRef(0);*/
+/*                     xmin = trackRef->LocalX(); xmax = trackRef->LocalX();
+      ymin = trackRef->LocalY(); ymax = trackRef->LocalY();
+      mom = trackRef->P();*/
+      Int_t *sector = new Int_t[trkInf->GetNTrackRefs()];
+      for(Int_t itr = 0; itr < trkInf->GetNTrackRefs(); itr++){
+        trackRef = trkInf->GetTrackRef(itr);
+                               if(fDebugLevel>10) printf("%d. x[%f], y[%f], z[%f]\n", itr, trackRef->LocalX(), trackRef->LocalY(), trackRef->Z());
+                               x = trackRef->LocalX(); 
+                               if(x < 250. || x > 370.) continue;      // Be Sure that we are inside TRD
+        sector[itr] = Int_t(trackRef->Alpha()/kAlpha);
+        if(x < xmin){
+          xmin = trackRef->LocalX();
+          ymin = trackRef->LocalY();
+          zmin = trackRef->Z();
+          mom = trackRef->P();
+        }
+        if(x > xmax){
+          xmax = trackRef->LocalX();
+          ymax = trackRef->LocalY();
+          zmax = trackRef->Z();
+        }
+                               if(itr > 0){
+                                       Float_t dist = TMath::Abs(x - lastx);
+                                       if(fDebugLevel>10) printf("x = %f, lastx = %f, dist = %f\n", x, lastx, dist);
+                                       if(TMath::Abs(dist - 3.7) < 0.1) nLayers++;     // ref(i+1) has to be larger than ref(i)
+                               }
+                               lastx = x;
+      }
+      // Apply cuts
+      Bool_t findable = kTRUE;
+      if(trkInf->GetNTrackRefs() > 2 && xmax > xmin){
+        if(mom < 0.55) findable = kFALSE;                                                                      // momentum cut at 0.6
+                               Double_t yangle = (ymax -ymin)/(xmax - xmin);
+                               Double_t zangle = (zmax -zmin)/(xmax - xmin);
+                               if(fDebugLevel>10) printf("track: y-Angle = %f, z-Angle = %f\n", yangle, zangle);
+                               if(fDebugLevel>10) printf("nLayers = %d\n", nLayers);
+                               if(TMath::ATan(TMath::Abs(yangle)) > 45.) findable = kFALSE;
+                               if(TMath::ATan(TMath::Abs(zangle)) > 45.) findable = kFALSE;
+                               if(nLayers < 4) findable = kFALSE;
+        if(!trkInf->IsPrimary()) findable = kFALSE;
+        Bool_t samesec = kTRUE;
+        for(Int_t iref = 1; iref < trkInf->GetNTrackRefs(); iref++)
+          if(sector[iref] != sector[0]) samesec = kFALSE;
+        if(!samesec) findable = kFALSE;                // Discard all tracks which are passing more than one sector
+        if(fDebugLevel){
+          Double_t trackAngle = TMath::ATan(yangle);
+          Bool_t primary = trkInf->IsPrimary();
+          (*fDebugStream) << "NotFoundTrack"
+            << "Momentum="     << mom
+            << "trackAngle="<< trackAngle
+            << "NLayers="      << nLayers
+            << "Primary="      << primary
+            << "\n";
+        }
+      }
+      else
+        findable = kFALSE;
+      delete[] sector;
+      if(findable){
+                               momrej[nrejected] = mom;
+                               labelsrej[nrejected++] = trkInf->GetLabel();
+/*                       printf("Not Reconstructed: event %3d Tracks: MC[%d] ESD[%d] NRefs[%d]\n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), trkInf->GetLabel(),  trkInf->GetTrackId(), trkInf->GetNTrackRefs());*/
+      }
+    }
+  }
+       for(Int_t itk = 0; itk < naccepted; itk++){
+               if(fDebugLevel > 2)printf("Accepted MC track: %d\n", labelsacc[itk]);
+               efficiency->Fill(momacc[itk], 1);
+               contamination->Fill(momacc[itk], 0);
+       }
+       Int_t nall = naccepted;
+       for(Int_t imis = 0; imis < nrejected; imis++){
+               Bool_t found = kFALSE;
+               for(Int_t ifound = 0; ifound < naccepted; ifound++){
+                       if(labelsacc[ifound] == labelsrej[imis]){
+                               found = kTRUE;
+                               break;
+                       }
+               }
+               if(!found){
+                       efficiency->Fill(momrej[imis], 0);
+                       contamination->Fill(momrej[imis], 0);
+                       if(fDebugLevel > 2)printf("Rejected MC track: %d\n", labelsrej[imis]);
+                       nall++;
+               }
+       }
+  //if(fDebugLevel>=1)
+  printf("%3d Tracks: MC[%3d] TRD[%3d | %5.2f%%] \n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), nall, naccepted, 1.E2*Float_t(naccepted)/Float_t(nall));
+  printf("%3d Tracks: ALL[%3d] DoubleCounted[%3d | %5.2f%%] \n", (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), nall + ndoublecounted, ndoublecounted, 1.E2*Float_t(ndoublecounted)/Float_t(nall + ndoublecounted));
+
+  PostData(0, fObjectContainer);
+}
+
+//_____________________________________________________________________________
+void AliTRDtrackingEfficiencyCombined::Terminate(Option_t *)
+{
+  //
+  // Termination
+  //
+
+/*  fObjectContainer = dynamic_cast<TObjArray*>(GetOutputData(0));
+  if (!fObjectContainer) {
+    Printf("ERROR: list not available");
+    return;
+  }
+
+  TProfile *hEff = (TProfile*)fObjectContainer->At(0);
+  TProfile *hEffCont = (TProfile*)fObjectContainer->At(1);
+  Printf("Eff[%p] EffCont[%p]\n", (void*)hEff, (void*)hEffCont);
+
+
+  TCanvas *c2 = new TCanvas("c2","",800,400);
+  c2->Divide(2,1);
+
+  c2->cd(1);
+  hEff->DrawCopy("e1");
+  c2->cd(2);
+  hEffCont->DrawCopy("e1");*/
+  if(fDebugStream) delete fDebugStream;
+}
diff --git a/TRD/qaRec/AliTRDtrackingEfficiencyCombined.h b/TRD/qaRec/AliTRDtrackingEfficiencyCombined.h
new file mode 100644 (file)
index 0000000..c603ba1
--- /dev/null
@@ -0,0 +1,44 @@
+#ifndef ALITRDTRACKINGEFFICIENCYCOMBINED_H
+#define ALITRDTRACKINGEFFICIENCYCOMBINED_H
+/* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
+ * See cxx source for full Copyright notice                               */
+
+/* $Id: AliTRDtrackingEfficiencyCombined.h 27496 2008-07-22 08:35:45Z cblume $ */
+
+////////////////////////////////////////////////////////////////////////////
+//                                                                        //
+//  Reconstruction QA                                                     //
+//                                                                        //
+////////////////////////////////////////////////////////////////////////////
+
+#include "AliAnalysisTask.h"
+
+class TObjArray;
+class TTreeSRedirector;
+
+class AliTRDtrackingEfficiencyCombined : public AliAnalysisTask{
+       public:
+               AliTRDtrackingEfficiencyCombined(const char *name = "combined tracking efficiency");
+               ~AliTRDtrackingEfficiencyCombined(){};
+               
+               virtual void ConnectInputData(Option_t *);
+               virtual void CreateOutputObjects();
+               virtual void Exec(Option_t *);
+               virtual void Terminate(Option_t *);
+               
+               void SetDebugLevel(Int_t debugLevel) { fDebugLevel = debugLevel;}
+               int GetDebugLevel() const { return fDebugLevel; }
+               
+       private:
+               AliTRDtrackingEfficiencyCombined(const AliTRDtrackingEfficiencyCombined &);
+               AliTRDtrackingEfficiencyCombined& operator=(const AliTRDtrackingEfficiencyCombined &);
+               
+               TObjArray *fObjectContainer;                   //! Container for output histograms
+               TObjArray *fTrackInfos;                        //! Input Container
+               Int_t fDebugLevel;                             // Debug level
+               TTreeSRedirector *fDebugStream;                // Debug streamer
+               
+               ClassDef(AliTRDtrackingEfficiencyCombined, 1); // Combined tracking efficiency
+};
+               
+#endif
diff --git a/TRD/qaRec/CreateESDChain.C b/TRD/qaRec/CreateESDChain.C
new file mode 100644 (file)
index 0000000..fbeaa6f
--- /dev/null
@@ -0,0 +1,22 @@
+TChain* CreateESDChain(const char* filename = "ESDfiles.txt", Int_t nfiles=-1 )
+{
+  // Create the chain
+  TChain* chain = new TChain("esdTree");
+
+  // Open the input stream
+  ifstream in;
+  in.open(filename);
+
+  // Read the input list of files and add them to the chain
+  TString esdfile;
+  while(in.good() && (nfiles--) ) {
+    in >> esdfile;
+    if (!esdfile.Contains("root")) continue; // protection
+
+    chain->Add(esdfile.Data());
+  }
+
+  in.close();
+
+  return chain;
+}
diff --git a/TRD/qaRec/run.C b/TRD/qaRec/run.C
new file mode 100644 (file)
index 0000000..bef6ada
--- /dev/null
@@ -0,0 +1,84 @@
+Char_t *source = "/d/alice12/mfasel/train"; // To Be Changed
+
+void run(const Char_t *list="list.txt", Int_t nmax=-1) {
+  TStopwatch timer;
+  timer.Start();
+  Char_t *source = Form("%s/TRD/qaRec",gSystem->ExpandPathName("$ALICE_ROOT"));
+
+/*  gSystem->Load("libSTEERBase.so");               
+  gSystem->Load("libESD.so");
+  gSystem->Load("libVMC.so");
+  gSystem->Load("libAOD.so");*/
+  gSystem->Load("libANALYSIS.so");
+  gSystem->Load(Form("%s/AliTRDtrackInfo/libTRDtrackInfo.so",source));
+       
+  gROOT->LoadMacro(Form("%s/AliTRDtrackInfoGen.cxx+", source));
+  gROOT->LoadMacro(Form("%s/AliTRDtrackingEfficiency.cxx+", source));
+  gROOT->LoadMacro(Form("%s/AliTRDtrackingEfficiencyCombined.cxx+", source));
+
+  //____________________________________________//
+
+  TChain* chain = 0x0;
+  gROOT->LoadMacro(Form("%s/CreateESDChain.C", source));
+  chain = CreateESDChain(list,nmax);
+  //chain->SetBranchStatus("*", 0);
+  chain->SetBranchStatus("*FMD*",0);
+  chain->SetBranchStatus("*Calo*",0);
+  chain->SetBranchStatus("Tracks", 1);
+  chain->SetBranchStatus("ESDfriend*",1);
+  chain->Lookup();
+  chain->GetListOfFiles()->Print();
+  printf("\n ----> CHAIN HAS %d ENTRIES <----\n\n", (Int_t)chain->GetEntries());
+  
+  AliLog::SetGlobalLogLevel(AliLog::kError);
+
+  //____________________________________________
+  // Make the analysis manager
+  AliAnalysisManager *mgr = new AliAnalysisManager("TRD Track Info Manager");
+  //mgr->SetSpecialOutputLocation(source); // To Be Changed
+  AliVEventHandler* esdH = new AliESDInputHandler;
+  mgr->SetInputEventHandler(esdH);  
+  AliMCEventHandler *mc = new AliMCEventHandler();
+  mgr->SetMCtruthEventHandler(mc);
+  //mgr->SetDebugLevel(10);
+
+  //____________________________________________
+  // TRD track summary generator
+  AliTRDtrackInfoGen *task1 = new AliTRDtrackInfoGen();
+  task1->SetDebugLevel(1);
+  mgr->AddTask(task1);
+  // Create containers for input/output
+  AliAnalysisDataContainer *cinput1 = mgr->CreateContainer("data", TChain::Class(), AliAnalysisManager::kInputContainer);
+  AliAnalysisDataContainer *coutput1 = mgr->CreateContainer("TrackInfoList", TObjArray::Class(), AliAnalysisManager::kExchangeContainer);
+  mgr->ConnectInput( task1, 0, cinput1);
+  mgr->ConnectOutput(task1, 0, coutput1);
+
+  //____________________________________________
+  // TRD barrel tracking efficiency
+  AliTRDtrackingEfficiency *task2 = new AliTRDtrackingEfficiency();
+  task2->SetDebugLevel(1);
+  mgr->AddTask(task2);
+  //Create containers for input/output
+  AliAnalysisDataContainer *coutput2 = mgr->CreateContainer("TrackingEfficiency", TList::Class(), AliAnalysisManager::kOutputContainer, "TRD.TrackingEfficiency.root");
+  mgr->ConnectInput( task2, 0, coutput1);
+  mgr->ConnectOutput(task2, 0, coutput2);
+
+  //____________________________________________
+  // TRD combined tracking efficiency
+  AliTRDtrackingEfficiencyCombined *task3 = new AliTRDtrackingEfficiencyCombined();
+  task3->SetDebugLevel(0);
+  mgr->AddTask(task3);
+  // Create containers for input/output
+  AliAnalysisDataContainer *coutput3 = mgr->CreateContainer("TrackingEfficiencyCombined", TObjArray::Class(), AliAnalysisManager::kOutputContainer, "TRD.TrackingEfficiencyCombined.root");
+  mgr->ConnectInput( task3, 0, coutput1);
+  mgr->ConnectOutput(task3, 0, coutput3);
+
+
+
+  if (!mgr->InitAnalysis()) return;
+  mgr->PrintStatus();
+  mgr->StartAnalysis("local",chain);
+
+  timer.Stop();
+  timer.Print();
+}