-/**************************************************************************
-* 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: AliTRDinfoGen.cxx 27496 2008-07-22 08:35:45Z cblume $ */
-
-////////////////////////////////////////////////////////////////////////////
-//
-// Tender wagon for TRD performance/calibration train
-//
-// In this wagon the information from
-// - ESD
-// - Friends [if available]
-// - MC [if available]
-// are grouped into AliTRDtrackInfo objects and fed to worker tasks
-//
-// Authors:
-// Markus Fasel <M.Fasel@gsi.de>
-// Alexandru Bercuci <A.Bercuci@gsi.de>
-//
-////////////////////////////////////////////////////////////////////////////
-
-#include <TClonesArray.h>
-#include <TObjArray.h>
-#include <TObject.h>
-#include <TString.h>
-#include <TH1S.h>
-#include <TPad.h>
-#include <TFile.h>
-#include <TTree.h>
-#include <TROOT.h>
-#include <TChain.h>
-#include <TParticle.h>
-
-#include "AliLog.h"
-#include "AliAnalysisManager.h"
-#include "AliGeomManager.h"
-#include "AliCDBManager.h"
-#include "AliCDBEntry.h"
-#include "AliCDBPath.h"
-#include "AliESDEvent.h"
-#include "AliMCEvent.h"
-#include "AliESDInputHandler.h"
-#include "AliMCEventHandler.h"
-
-#include "AliESDfriend.h"
-#include "AliESDfriendTrack.h"
-#include "AliESDHeader.h"
-#include "AliESDtrack.h"
-#include "AliESDv0.h"
-#include "AliESDtrackCuts.h"
-#include "AliMCParticle.h"
-#include "AliPID.h"
-#include "AliStack.h"
-#include "AliTrackReference.h"
-#include "TTreeStream.h"
-
-#include <cstdio>
-#include <climits>
-#include <cstring>
-#include <iostream>
-#include <memory>
-
-#include "AliTRDReconstructor.h"
-#include "AliTRDrecoParam.h"
-#include "AliTRDcalibDB.h"
-#include "AliTRDtrackerV1.h"
-#include "AliTRDgeometry.h"
-#include "AliTRDtrackV1.h"
-#include "AliTRDseedV1.h"
-#include "AliTRDcluster.h"
-#include "AliTRDinfoGen.h"
-#include "AliTRDpwg1Helper.h"
-#include "info/AliTRDtrackInfo.h"
-#include "info/AliTRDeventInfo.h"
-#include "info/AliTRDv0Info.h"
-#include "info/AliTRDeventCuts.h"
-
-ClassImp(AliTRDinfoGen)
-
-const Float_t AliTRDinfoGen::fgkITS = 100.; // to be checked
-const Float_t AliTRDinfoGen::fgkTPC = 290.;
-const Float_t AliTRDinfoGen::fgkTRD = 365.;
-
-const Float_t AliTRDinfoGen::fgkEvVertexZ = 15.;
-const Int_t AliTRDinfoGen::fgkEvVertexN = 1;
-
-const Float_t AliTRDinfoGen::fgkTrkDCAxy = 3.;
-const Float_t AliTRDinfoGen::fgkTrkDCAz = 10.;
-const Int_t AliTRDinfoGen::fgkNclTPC = 70;
-const Float_t AliTRDinfoGen::fgkPt = 0.2;
-const Float_t AliTRDinfoGen::fgkEta = 0.9;
-AliTRDReconstructor* AliTRDinfoGen::fgReconstructor(NULL);
-AliTRDgeometry* AliTRDinfoGen::fgGeo(NULL);
-//____________________________________________________________________
-AliTRDinfoGen::AliTRDinfoGen()
- :AliAnalysisTaskSE()
- ,fEvTrigger(NULL)
- ,fESDev(NULL)
- ,fMCev(NULL)
- ,fEventCut(NULL)
- ,fTrackCut(NULL)
- ,fV0Cut(NULL)
- ,fOCDB("local://$ALICE_ROOT/OCDB")
- ,fTrackInfo(NULL)
- ,fEventInfo(NULL)
- ,fV0Info(NULL)
- ,fTracksBarrel(NULL)
- ,fTracksSA(NULL)
- ,fTracksKink(NULL)
- ,fV0List(NULL)
- ,fContainer(NULL)
- ,fDebugStream(NULL)
-{
- //
- // Default constructor
- //
- SetNameTitle("TRDinfoGen", "MC-REC TRD-track list generator");
-}
-
-//____________________________________________________________________
-AliTRDinfoGen::AliTRDinfoGen(char* name)
- :AliAnalysisTaskSE(name)
- ,fEvTrigger(NULL)
- ,fESDev(NULL)
- ,fMCev(NULL)
- ,fEventCut(NULL)
- ,fTrackCut(NULL)
- ,fV0Cut(NULL)
- ,fOCDB("local://$ALICE_ROOT/OCDB")
- ,fTrackInfo(NULL)
- ,fEventInfo(NULL)
- ,fV0Info(NULL)
- ,fTracksBarrel(NULL)
- ,fTracksSA(NULL)
- ,fTracksKink(NULL)
- ,fV0List(NULL)
- ,fContainer(NULL)
- ,fDebugStream(NULL)
-{
- //
- // Default constructor
- //
- SetTitle("MC-REC TRD-track list generator");
- DefineOutput(AliTRDpwg1Helper::kTracksBarrel, TObjArray::Class());
- DefineOutput(AliTRDpwg1Helper::kTracksSA, TObjArray::Class());
- DefineOutput(AliTRDpwg1Helper::kTracksKink, TObjArray::Class());
- DefineOutput(AliTRDpwg1Helper::kEventInfo, AliTRDeventInfo::Class());
- DefineOutput(AliTRDpwg1Helper::kV0List, TObjArray::Class());
- DefineOutput(AliTRDpwg1Helper::kMonitor, TObjArray::Class()); // histogram list
-}
-
-//____________________________________________________________________
-AliTRDinfoGen::~AliTRDinfoGen()
-{
-// Destructor
- if(fgGeo) delete fgGeo;
- if(fgReconstructor) delete fgReconstructor;
- if(fDebugStream) delete fDebugStream;
- if(fEvTrigger) delete fEvTrigger;
- if(fV0Cut) delete fV0Cut;
- if(fTrackCut) delete fTrackCut;
- if(fEventCut) delete fEventCut;
- if(fTrackInfo) delete fTrackInfo; fTrackInfo = NULL;
- if(fEventInfo) delete fEventInfo; fEventInfo = NULL;
- if(fV0Info) delete fV0Info; fV0Info = NULL;
- if(fTracksBarrel){
- fTracksBarrel->Delete(); delete fTracksBarrel;
- fTracksBarrel = NULL;
- }
- if(fTracksSA){
- fTracksSA->Delete(); delete fTracksSA;
- fTracksSA = NULL;
- }
- if(fTracksKink){
- fTracksKink->Delete(); delete fTracksKink;
- fTracksKink = NULL;
- }
- if(fV0List){
- fV0List->Delete();
- delete fV0List;
- fV0List = NULL;
- }
- if(fContainer){
- fContainer->Delete();
- delete fContainer;
- fContainer = NULL;
- }
-}
-
-//____________________________________________________________________
-Bool_t AliTRDinfoGen::GetRefFigure(Int_t)
-{
-// General graphs for PWG1/TRD train
- if(!gPad){
- AliWarning("Please provide a canvas to draw results.");
- return kFALSE;
- }
- fContainer->At(0)->Draw("bar");
- return kTRUE;
-}
-
-//____________________________________________________________________
-void AliTRDinfoGen::UserCreateOutputObjects()
-{
- //
- // Create Output Containers (TObjectArray containing 1D histograms)
- //
-
- fTrackInfo = new AliTRDtrackInfo();
- fEventInfo = new AliTRDeventInfo();
- fV0Info = new AliTRDv0Info();
- fTracksBarrel = new TObjArray(200); fTracksBarrel->SetOwner(kTRUE);
- fTracksSA = new TObjArray(20); fTracksSA->SetOwner(kTRUE);
- fTracksKink = new TObjArray(20); fTracksKink->SetOwner(kTRUE);
- fV0List = new TObjArray(10); fV0List->SetOwner(kTRUE);
-
- // define general monitor
- fContainer = new TObjArray(1); fContainer->SetOwner(kTRUE);
- TH1 *h=new TH1I("hStat", "Run statistics;Observable;Entries", Int_t(kNObjects), -0.5, Float_t(kNObjects)-0.5);
- TAxis *ax(h->GetXaxis());
- ax->SetBinLabel(Int_t(kTracksESD) + 1, "ESD");
- ax->SetBinLabel(Int_t(kTracksMC) + 1, "MC");
- ax->SetBinLabel(Int_t(kV0) + 1, "V0");
- ax->SetBinLabel(Int_t(kTPC) + 1, "TPC");
- ax->SetBinLabel(Int_t(kTRDin) + 1, "TRDin");
- ax->SetBinLabel(Int_t(kTRDout) + 1, "TRDout");
- ax->SetBinLabel(Int_t(kBarrel) + 1, "Barrel");
- ax->SetBinLabel(Int_t(kBarrelMC) + 1, "BarrelMC");
- ax->SetBinLabel(Int_t(kSA) + 1, "SA");
- ax->SetBinLabel(Int_t(kSAMC) + 1, "SAMC");
- ax->SetBinLabel(Int_t(kKink) + 1, "Kink");
- ax->SetBinLabel(Int_t(kKinkMC) + 1, "KinkMC");
- ax->SetBinLabel(Int_t(kBarrelFriend) + 1, "BFriend");
- ax->SetBinLabel(Int_t(kSAFriend) + 1, "SFriend");
- fContainer->AddAt(h, 0);
- PostData(AliTRDpwg1Helper::kMonitor, fContainer);
-}
-
-//____________________________________________________________________
-Bool_t AliTRDinfoGen::Load(const Char_t *file, const Char_t *dir, const Char_t *name)
-{
-// Load data from performance file
-
- if(!TFile::Open(file)){
- AliWarning(Form("Couldn't open file %s.", file));
- return kFALSE;
- }
- if(dir){
- if(!gFile->cd(dir)){
- AliWarning(Form("Couldn't cd to %s in %s.", dir, file));
- return kFALSE;
- }
- }
- TObjArray *o(NULL);
- const Char_t *tn=(name ? name : GetName());
- if(!(o = (TObjArray*)gDirectory->Get(tn))){
- AliWarning(Form("Missing histogram container %s.", tn));
- return kFALSE;
- }
- fContainer = (TObjArray*)o->Clone(GetName());
- gFile->Close();
- return kTRUE;
-}
-
-//____________________________________________________________________
-void AliTRDinfoGen::UserExec(Option_t *){
- //
- // Run the Analysis
- //
-
- fTracksBarrel->Delete();
- fTracksSA->Delete();
- fTracksKink->Delete();
- fV0List->Delete();
- fEventInfo->Delete("");
-
- fESDev = dynamic_cast<AliESDEvent*>(InputEvent());
- if(!fESDev){
- AliError("Failed retrieving ESD event");
- return;
- }
- // WARNING
- // This part may conflict with other detectors !!
- if(!IsInitOCDB()){
- AliInfo("Initializing OCDB ...");
- // prepare OCDB access
- AliCDBManager* ocdb = AliCDBManager::Instance();
- ocdb->SetDefaultStorage(fOCDB.Data());
- ocdb->SetRun(fESDev->GetRunNumber());
- // create geo manager
- AliCDBEntry* obj(NULL);
- if(!(obj = ocdb->Get(AliCDBPath("GRP", "Geometry", "Data")))){
- AliError("GEOMETRY failed initialization.");
- } else {
- AliGeomManager::SetGeometry((TGeoManager*)obj->GetObject());
- AliGeomManager::GetNalignable("TRD");
- AliGeomManager::ApplyAlignObjsFromCDB("TRD");
- }
- fgGeo = new AliTRDgeometry;
- fgGeo->CreateClusterMatrixArray();
- //init magnetic field
- if(!TGeoGlobalMagField::Instance()->IsLocked() &&
- !fESDev->InitMagneticField()){
- AliError("MAGNETIC FIELD failed initialization.");
- }
-
- // set no of time bins
- AliTRDtrackerV1::SetNTimeBins(AliTRDcalibDB::Instance()->GetNumberOfTimeBinsDCS());
- AliInfo(Form("OCDB : Loc[%s] Run[%d] TB[%d]", fOCDB.Data(), ocdb->GetRun(), AliTRDtrackerV1::GetNTimeBins()));
-
- // set reco param valid for this run/event
- AliInfo(Form("Initializing TRD reco params for EventSpecie[%d]...",
- fESDev->GetEventSpecie()));
- fgReconstructor = new AliTRDReconstructor();
- if(!(obj = ocdb->Get(AliCDBPath("TRD", "Calib", "RecoParam")))){
- AliError("RECO PARAM failed initialization.");
- fgReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());
- } else {
- obj->PrintMetaData();
- TObjArray *recos((TObjArray*)obj->GetObject());
- for(Int_t ireco(0); ireco<recos->GetEntriesFast(); ireco++){
- AliTRDrecoParam *reco((AliTRDrecoParam*)recos->At(ireco));
- Int_t es(reco->GetEventSpecie());
- if(!(es&fESDev->GetEventSpecie())) continue;
- fgReconstructor->SetRecoParam(reco);
- TString s;
- if(es&AliRecoParam::kLowMult) s="LowMult";
- else if(es&AliRecoParam::kHighMult) s="HighMult";
- else if(es&AliRecoParam::kCosmic) s="Cosmic";
- else if(es&AliRecoParam::kCalib) s="Calib";
- else s="Unknown";
- AliInfo(Form("Using reco params for %s", s.Data()));
- break;
- }
- }
- SetInitOCDB();
- }
-
- // link MC if available
- fMCev = MCEvent();
-
- // event selection : trigger cut
- if(UseLocalEvSelection() && fEvTrigger){
- Bool_t kTRIGGERED(kFALSE);
- std::auto_ptr<TObjArray> trig(fEvTrigger->Tokenize(" "));
- for(Int_t itrig=trig->GetEntriesFast(); itrig--;){
- const Char_t *trigClass(((TObjString*)(*trig)[itrig])->GetName());
- if(fESDev->IsTriggerClassFired(trigClass)) {
- AliDebug(2, Form("Ev[%4d] Trigger[%s]", fESDev->GetEventNumberInFile(), trigClass));
- kTRIGGERED = kTRUE;
- break;
- }
- }
- if(!kTRIGGERED){
- AliDebug(2, Form("Reject Ev[%4d] Trigger", fESDev->GetEventNumberInFile()));
- return;
- }
- // select only physical events
- if(fESDev->GetEventType() != 7){
- AliDebug(2, Form("Reject Ev[%4d] EvType[%d]", fESDev->GetEventNumberInFile(), fESDev->GetEventType()));
- return;
- }
- }
-
- // if the required trigger is a collision trigger then apply event vertex cut
- if(UseLocalEvSelection() && IsCollision()){
- const AliESDVertex *vertex = fESDev->GetPrimaryVertex();
- if(TMath::Abs(vertex->GetZv())<1.e-10 ||
- TMath::Abs(vertex->GetZv())>fgkEvVertexZ ||
- vertex->GetNContributors()<fgkEvVertexN) {
- AliDebug(2, Form("Reject Ev[%4d] Vertex Zv[%f] Nv[%d]", fESDev->GetEventNumberInFile(), TMath::Abs(vertex->GetZv()), vertex->GetNContributors()));
- return;
- }
- }
-
- if(fEventCut && !fEventCut->IsSelected(fESDev, IsCollision())) return;
-
- if(!fESDfriend){
- AliError("Failed retrieving ESD friend event");
- return;
- }
- if(HasMCdata() && !fMCev){
- AliError("Failed retrieving MC event");
- return;
- }
-
- new(fEventInfo)AliTRDeventInfo(fESDev->GetHeader(), const_cast<AliESDRun *>(fESDev->GetESDRun()));
-
- Bool_t *trackMap(NULL);
- AliStack * mStack(NULL);
- Int_t nTracksMC = HasMCdata() ? fMCev->GetNumberOfTracks() : 0, nTracksESD = fESDev->GetNumberOfTracks();
- if(HasMCdata()){
- mStack = fMCev->Stack();
- if(!mStack){
- AliError("Failed retrieving MC Stack");
- return;
- }
- trackMap = new Bool_t[nTracksMC];
- memset(trackMap, 0, sizeof(Bool_t) * nTracksMC);
- }
-
- Double32_t dedx[100]; Int_t nSlices(0);
- Int_t nTRDout(0), nTRDin(0), nTPC(0)
- ,nclsTrklt
- ,nBarrel(0), nSA(0), nKink(0)
- ,nBarrelFriend(0), nSAFriend(0)
- ,nBarrelMC(0), nSAMC(0), nKinkMC(0);
- AliESDtrack *esdTrack = NULL;
- AliESDfriendTrack *esdFriendTrack = NULL;
- TObject *calObject = NULL;
- AliTRDtrackV1 *track = NULL;
- AliTRDseedV1 *tracklet = NULL;
- AliTRDcluster *cl = NULL;
-
-
- // LOOP 0 - over ESD v0s
- Float_t bField(fESDev->GetMagneticField());
- AliESDv0 *v0(NULL);
- Int_t v0pid[AliPID::kSPECIES];
- for(Int_t iv0(0); iv0<fESDev->GetNumberOfV0s(); iv0++){
- if(!(v0 = fESDev->GetV0(iv0))) continue;
- // register v0
- if(fV0Cut) new(fV0Info) AliTRDv0Info(*fV0Cut);
- else new(fV0Info) AliTRDv0Info();
- fV0Info->SetMagField(bField);
- fV0Info->SetV0tracks(fESDev->GetTrack(v0->GetPindex()), fESDev->GetTrack(v0->GetNindex()));
- fV0Info->SetV0Info(v0);
- fV0List->Add(new AliTRDv0Info(*fV0Info));// kFOUND=kFALSE;
- }
-
-
- // LOOP 1 - over ESD tracks
- AliTRDv0Info *v0info=NULL;
- for(Int_t itrk = 0; itrk < nTracksESD; itrk++){
- new(fTrackInfo) AliTRDtrackInfo();
- esdTrack = fESDev->GetTrack(itrk);
- AliDebug(3, Form("\n%3d ITS[%d] TPC[%d] TRD[%d]\n", itrk, esdTrack->GetNcls(0), esdTrack->GetNcls(1), esdTrack->GetNcls(2)));
-
- if(esdTrack->GetStatus()&AliESDtrack::kTPCout) nTPC++;
- if(esdTrack->GetStatus()&AliESDtrack::kTRDout) nTRDout++;
- if(esdTrack->GetStatus()&AliESDtrack::kTRDin) nTRDin++;
-
- // look at external track param
- const AliExternalTrackParam *op = esdTrack->GetOuterParam();
- Double_t xyz[3];
- if(op){
- op->GetXYZ(xyz);
- op->Global2LocalPosition(xyz, op->GetAlpha());
- AliDebug(3, Form("op @ X[%7.3f]\n", xyz[0]));
- }
-
- // read MC info
- Int_t fPdg = -1;
- Int_t label = -1; UInt_t alab=UINT_MAX;
- if(HasMCdata()){
- label = esdTrack->GetLabel();
- alab = TMath::Abs(label);
- // register the track
- if(alab < UInt_t(nTracksMC)){
- trackMap[alab] = kTRUE;
- } else {
- AliError(Form("MC label[%d] outside scope for Ev[%d] Trk[%d].", label, (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), itrk));
- continue;
- }
- AliMCParticle *mcParticle = NULL;
- if(!(mcParticle = (AliMCParticle*) fMCev->GetTrack(alab))){
- AliError(Form("MC particle label[%d] missing for Ev[%d] Trk[%d].", label, (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), itrk));
- continue;
- }
- fPdg = mcParticle->Particle()->GetPdgCode();
- Int_t nRefs = mcParticle->GetNumberOfTrackReferences();
- Int_t iref = 0; AliTrackReference *ref = NULL;
- while(iref<nRefs){
- ref = mcParticle->GetTrackReference(iref);
- if(ref->LocalX() > fgkTPC) break;
- iref++;
- }
-
- fTrackInfo->SetMC();
- fTrackInfo->SetPDG(fPdg);
- fTrackInfo->SetPrimary(mcParticle->Particle()->IsPrimary());
- fTrackInfo->SetLabel(label);
- Int_t jref = iref;//, kref = 0;
- while(jref<nRefs){
- ref = mcParticle->GetTrackReference(jref);
- if(ref->LocalX() > fgkTRD) break;
- AliDebug(4, Form(" trackRef[%2d (%2d)] @ %7.3f OK", jref-iref, jref, ref->LocalX()));
- fTrackInfo->AddTrackRef(ref);
- jref++;
- }
- AliDebug(3, Form("NtrackRefs[%d(%d)]", fTrackInfo->GetNTrackRefs(), nRefs));
- }
-
- // copy some relevant info to TRD track info
- fTrackInfo->SetStatus(esdTrack->GetStatus());
- fTrackInfo->SetTrackId(esdTrack->GetID());
- Double_t p[AliPID::kSPECIES]; esdTrack->GetTRDpid(p);
- fTrackInfo->SetESDpid(p);
- fTrackInfo->SetESDpidQuality(esdTrack->GetTRDntrackletsPID());
- if(!nSlices) nSlices = esdTrack->GetNumberOfTRDslices();
- memset(dedx, 0, 100*sizeof(Double32_t));
- Int_t in(0);
- for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++)
- for(Int_t is=0; is<nSlices; is++)
- dedx[in++]=esdTrack->GetTRDslice(il, is);
- for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++) dedx[in++]=esdTrack->GetTRDmomentum(il);
- fTrackInfo->SetSlices(in, dedx);
- fTrackInfo->SetNumberOfClustersRefit(esdTrack->GetNcls(2));
- // some other Informations which we may wish to store in order to find problematic cases
- fTrackInfo->SetKinkIndex(esdTrack->GetKinkIndex(0));
- fTrackInfo->SetTPCncls(static_cast<UShort_t>(esdTrack->GetNcls(1)));
- nclsTrklt = 0;
-
- // set V0pid info
- for(Int_t iv(0); iv<fV0List->GetEntriesFast(); iv++){
- if(!(v0info = (AliTRDv0Info*)fV0List->At(iv))) continue;
- if(!v0info->GetV0Daughter(1) && !v0info->GetV0Daughter(-1)) continue;
- if(!v0info->HasTrack(fTrackInfo)) continue;
- memset(v0pid, 0, AliPID::kSPECIES*sizeof(Int_t));
- fTrackInfo->SetV0();
- for(Int_t is=AliPID::kSPECIES; is--;){v0pid[is] = v0info->GetPID(is, fTrackInfo);}
- fTrackInfo->SetV0pid(v0pid);
- fTrackInfo->SetV0();
- //const AliTRDtrackInfo::AliESDinfo *ei = fTrackInfo->GetESDinfo();
- break;
- }
-
- // 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
- if(!(track = dynamic_cast<AliTRDtrackV1*>(calObject))) break;
- AliDebug(4, Form("TRD track OK"));
- // Set the clusters to unused
- for(Int_t ipl = 0; ipl < AliTRDgeometry::kNlayer; ipl++){
- if(!(tracklet = track->GetTracklet(ipl))) continue;
- tracklet->ResetClusterIter();
- while((cl = tracklet->NextCluster())) cl->Use(0);
- }
- fTrackInfo->SetTrack(track);
- break;
- }
- AliDebug(3, Form("Ntracklets[%d]\n", fTrackInfo->GetNTracklets()));
- } else AliDebug(3, "No ESD friends");
- if(op) fTrackInfo->SetOuterParam(op);
-
- if(DebugLevel() >= 1){
- AliTRDtrackInfo info(*fTrackInfo);
- (*DebugStream()) << "trackInfo"
- << "TrackInfo.=" << &info
- << "\n";
- info.Delete("");
- }
-
- ULong_t status(esdTrack->GetStatus());
- if((status&AliESDtrack::kTPCout)){
- if(!esdTrack->GetKinkIndex(0)){ // Barrel Track Selection
- Bool_t selected(kTRUE);
- if(UseLocalTrkSelection()){
- if(esdTrack->Pt() < fgkPt){
- AliDebug(2, Form("Reject Ev[%4d] Trk[%3d] Pt[%5.2f]", fESDev->GetEventNumberInFile(), itrk, esdTrack->Pt()));
- selected = kFALSE;
- }
- if(selected && TMath::Abs(esdTrack->Eta()) > fgkEta){
- AliDebug(2, Form("Reject Ev[%4d] Trk[%3d] Eta[%5.2f]", fESDev->GetEventNumberInFile(), itrk, TMath::Abs(esdTrack->Eta())));
- selected = kFALSE;
- }
- if(selected && esdTrack->GetTPCNcls() < fgkNclTPC){
- AliDebug(2, Form("Reject Ev[%4d] Trk[%3d] NclTPC[%d]", fESDev->GetEventNumberInFile(), itrk, esdTrack->GetTPCNcls()));
- selected = kFALSE;
- }
- Float_t par[2], cov[3];
- esdTrack->GetImpactParameters(par, cov);
- if(IsCollision()){ // cuts on DCA
- if(selected && TMath::Abs(par[0]) > fgkTrkDCAxy){
- AliDebug(2, Form("Reject Ev[%4d] Trk[%3d] DCAxy[%f]", fESDev->GetEventNumberInFile(), itrk, TMath::Abs(par[0])));
- selected = kFALSE;
- }
- if(selected && TMath::Abs(par[1]) > fgkTrkDCAz){
- AliDebug(2, Form("Reject Ev[%4d] Trk[%3d] DCAz[%f]", fESDev->GetEventNumberInFile(), itrk, TMath::Abs(par[1])));
- selected = kFALSE;
- }
- } else if(selected && fMCev && !fMCev->IsPhysicalPrimary(alab)){;
- AliDebug(2, Form("Reject Ev[%4d] Trk[%3d] Primary", fESDev->GetEventNumberInFile(), itrk));
- selected = kFALSE;
- }
- }
- if(fTrackCut && !fTrackCut->IsSelected(esdTrack)) selected = kFALSE;
- if(selected){
- fTracksBarrel->Add(new AliTRDtrackInfo(*fTrackInfo));
- nBarrel++;
- if(fTrackInfo->GetTrack())
- nBarrelFriend++;
- }
- } else {
- fTracksKink->Add(new AliTRDtrackInfo(*fTrackInfo));
- nKink++;
- }
- } else if((status&AliESDtrack::kTRDout) && !(status&AliESDtrack::kTRDin)){
- fTracksSA->Add(new AliTRDtrackInfo(*fTrackInfo));
- nSA++;
- if(fTrackInfo->GetTrack())
- nSAFriend++;
- }
- fTrackInfo->Delete("");
- }
-
- // LOOP 2 - over MC tracks which are passing TRD where the track is not reconstructed
- if(HasMCdata()){
- AliDebug(10, "Output of the MC track map:");
- for(Int_t itk = 0; itk < nTracksMC; itk++) AliDebug(10, Form("trackMap[%d] = %s", itk, trackMap[itk] == kTRUE ? "TRUE" : "kFALSE"));
-
- for(Int_t itk = 0; itk < nTracksMC; itk++){
- if(trackMap[itk]) continue;
- AliMCParticle *mcParticle = (AliMCParticle*) fMCev->GetTrack(TMath::Abs(itk));
- Int_t fPdg = mcParticle->Particle()->GetPdgCode();
- Int_t nRefs = mcParticle->GetNumberOfTrackReferences();
- Int_t iref = 0; AliTrackReference *ref = NULL;
- Int_t nRefsTRD = 0;
- new(fTrackInfo) AliTRDtrackInfo();
- fTrackInfo->SetMC();
- fTrackInfo->SetPDG(fPdg);
- while(iref<nRefs){ // count TRD TR
- Bool_t kIN(kFALSE);
- ref = mcParticle->GetTrackReference(iref);
- if(ref->LocalX() > fgkTPC && ref->LocalX() < fgkTRD){
- fTrackInfo->AddTrackRef(ref);
- nRefsTRD++;kIN=kTRUE;
- }
- AliDebug(4, Form(" trackRef[%2d] @ x[%7.3f] %s", iref, ref->LocalX(), kIN?"IN":"OUT"));
- 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(DebugLevel() >= 1){
- AliTRDtrackInfo info(*fTrackInfo);
- (*DebugStream()) << "trackInfo"
- << "TrackInfo.=" << &info
- << "\n";
- info.Delete("");
- }
- AliDebug(3, Form("Add MC track @ label[%d] nTRDrefs[%d].", itk, nRefsTRD));
- // check where the track starts
- ref = mcParticle->GetTrackReference(0);
- if(ref->LocalX() < fgkITS){
- fTracksBarrel->Add(new AliTRDtrackInfo(*fTrackInfo));
- nBarrelMC++;
- } else if(ref->LocalX() < fgkTPC) {
- fTracksKink->Add(new AliTRDtrackInfo(*fTrackInfo));
- nKinkMC++;
- } else if(nRefsTRD>6){
- fTracksSA->Add(new AliTRDtrackInfo(*fTrackInfo));
- nSAMC++;
- }
- fTrackInfo->Delete("");
- }
- delete[] trackMap;
- }
- AliDebug(1, Form(
- "\nEv[%3d] Tracks: ESD[%d] MC[%d] V0[%d]\n"
- " TPCout[%d] TRDin[%d] TRDout[%d]\n"
- " Barrel[%3d+%3d=%3d] SA[%2d+%2d=%2d] Kink[%2d+%2d=%2d]"
- ,(Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), nTracksESD, nTracksMC, fV0List->GetEntries()
- , nTPC, nTRDin, nTRDout
- ,nBarrel, nBarrelMC, fTracksBarrel->GetEntries()
- ,nSA, nSAMC, fTracksSA->GetEntries()
- ,nKink, nKinkMC, fTracksKink->GetEntries()
- ));
- // save track statistics
- TH1 *h((TH1S*)fContainer->At(0));
- h->Fill(Float_t(kTracksESD), nTracksESD);
- h->Fill(Float_t(kTracksMC), nTracksMC);
- h->Fill(Float_t(kV0), fV0List->GetEntries());
- h->Fill(Float_t(kTPC), nTPC);
- h->Fill(Float_t(kTRDin), nTRDin);
- h->Fill(Float_t(kTRDout), nTRDout);
- h->Fill(Float_t(kBarrel), nBarrel);
- h->Fill(Float_t(kBarrelMC), nBarrelMC);
- h->Fill(Float_t(kSA), nSA);
- h->Fill(Float_t(kSAMC), nSAMC);
- h->Fill(Float_t(kKink), nKink);
- h->Fill(Float_t(kKinkMC), nKinkMC);
- h->Fill(Float_t(kBarrelFriend), nBarrelFriend);
- h->Fill(Float_t(kSAFriend), nSAFriend);
-
- PostData(AliTRDpwg1Helper::kTracksBarrel, fTracksBarrel);
- PostData(AliTRDpwg1Helper::kTracksSA, fTracksSA);
- PostData(AliTRDpwg1Helper::kTracksKink, fTracksKink);
- PostData(AliTRDpwg1Helper::kEventInfo, fEventInfo);
- PostData(AliTRDpwg1Helper::kV0List, fV0List);
-}
-
-//____________________________________________________________________
-void AliTRDinfoGen::SetLocalV0Selection(const AliTRDv0Info *v0)
-{
-// Set V0 cuts from outside
-
- if(!fV0Cut) fV0Cut = new AliTRDv0Info(*v0);
- else new(fV0Cut) AliTRDv0Info(*v0);
-}
-
-//____________________________________________________________________
-void AliTRDinfoGen::SetTrigger(const Char_t *trigger)
-{
- if(!fEvTrigger) fEvTrigger = new TString(trigger);
- else (*fEvTrigger) = trigger;
-}
-
-//____________________________________________________________________
-TTreeSRedirector* AliTRDinfoGen::DebugStream()
-{
-// Manage debug stream for task
- if(!fDebugStream){
- TDirectory *savedir = gDirectory;
- fDebugStream = new TTreeSRedirector("TRD.DebugInfoGen.root");
- savedir->cd();
- }
- return fDebugStream;
-}
-
-
-
+/**************************************************************************\r
+* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+* *\r
+* Author: The ALICE Off-line Project. *\r
+* Contributors are mentioned in the code where appropriate. *\r
+* *\r
+* Permission to use, copy, modify and distribute this software and its *\r
+* documentation strictly for non-commercial purposes is hereby granted *\r
+* without fee, provided that the above copyright notice appears in all *\r
+* copies and that both the copyright notice and this permission notice *\r
+* appear in the supporting documentation. The authors make no claims *\r
+* about the suitability of this software for any purpose. It is *\r
+* provided "as is" without express or implied warranty. *\r
+**************************************************************************/\r
+\r
+/* $Id: AliTRDinfoGen.cxx 27496 2008-07-22 08:35:45Z cblume $ */\r
+\r
+////////////////////////////////////////////////////////////////////////////\r
+//\r
+// Tender wagon for TRD performance/calibration train\r
+//\r
+// In this wagon the information from\r
+// - ESD\r
+// - Friends [if available]\r
+// - MC [if available]\r
+// are grouped into AliTRDtrackInfo objects and fed to worker tasks\r
+//\r
+// Authors:\r
+// Markus Fasel <M.Fasel@gsi.de>\r
+// Alexandru Bercuci <A.Bercuci@gsi.de>\r
+//\r
+////////////////////////////////////////////////////////////////////////////\r
+\r
+#include <TClonesArray.h>\r
+#include <TObjArray.h>\r
+#include <TObject.h>\r
+#include <TString.h>\r
+#include <TVectorT.h>\r
+#include <TH1.h>\r
+#include <TH2.h>\r
+#include <TCanvas.h>\r
+#include <TPad.h>\r
+#include <TFile.h>\r
+#include <TTree.h>\r
+#include <TROOT.h>\r
+#include <TChain.h>\r
+#include <TParticle.h>\r
+\r
+#include "AliLog.h"\r
+#include "AliAnalysisManager.h"\r
+#include "AliGeomManager.h"\r
+#include "AliCDBManager.h"\r
+#include "AliCDBStorage.h"\r
+#include "AliCDBEntry.h"\r
+#include "AliCDBPath.h"\r
+#include "AliESDEvent.h"\r
+#include "AliMCEvent.h"\r
+#include "AliESDInputHandler.h"\r
+#include "AliMCEventHandler.h"\r
+\r
+#include "AliESDfriend.h"\r
+#include "AliESDfriendTrack.h"\r
+#include "AliESDHeader.h"\r
+#include "AliESDRun.h"\r
+#include "AliESDtrack.h"\r
+#include "AliESDv0.h"\r
+#include "AliESDtrackCuts.h"\r
+#include "AliMCParticle.h"\r
+#include "AliMultiplicity.h"\r
+#include "AliPID.h"\r
+#include "AliStack.h"\r
+#include "AliTrackReference.h"\r
+#include "TTreeStream.h"\r
+\r
+#include <cstdio>\r
+#include <climits>\r
+#include <cstring>\r
+#include <iostream>\r
+#include <memory>\r
+\r
+#include "AliTRDReconstructor.h"\r
+#include "AliTRDrecoParam.h"\r
+#include "AliTRDcalibDB.h"\r
+#include "AliTRDtrackerV1.h"\r
+#include "AliTRDpadPlane.h"\r
+#include "AliTRDgeometry.h"\r
+#include "AliTRDtrackV1.h"\r
+#include "AliTRDseedV1.h"\r
+#include "AliTRDcluster.h"\r
+#include "AliTRDinfoGen.h"\r
+#include "AliTRDpwg1Helper.h"\r
+#include "info/AliTRDtrackInfo.h"\r
+#include "info/AliTRDeventInfo.h"\r
+#include "info/AliTRDv0Info.h"\r
+#include "info/AliTRDeventCuts.h"\r
+\r
+ClassImp(AliTRDinfoGen)\r
+\r
+const Float_t AliTRDinfoGen::fgkITS = 100.; // to be checked\r
+const Float_t AliTRDinfoGen::fgkTPC = 290.;\r
+const Float_t AliTRDinfoGen::fgkTRD = 365.;\r
+\r
+const Float_t AliTRDinfoGen::fgkTrkDCAxy = 3.;\r
+const Float_t AliTRDinfoGen::fgkTrkDCAz = 10.;\r
+const Int_t AliTRDinfoGen::fgkNclTPC = 70;\r
+const Float_t AliTRDinfoGen::fgkPt = 0.2;\r
+const Float_t AliTRDinfoGen::fgkEta = 0.9;\r
+AliTRDReconstructor* AliTRDinfoGen::fgReconstructor(NULL);\r
+AliTRDgeometry* AliTRDinfoGen::fgGeo(NULL);\r
+//____________________________________________________________________\r
+AliTRDinfoGen::AliTRDinfoGen()\r
+ :AliAnalysisTaskSE()\r
+ ,fESDev(NULL)\r
+ ,fMCev(NULL)\r
+ ,fEventCut(NULL)\r
+ ,fTrackCut(NULL)\r
+ ,fV0Cut(NULL)\r
+ ,fOCDB("local://$ALICE_ROOT/OCDB")\r
+ ,fTrackInfo(NULL)\r
+ ,fEventInfo(NULL)\r
+ ,fV0Info(NULL)\r
+ ,fTracksBarrel(NULL)\r
+ ,fTracksSA(NULL)\r
+ ,fTracksKink(NULL)\r
+ ,fV0List(NULL)\r
+ ,fContainer(NULL)\r
+ ,fRecos(NULL)\r
+ ,fDebugStream(NULL)\r
+{\r
+ //\r
+ // Default constructor\r
+ //\r
+ SetNameTitle("TRDinfoGen", "MC-REC TRD-track list generator");\r
+}\r
+\r
+//____________________________________________________________________\r
+AliTRDinfoGen::AliTRDinfoGen(char* name)\r
+ :AliAnalysisTaskSE(name)\r
+ ,fESDev(NULL)\r
+ ,fMCev(NULL)\r
+ ,fEventCut(NULL)\r
+ ,fTrackCut(NULL)\r
+ ,fV0Cut(NULL)\r
+ ,fOCDB("local://$ALICE_ROOT/OCDB")\r
+ ,fTrackInfo(NULL)\r
+ ,fEventInfo(NULL)\r
+ ,fV0Info(NULL)\r
+ ,fTracksBarrel(NULL)\r
+ ,fTracksSA(NULL)\r
+ ,fTracksKink(NULL)\r
+ ,fV0List(NULL)\r
+ ,fContainer(NULL)\r
+ ,fRecos(NULL)\r
+ ,fDebugStream(NULL)\r
+{\r
+ //\r
+ // Default constructor\r
+ //\r
+ SetTitle("MC-REC TRD-track list generator");\r
+ DefineOutput(AliTRDpwg1Helper::kTracksBarrel, TObjArray::Class());\r
+ DefineOutput(AliTRDpwg1Helper::kTracksSA, TObjArray::Class());\r
+ DefineOutput(AliTRDpwg1Helper::kTracksKink, TObjArray::Class());\r
+ DefineOutput(AliTRDpwg1Helper::kEventInfo, AliTRDeventInfo::Class());\r
+ DefineOutput(AliTRDpwg1Helper::kV0List, TObjArray::Class());\r
+ DefineOutput(AliTRDpwg1Helper::kMonitor, TObjArray::Class()); // histogram list\r
+}\r
+\r
+//____________________________________________________________________\r
+AliTRDinfoGen::~AliTRDinfoGen()\r
+{\r
+// Destructor\r
+ if(fgGeo) delete fgGeo;\r
+ if(fgReconstructor) delete fgReconstructor;\r
+ if(fDebugStream) delete fDebugStream;\r
+ if(fV0Cut) delete fV0Cut;\r
+ if(fTrackCut) delete fTrackCut;\r
+ if(fEventCut) delete fEventCut;\r
+ if(fTrackInfo) delete fTrackInfo; fTrackInfo = NULL;\r
+ if(fEventInfo) delete fEventInfo; fEventInfo = NULL;\r
+ if(fV0Info) delete fV0Info; fV0Info = NULL;\r
+ if(fTracksBarrel){ \r
+ fTracksBarrel->Delete(); delete fTracksBarrel;\r
+ fTracksBarrel = NULL;\r
+ }\r
+ if(fTracksSA){ \r
+ fTracksSA->Delete(); delete fTracksSA;\r
+ fTracksSA = NULL;\r
+ }\r
+ if(fTracksKink){ \r
+ fTracksKink->Delete(); delete fTracksKink;\r
+ fTracksKink = NULL;\r
+ }\r
+ if(fV0List){ \r
+ fV0List->Delete(); \r
+ delete fV0List;\r
+ fV0List = NULL;\r
+ }\r
+ if(fContainer && !(AliAnalysisManager::GetAnalysisManager() && AliAnalysisManager::GetAnalysisManager()->IsProofMode())){\r
+ fContainer->Delete(); \r
+ delete fContainer;\r
+ fContainer = NULL;\r
+ }\r
+}\r
+\r
+//____________________________________________________________________\r
+Bool_t AliTRDinfoGen::GetRefFigure(Int_t)\r
+{\r
+// General graphs for PWG1/TRD train\r
+ if(!gPad){\r
+ AliWarning("Please provide a canvas to draw results.");\r
+ return kFALSE;\r
+ }\r
+ fContainer->At(kStatTrk)->Draw("bar");\r
+ return kTRUE;\r
+}\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::UserCreateOutputObjects()\r
+{ \r
+ //\r
+ // Create Output Containers (TObjectArray containing 1D histograms)\r
+ //\r
+ \r
+ fTrackInfo = new AliTRDtrackInfo();\r
+ fEventInfo = new AliTRDeventInfo();\r
+ fV0Info = new AliTRDv0Info();\r
+ fTracksBarrel = new TObjArray(200); fTracksBarrel->SetOwner(kTRUE);\r
+ fTracksSA = new TObjArray(20); fTracksSA->SetOwner(kTRUE);\r
+ fTracksKink = new TObjArray(20); fTracksKink->SetOwner(kTRUE);\r
+ fV0List = new TObjArray(10); fV0List->SetOwner(kTRUE);\r
+\r
+ // define general monitor\r
+ fContainer = new TObjArray(kNclasses); fContainer->SetOwner(kTRUE);\r
+ TH1 *h=new TH1I("hStat", "Run statistics;Observable;Entries", Int_t(kNObjects), -0.5, Float_t(kNObjects)-0.5);\r
+ TAxis *ax(h->GetXaxis());\r
+ ax->SetBinLabel(Int_t(kTracksESD) + 1, "ESD");\r
+ ax->SetBinLabel(Int_t(kTracksMC) + 1, "MC");\r
+ ax->SetBinLabel(Int_t(kV0) + 1, "V0");\r
+ ax->SetBinLabel(Int_t(kTPC) + 1, "TPC");\r
+ ax->SetBinLabel(Int_t(kTRDin) + 1, "TRDin");\r
+ ax->SetBinLabel(Int_t(kTRDout) + 1, "TRDout");\r
+ ax->SetBinLabel(Int_t(kBarrel) + 1, "Barrel");\r
+ ax->SetBinLabel(Int_t(kBarrelMC) + 1, "BarrelMC");\r
+ ax->SetBinLabel(Int_t(kSA) + 1, "SA");\r
+ ax->SetBinLabel(Int_t(kSAMC) + 1, "SAMC");\r
+ ax->SetBinLabel(Int_t(kKink) + 1, "Kink");\r
+ ax->SetBinLabel(Int_t(kKinkMC) + 1, "KinkMC");\r
+ ax->SetBinLabel(Int_t(kBarrelFriend) + 1, "BFriend");\r
+ ax->SetBinLabel(Int_t(kSAFriend) + 1, "SFriend");\r
+ fContainer->AddAt(h, kStatTrk);\r
+ h=new TH1I("hEv", "Run statistics;Event Class;Entries", 4, -0.5, 3.5);\r
+ ax = h->GetXaxis();\r
+ ax->SetBinLabel(1, "Low");\r
+ ax->SetBinLabel(2, "High");\r
+ ax->SetBinLabel(3, "Cosmic");\r
+ ax->SetBinLabel(4, "Calib");\r
+ fContainer->AddAt(h, kEvType);\r
+ TH2I* h2=new TH2I("hBC", "Bunch statistics;Fill Bunch;TOF BC;Entries", 3500, -0.5, 3499.5, 31, -10.5, 20.5);\r
+ fContainer->AddAt(h2, kBC);\r
+ h=new TH1I("hTriggers", "Triggers statistics;;Entries", 21, -0.5, 20.5);\r
+ fContainer->AddAt(h, kTrigger);\r
+ TObjArray *chmb = new TObjArray(AliTRDgeometry::kNdet);\r
+ chmb->SetName("Chambers"); chmb->SetOwner();\r
+ fContainer->AddAt(chmb, kChmb);\r
+ PostData(AliTRDpwg1Helper::kMonitor, fContainer);\r
+}\r
+\r
+//____________________________________________________________________\r
+Bool_t AliTRDinfoGen::Load(const Char_t *file, const Char_t *dir, const Char_t *name)\r
+{\r
+// Load data from performance file\r
+\r
+ if(!TFile::Open(file)){\r
+ AliWarning(Form("Couldn't open file %s.", file));\r
+ return kFALSE;\r
+ }\r
+ if(dir){\r
+ if(!gFile->cd(dir)){\r
+ AliWarning(Form("Couldn't cd to %s in %s.", dir, file));\r
+ return kFALSE;\r
+ }\r
+ }\r
+ TObjArray *o(NULL);\r
+ const Char_t *tn=(name ? name : GetName());\r
+ if(!(o = (TObjArray*)gDirectory->Get(tn))){\r
+ AliWarning(Form("Missing histogram container %s.", tn));\r
+ return kFALSE;\r
+ }\r
+ fContainer = (TObjArray*)o->Clone(GetName());\r
+ gFile->Close();\r
+ return kTRUE;\r
+}\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::UserExec(Option_t *){\r
+ //\r
+ // Run the Analysis\r
+ //\r
+\r
+ fTracksBarrel->Delete();\r
+ fTracksSA->Delete();\r
+ fTracksKink->Delete();\r
+ fV0List->Delete();\r
+ fEventInfo->Delete("");\r
+\r
+ fESDev = dynamic_cast<AliESDEvent*>(InputEvent());\r
+ if(!fESDev){\r
+ AliError("Failed retrieving ESD event");\r
+ return;\r
+ }\r
+ // WARNING\r
+ // This part may conflict with other detectors !!\r
+ if(!IsInitOCDB()){\r
+ AliCDBEntry* obj(NULL);\r
+ AliCDBManager* ocdb = AliCDBManager::Instance();\r
+ if(ocdb->IsDefaultStorageSet()){\r
+ AliInfo("OCDB :: initialized externally.");\r
+ } else {\r
+ AliInfo("OCDB :: initializing locally ...");\r
+ // prepare OCDB access\r
+ ocdb->SetDefaultStorage(fOCDB.Data());\r
+ ocdb->SetRun(fESDev->GetRunNumber());\r
+ // create geo manager\r
+ if(!(obj = ocdb->Get(AliCDBPath("GRP", "Geometry", "Data")))){\r
+ AliError("GEOMETRY failed initialization.");\r
+ } else {\r
+ AliGeomManager::SetGeometry((TGeoManager*)obj->GetObject());\r
+ AliGeomManager::GetNalignable("TRD");\r
+ AliGeomManager::ApplyAlignObjsFromCDB("TRD");\r
+ }\r
+ //init magnetic field\r
+ if(!TGeoGlobalMagField::Instance()->IsLocked() &&\r
+ !fESDev->InitMagneticField()){\r
+ AliError("MAGNETIC FIELD failed initialization.");\r
+ }\r
+ // set no of time bins\r
+ AliTRDtrackerV1::SetNTimeBins(AliTRDcalibDB::Instance()->GetNumberOfTimeBinsDCS());\r
+ }\r
+ AliInfo(Form("OCDB : Loc[%s] Run[%d] TB[%d]", ocdb->GetDefaultStorage()->GetURI().Data(), ocdb->GetRun(), AliTRDtrackerV1::GetNTimeBins()));\r
+\r
+ // load misalignment\r
+ fgGeo = new AliTRDgeometry;\r
+ fgGeo->CreateClusterMatrixArray();\r
+ MakeChambers();\r
+ // load reco param list from OCDB\r
+ AliInfo("Initializing TRD reco params ...");\r
+ fgReconstructor = new AliTRDReconstructor();\r
+ if(!(obj = ocdb->Get(AliCDBPath("TRD", "Calib", "RecoParam")))){\r
+ AliError("RECO PARAM failed initialization.");\r
+ } else {\r
+ obj->PrintMetaData();\r
+ fRecos = (TObjArray*)obj->GetObject();\r
+ }\r
+ SetInitOCDB();\r
+ }\r
+\r
+ AliDebug(2, "+++++++++++++++++++++++++++++++++++++++");\r
+ AliDebug(2, Form("Analyse Ev[%d]", fESDev->GetEventNumberInFile()));\r
+\r
+ // set reco param valid for this event\r
+ TH1 *h = (TH1I*)fContainer->At(kEvType);\r
+ if(!fRecos){\r
+ fgReconstructor->SetRecoParam(AliTRDrecoParam::GetLowFluxParam());\r
+ h->Fill(0);\r
+ } else {\r
+ for(Int_t ireco(0); ireco<fRecos->GetEntriesFast(); ireco++){\r
+ AliTRDrecoParam *reco((AliTRDrecoParam*)fRecos->At(ireco));\r
+ Int_t es(reco->GetEventSpecie());\r
+ if(!(es&fESDev->GetEventSpecie())) continue;\r
+ fgReconstructor->SetRecoParam(reco);\r
+ if(AliLog::GetDebugLevel("PWG1/TRD", "AliTRDinfoGen")>5) reco->Dump();\r
+ TString s;\r
+ if(es&AliRecoParam::kLowMult){ s="LowMult"; h->Fill(0);}\r
+ else if(es&AliRecoParam::kHighMult){ s="HighMult"; h->Fill(1);}\r
+ else if(es&AliRecoParam::kCosmic){ s="Cosmic"; h->Fill(2);}\r
+ else if(es&AliRecoParam::kCalib){ s="Calib"; h->Fill(3);}\r
+ else s="Unknown";\r
+ AliDebug(2, Form(" Using reco param \"%s\" for event %d.", s.Data(), fESDev->GetEventNumberInFile()));\r
+ break;\r
+ }\r
+ }\r
+\r
+ // link MC if available\r
+ fMCev = MCEvent();\r
+ \r
+ // trigger monitor\r
+ h = (TH1I*)fContainer->At(kTrigger);\r
+ TAxis *ax(h->GetXaxis());\r
+ TObjArray *evTriggers = fESDev->GetFiredTriggerClasses().Tokenize(" ");\r
+ for(Int_t iet(evTriggers->GetEntriesFast()); iet--;){\r
+ Int_t ix(1);\r
+ for(; ix<=ax->GetNbins(); ix++){\r
+ if(!Int_t(h->GetBinContent(ix))){\r
+ ax->SetBinLabel(ix, (*evTriggers)[iet]->GetName());\r
+ break;\r
+ }\r
+ if(strcmp((*evTriggers)[iet]->GetName(), ax->GetBinLabel(ix))==0) break;\r
+ }\r
+ h->AddBinContent(ix);\r
+ }\r
+\r
+ // event selection based on vertex cuts and trigger\r
+ if(UseLocalEvSelection() && !fEventCut->IsSelected(fESDev, IsCollision())) return;\r
+\r
+ if(!fESDfriend){\r
+ AliError("Failed retrieving ESD friend event");\r
+ return;\r
+ }\r
+ if(HasMCdata() && !fMCev){\r
+ AliError("Failed retrieving MC event");\r
+ return;\r
+ }\r
+\r
+ new(fEventInfo)AliTRDeventInfo(fESDev->GetHeader(), const_cast<AliESDRun *>(fESDev->GetESDRun()));\r
+ // Determine centrality\r
+ // Author: Ionut Arsene <I.C.Arsene@gsi.de>\r
+ Int_t centralityBin = -1;\r
+ AliDebug(2, Form(" Beam Type: %s", fESDev->GetESDRun()->GetBeamType()));\r
+ TString beamtype = fESDev->GetESDRun()->GetBeamType();\r
+ if(beamtype.Contains("Pb-Pb") || beamtype.Contains("A-A")){\r
+ centralityBin = 4;\r
+ const AliMultiplicity *mult = fESDev->GetMultiplicity();\r
+ Double_t zdcNeutronEnergy = fESDev->GetZDCN1Energy()+fESDev->GetZDCN2Energy();\r
+ Double_t itsNTracklets = mult->GetNumberOfTracklets();\r
+ Double_t centralitySlopes[6] = {0.0, 4.0, 8.0, 20.0, 50.0, 1000000.};\r
+ AliDebug(1, Form("zdcNeutronEnergy: %f, itsNTracklets: %f\n", zdcNeutronEnergy, itsNTracklets));\r
+ for(Int_t iCent=1; iCent<=5; ++iCent) {\r
+ if(zdcNeutronEnergy>centralitySlopes[iCent-1]*itsNTracklets && zdcNeutronEnergy<centralitySlopes[iCent]*itsNTracklets)\r
+ centralityBin=iCent - 1;\r
+ }\r
+ AliDebug(2, Form(" Centrality Class: %d", centralityBin));\r
+ }\r
+ fEventInfo->SetCentrality(centralityBin);\r
+ UShort_t evBC(fESDev->GetBunchCrossNumber());\r
+\r
+ Bool_t *trackMap(NULL);\r
+ AliStack * mStack(NULL);\r
+ Int_t nTracksMC = HasMCdata() ? fMCev->GetNumberOfTracks() : 0, nTracksESD = fESDev->GetNumberOfTracks();\r
+ if(HasMCdata()){\r
+ mStack = fMCev->Stack();\r
+ if(!mStack){\r
+ AliError("Failed retrieving MC Stack");\r
+ return;\r
+ }\r
+ trackMap = new Bool_t[nTracksMC];\r
+ memset(trackMap, 0, sizeof(Bool_t) * nTracksMC);\r
+ }\r
+ \r
+ Double32_t dedx[100]; Int_t nSlices(0);\r
+ Int_t nTRDout(0), nTRDin(0), nTPC(0)\r
+ ,nclsTrklt\r
+ ,nBarrel(0), nSA(0), nKink(0)\r
+ ,nBarrelFriend(0), nSAFriend(0)\r
+ ,nBarrelMC(0), nSAMC(0), nKinkMC(0);\r
+ AliESDtrack *esdTrack = NULL;\r
+ AliESDfriendTrack *esdFriendTrack = NULL;\r
+ TObject *calObject = NULL;\r
+ AliTRDtrackV1 *track = NULL;\r
+ AliTRDseedV1 *tracklet = NULL;\r
+ AliTRDcluster *cl = NULL;\r
+\r
+\r
+ // LOOP 0 - over ESD v0s\r
+ Float_t bField(fESDev->GetMagneticField());\r
+ AliESDv0 *v0(NULL);\r
+ Int_t v0pid[AliPID::kSPECIES];\r
+ for(Int_t iv0(0); iv0<fESDev->GetNumberOfV0s(); iv0++){\r
+ if(!(v0 = fESDev->GetV0(iv0))) continue;\r
+ // register v0\r
+ if(fV0Cut) new(fV0Info) AliTRDv0Info(*fV0Cut);\r
+ else new(fV0Info) AliTRDv0Info();\r
+ fV0Info->SetMagField(bField);\r
+ fV0Info->SetV0tracks(fESDev->GetTrack(v0->GetPindex()), fESDev->GetTrack(v0->GetNindex()));\r
+ fV0Info->SetV0Info(v0);\r
+ fV0List->Add(new AliTRDv0Info(*fV0Info));// kFOUND=kFALSE;\r
+ }\r
+\r
+\r
+ // LOOP 1 - over ESD tracks\r
+ AliTRDv0Info *v0info=NULL;\r
+ for(Int_t itrk = 0; itrk < nTracksESD; itrk++){\r
+ new(fTrackInfo) AliTRDtrackInfo();\r
+ esdTrack = fESDev->GetTrack(itrk);\r
+ AliDebug(3, Form("\n%3d ITS[%d] TPC[%d] TRD[%d] TOF-BC[%d]\n", itrk, esdTrack->GetNcls(0), esdTrack->GetNcls(1), esdTrack->GetNcls(2), esdTrack->GetTOFBunchCrossing()));\r
+\r
+ if(esdTrack->GetStatus()&AliESDtrack::kTPCout) nTPC++;\r
+ if(esdTrack->GetStatus()&AliESDtrack::kTRDout) nTRDout++;\r
+ if(esdTrack->GetStatus()&AliESDtrack::kTRDin) nTRDin++;\r
+\r
+ // look at external track param\r
+ const AliExternalTrackParam *op = esdTrack->GetOuterParam();\r
+ Double_t xyz[3];\r
+ if(op){\r
+ op->GetXYZ(xyz);\r
+ op->Global2LocalPosition(xyz, op->GetAlpha());\r
+ AliDebug(3, Form("op @ X[%7.3f]\n", xyz[0]));\r
+ }\r
+\r
+ // read MC info\r
+ Int_t fPdg = -1;\r
+ Int_t label = -1; UInt_t alab=UINT_MAX;\r
+ if(HasMCdata()){\r
+ label = esdTrack->GetLabel(); \r
+ alab = TMath::Abs(label);\r
+ // register the track\r
+ if(alab < UInt_t(nTracksMC)){ \r
+ trackMap[alab] = kTRUE; \r
+ } else { \r
+ AliError(Form("MC label[%d] outside scope for Ev[%d] Trk[%d].", label, (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), itrk));\r
+ continue; \r
+ }\r
+ AliMCParticle *mcParticle = NULL; \r
+ if(!(mcParticle = (AliMCParticle*) fMCev->GetTrack(alab))){\r
+ AliError(Form("MC particle label[%d] missing for Ev[%d] Trk[%d].", label, (Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), itrk));\r
+ continue;\r
+ }\r
+ fPdg = mcParticle->Particle()->GetPdgCode();\r
+ Int_t nRefs = mcParticle->GetNumberOfTrackReferences();\r
+ Int_t iref = 0; AliTrackReference *ref = NULL; \r
+ while(iref<nRefs){\r
+ ref = mcParticle->GetTrackReference(iref);\r
+ if(ref->LocalX() > fgkTPC) break;\r
+ iref++;\r
+ }\r
+\r
+ fTrackInfo->SetMC();\r
+ fTrackInfo->SetPDG(fPdg);\r
+ fTrackInfo->SetPrimary(mcParticle->Particle()->IsPrimary());\r
+ fTrackInfo->SetLabel(label);\r
+ Int_t jref = iref;//, kref = 0;\r
+ while(jref<nRefs){\r
+ ref = mcParticle->GetTrackReference(jref);\r
+ if(ref->LocalX() > fgkTRD) break;\r
+ AliDebug(4, Form(" trackRef[%2d (%2d)] @ %7.3f OK", jref-iref, jref, ref->LocalX()));\r
+ fTrackInfo->AddTrackRef(ref);\r
+ jref++;\r
+ }\r
+ AliDebug(3, Form("NtrackRefs[%d(%d)]", fTrackInfo->GetNTrackRefs(), nRefs));\r
+ }\r
+\r
+ // copy some relevant info to TRD track info\r
+ fTrackInfo->SetStatus(esdTrack->GetStatus());\r
+ fTrackInfo->SetTrackId(esdTrack->GetID());\r
+ Double_t p[AliPID::kSPECIES]; esdTrack->GetTRDpid(p);\r
+ fTrackInfo->SetESDpid(p);\r
+ fTrackInfo->SetESDpidQuality(esdTrack->GetTRDntrackletsPID());\r
+ if(!nSlices) nSlices = esdTrack->GetNumberOfTRDslices();\r
+ memset(dedx, 0, 100*sizeof(Double32_t));\r
+ Int_t in(0);\r
+ for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++)\r
+ for(Int_t is=0; is<nSlices; is++) \r
+ dedx[in++]=esdTrack->GetTRDslice(il, is);\r
+ for(Int_t il=0; il<AliTRDgeometry::kNlayer; il++) dedx[in++]=esdTrack->GetTRDmomentum(il);\r
+ fTrackInfo->SetSlices(in, dedx);\r
+ fTrackInfo->SetNumberOfClustersRefit(esdTrack->GetNcls(2));\r
+ // some other Informations which we may wish to store in order to find problematic cases\r
+ fTrackInfo->SetKinkIndex(esdTrack->GetKinkIndex(0));\r
+ fTrackInfo->SetTPCncls(static_cast<UShort_t>(esdTrack->GetNcls(1)));\r
+ fTrackInfo->SetTOFbc(esdTrack->GetTOFBunchCrossing()==AliVTrack::kTOFBCNA?0:esdTrack->GetTOFBunchCrossing());\r
+ nclsTrklt = 0;\r
+ \r
+ // set V0pid info\r
+ for(Int_t iv(0); iv<fV0List->GetEntriesFast(); iv++){\r
+ if(!(v0info = (AliTRDv0Info*)fV0List->At(iv))) continue;\r
+ if(!v0info->GetV0Daughter(1) && !v0info->GetV0Daughter(-1)) continue;\r
+ if(!v0info->HasTrack(fTrackInfo)) continue;\r
+ memset(v0pid, 0, AliPID::kSPECIES*sizeof(Int_t));\r
+ fTrackInfo->SetV0();\r
+ for(Int_t is=AliPID::kSPECIES; is--;){v0pid[is] = v0info->GetPID(is, fTrackInfo);}\r
+ fTrackInfo->SetV0pid(v0pid);\r
+ fTrackInfo->SetV0();\r
+ //const AliTRDtrackInfo::AliESDinfo *ei = fTrackInfo->GetESDinfo();\r
+ break;\r
+ }\r
+\r
+ // read REC info\r
+ esdFriendTrack = (fESDfriend->GetNumberOfTracks() > itrk) ? fESDfriend->GetTrack(itrk): NULL;\r
+\r
+ if(esdFriendTrack){\r
+ fTrackInfo->SetTPCoutParam(esdFriendTrack->GetTPCOut());\r
+ Int_t icalib = 0;\r
+ while((calObject = esdFriendTrack->GetCalibObject(icalib++))){\r
+ if(strcmp(calObject->IsA()->GetName(),"AliTRDtrackV1") != 0) continue; // Look for the TRDtrack\r
+ if(!(track = dynamic_cast<AliTRDtrackV1*>(calObject))) break;\r
+ AliDebug(4, Form("TRD track OK"));\r
+ // Set the clusters to unused\r
+ for(Int_t ipl = 0; ipl < AliTRDgeometry::kNlayer; ipl++){\r
+ if(!(tracklet = track->GetTracklet(ipl))) continue;\r
+ tracklet->ResetClusterIter();\r
+ while((cl = tracklet->NextCluster())) cl->Use(0);\r
+ }\r
+ fTrackInfo->SetTrack(track);\r
+ ((TH2I*)fContainer->At(kBC))->Fill(evBC, esdTrack->GetTOFBunchCrossing());\r
+ break;\r
+ }\r
+ AliDebug(3, Form("Ntracklets[%d]\n", fTrackInfo->GetNTracklets()));\r
+ } else AliDebug(3, "No ESD friends");\r
+ if(op) fTrackInfo->SetOuterParam(op);\r
+\r
+ if(DebugLevel() >= 1){\r
+ AliTRDtrackInfo info(*fTrackInfo);\r
+ (*DebugStream()) << "trackInfo"\r
+ << "TrackInfo.=" << &info\r
+ << "\n";\r
+ info.Delete("");\r
+ }\r
+\r
+ ULong_t status(esdTrack->GetStatus());\r
+ if((status&AliESDtrack::kTPCout)){\r
+ if(!esdTrack->GetKinkIndex(0)){ // Barrel Track Selection\r
+ Bool_t selected(kTRUE);\r
+ if(UseLocalTrkSelection()){\r
+ if(esdTrack->Pt() < fgkPt){ \r
+ AliDebug(3, Form("Reject Trk[%3d] Ev[%4d] Pt[%5.2f]", itrk, fESDev->GetEventNumberInFile(), esdTrack->Pt()));\r
+ selected = kFALSE;\r
+ }\r
+ if(selected && TMath::Abs(esdTrack->Eta()) > fgkEta){\r
+ AliDebug(3, Form("Reject Trk[%3d] Ev[%4d] Eta[%5.2f]", itrk, fESDev->GetEventNumberInFile(), TMath::Abs(esdTrack->Eta())));\r
+ selected = kFALSE;\r
+ }\r
+ if(selected && esdTrack->GetTPCNcls() < fgkNclTPC){ \r
+ AliDebug(3, Form("Reject Trk[%3d] Ev[%4d] NclTPC[%d]", itrk, fESDev->GetEventNumberInFile(), esdTrack->GetTPCNcls()));\r
+ selected = kFALSE;\r
+ }\r
+ Float_t par[2], cov[3];\r
+ esdTrack->GetImpactParameters(par, cov);\r
+ if(IsCollision()){ // cuts on DCA\r
+ if(selected && TMath::Abs(par[0]) > fgkTrkDCAxy){ \r
+ AliDebug(3, Form("Reject Trk[%3d] Ev[%4d] DCAxy[%f]", itrk, fESDev->GetEventNumberInFile(), TMath::Abs(par[0])));\r
+ selected = kFALSE;\r
+ }\r
+ if(selected && TMath::Abs(par[1]) > fgkTrkDCAz){ \r
+ AliDebug(3, Form("Reject Trk[%3d] Ev[%4d] DCAz[%f]", itrk, fESDev->GetEventNumberInFile(), TMath::Abs(par[1])));\r
+ selected = kFALSE;\r
+ }\r
+ } else if(selected && fMCev && !fMCev->IsPhysicalPrimary(alab)){;\r
+ AliDebug(3, Form("Reject Trk[%3d] Ev[%4d] Primary", itrk, fESDev->GetEventNumberInFile()));\r
+ selected = kFALSE;\r
+ }\r
+ }\r
+ if(fTrackCut && !fTrackCut->IsSelected(esdTrack)) selected = kFALSE;\r
+ if(selected){ \r
+ fTracksBarrel->Add(new AliTRDtrackInfo(*fTrackInfo));\r
+ nBarrel++;\r
+ if(fTrackInfo->GetTrack()) \r
+ nBarrelFriend++;\r
+ }\r
+ } else {\r
+ fTracksKink->Add(new AliTRDtrackInfo(*fTrackInfo));\r
+ nKink++;\r
+ }\r
+ } else if((status&AliESDtrack::kTRDout) && !(status&AliESDtrack::kTRDin)){ \r
+ fTracksSA->Add(new AliTRDtrackInfo(*fTrackInfo));\r
+ nSA++;\r
+ if(fTrackInfo->GetTrack()) \r
+ nSAFriend++;\r
+ }\r
+ fTrackInfo->Delete("");\r
+ }\r
+\r
+ // LOOP 2 - over MC tracks which are passing TRD where the track is not reconstructed\r
+ if(HasMCdata()){\r
+ AliDebug(10, "Output of the MC track map:");\r
+ for(Int_t itk = 0; itk < nTracksMC; itk++) AliDebug(10, Form("trackMap[%d] = %s", itk, trackMap[itk] == kTRUE ? "TRUE" : "kFALSE"));\r
+ \r
+ for(Int_t itk = 0; itk < nTracksMC; itk++){\r
+ if(trackMap[itk]) continue;\r
+ AliMCParticle *mcParticle = (AliMCParticle*) fMCev->GetTrack(TMath::Abs(itk));\r
+ Int_t fPdg = mcParticle->Particle()->GetPdgCode();\r
+ Int_t nRefs = mcParticle->GetNumberOfTrackReferences();\r
+ Int_t iref = 0; AliTrackReference *ref = NULL; \r
+ Int_t nRefsTRD = 0;\r
+ new(fTrackInfo) AliTRDtrackInfo();\r
+ fTrackInfo->SetMC();\r
+ fTrackInfo->SetPDG(fPdg);\r
+ while(iref<nRefs){ // count TRD TR\r
+ Bool_t kIN(kFALSE);\r
+ ref = mcParticle->GetTrackReference(iref);\r
+ if(ref->LocalX() > fgkTPC && ref->LocalX() < fgkTRD){\r
+ fTrackInfo->AddTrackRef(ref);\r
+ nRefsTRD++;kIN=kTRUE;\r
+ }\r
+ AliDebug(4, Form(" trackRef[%2d] @ x[%7.3f] %s", iref, ref->LocalX(), kIN?"IN":"OUT"));\r
+ iref++;\r
+ }\r
+ if(!nRefsTRD){\r
+ // In this stage we at least require 1 hit inside TRD. What will be done with this tracks is a task for the \r
+ // analysis job\r
+ fTrackInfo->Delete("");\r
+ continue;\r
+ }\r
+ fTrackInfo->SetPrimary(mcParticle->Particle()->IsPrimary());\r
+ fTrackInfo->SetLabel(itk);\r
+ if(DebugLevel() >= 1){\r
+ AliTRDtrackInfo info(*fTrackInfo);\r
+ (*DebugStream()) << "trackInfo"\r
+ << "TrackInfo.=" << &info\r
+ << "\n";\r
+ info.Delete("");\r
+ }\r
+ AliDebug(3, Form("Add MC track @ label[%d] nTRDrefs[%d].", itk, nRefsTRD));\r
+ // check where the track starts\r
+ ref = mcParticle->GetTrackReference(0);\r
+ if(ref->LocalX() < fgkITS){ \r
+ fTracksBarrel->Add(new AliTRDtrackInfo(*fTrackInfo));\r
+ nBarrelMC++;\r
+ } else if(ref->LocalX() < fgkTPC) {\r
+ fTracksKink->Add(new AliTRDtrackInfo(*fTrackInfo));\r
+ nKinkMC++;\r
+ } else if(nRefsTRD>6){\r
+ fTracksSA->Add(new AliTRDtrackInfo(*fTrackInfo));\r
+ nSAMC++;\r
+ }\r
+ fTrackInfo->Delete("");\r
+ }\r
+ delete[] trackMap;\r
+ }\r
+ AliDebug(1, Form(\r
+ "\nEv[%3d] Tracks: ESD[%d] MC[%d] V0[%d]\n"\r
+ " TPCout[%d] TRDin[%d] TRDout[%d]\n"\r
+ " Barrel[%3d+%3d=%3d] SA[%2d+%2d=%2d] Kink[%2d+%2d=%2d]"\r
+ ,(Int_t)AliAnalysisManager::GetAnalysisManager()->GetCurrentEntry(), nTracksESD, nTracksMC, fV0List->GetEntries()\r
+ , nTPC, nTRDin, nTRDout\r
+ ,nBarrel, nBarrelMC, fTracksBarrel->GetEntries()\r
+ ,nSA, nSAMC, fTracksSA->GetEntries()\r
+ ,nKink, nKinkMC, fTracksKink->GetEntries()\r
+ ));\r
+ // save track statistics\r
+ h = (TH1I*)fContainer->At(kStatTrk);\r
+ h->Fill(Float_t(kTracksESD), nTracksESD);\r
+ h->Fill(Float_t(kTracksMC), nTracksMC);\r
+ h->Fill(Float_t(kV0), fV0List->GetEntries());\r
+ h->Fill(Float_t(kTPC), nTPC);\r
+ h->Fill(Float_t(kTRDin), nTRDin);\r
+ h->Fill(Float_t(kTRDout), nTRDout);\r
+ h->Fill(Float_t(kBarrel), nBarrel);\r
+ h->Fill(Float_t(kBarrelMC), nBarrelMC);\r
+ h->Fill(Float_t(kSA), nSA);\r
+ h->Fill(Float_t(kSAMC), nSAMC);\r
+ h->Fill(Float_t(kKink), nKink);\r
+ h->Fill(Float_t(kKinkMC), nKinkMC);\r
+ h->Fill(Float_t(kBarrelFriend), nBarrelFriend);\r
+ h->Fill(Float_t(kSAFriend), nSAFriend);\r
+\r
+ PostData(AliTRDpwg1Helper::kTracksBarrel, fTracksBarrel);\r
+ PostData(AliTRDpwg1Helper::kTracksSA, fTracksSA);\r
+ PostData(AliTRDpwg1Helper::kTracksKink, fTracksKink);\r
+ PostData(AliTRDpwg1Helper::kEventInfo, fEventInfo);\r
+ PostData(AliTRDpwg1Helper::kV0List, fV0List);\r
+}\r
+\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::MakeChambers()\r
+{\r
+// Build chamber position and status\r
+ if(!fContainer){\r
+ AliError("Missing container");\r
+ return;\r
+ }\r
+ AliTRDcalibDB *calib(AliTRDcalibDB::Instance());\r
+ if(!calib){\r
+ AliError("No access to calibration data");\r
+ return;\r
+ }\r
+\r
+ TObjArray *chmb((TObjArray*)fContainer->At(kChmb));\r
+ Double_t alpha(0.), cs(-2.), sn(0.);\r
+ TVectorF pos(5);\r
+ for(Int_t isec(0); isec<AliTRDgeometry::kNsector; isec++){\r
+ alpha = (0.5+isec)*AliTRDgeometry::GetAlpha();\r
+ cs = TMath::Cos(alpha);\r
+ sn = TMath::Sin(alpha);\r
+\r
+ for(Int_t istk(0); istk<AliTRDgeometry::kNstack; istk++){\r
+ for(Int_t ilyr(0); ilyr<AliTRDgeometry::kNlayer; ilyr++){\r
+ Int_t idet(AliTRDgeometry::GetDetector(ilyr, istk, isec));\r
+ TGeoHMatrix *matrix(fgGeo->GetClusterMatrix(idet));\r
+ if(!matrix){\r
+ AliDebug(2, Form("Missing matrix for %03d [%02d_%d_%d]", idet, isec, istk, ilyr));\r
+ continue;\r
+ }\r
+ AliDebug(2, Form("Read info for %03d [%02d_%d_%d]", idet, isec, istk, ilyr));\r
+ AliTRDpadPlane *pp(fgGeo->GetPadPlane(ilyr, istk));\r
+ Double_t zm(0.5 * (pp->GetRow0() + pp->GetRowEnd())),\r
+ loc0[] = {AliTRDgeometry::AnodePos(), pp->GetCol0(), zm-pp->GetRow0()},\r
+ loc1[] = {AliTRDgeometry::AnodePos(), pp->GetColEnd(), zm-pp->GetRowEnd()},\r
+ glb[3];\r
+ matrix->LocalToMaster(loc0, glb);\r
+ Float_t phi = TMath::ATan2(glb[0]*sn + glb[1]*cs, glb[0]*cs - glb[1]*sn),\r
+ tgl = glb[2]/glb[0]/TMath::Sqrt(1.+glb[1]*glb[1]/glb[0]/glb[0]),\r
+ eta = -TMath::Log(TMath::Tan(0.5 * (0.5*TMath::Pi() - TMath::ATan(tgl))));\r
+ pos[0] = eta; pos[1] = phi;\r
+ matrix->LocalToMaster(loc1, glb);\r
+ phi = TMath::ATan2(glb[0]*sn + glb[1]*cs, glb[0]*cs - glb[1]*sn);\r
+ tgl = glb[2]/glb[0]/TMath::Sqrt(1.+glb[1]*glb[1]/glb[0]/glb[0]);\r
+ eta = -TMath::Log(TMath::Tan(0.5 * (0.5*TMath::Pi() - TMath::ATan(tgl))));\r
+ pos[2] = eta; pos[3] = phi;\r
+ pos[4] = calib->IsChamberGood(idet)?0.:1.;\r
+ chmb->AddAt(new TVectorF(pos), idet);\r
+ }\r
+ }\r
+ }\r
+}\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::MakeSummary()\r
+{\r
+// Build summary plots\r
+ if(!fContainer){\r
+ AliError("Missing results");\r
+ return;\r
+ }\r
+ TH1 *h1(NULL); TVirtualPad *p(NULL); TCanvas *cOut(NULL);\r
+\r
+ const Int_t nx(2048), ny(750);\r
+ cOut = new TCanvas(GetName(), "Run Statistics", nx, ny);\r
+ cOut->Divide(3,1, 1.e-5, 1.e-5);\r
+ //=========\r
+ p=cOut->cd(1);p->SetRightMargin(0.025);p->SetTopMargin(0.01);p->SetLogy();\r
+ h1 = (TH1*)fContainer->At(kStatTrk);\r
+ h1->SetBarOffset(0.06); h1->SetBarWidth(0.88); h1->SetFillColor(3);\r
+ h1->Draw("bar1");\r
+ //=========\r
+ p=cOut->cd(2);p->SetRightMargin(0.025);p->SetTopMargin(0.01);p->SetLogy();\r
+ h1 = (TH1*)fContainer->At(kEvType);\r
+ h1->SetBarOffset(0.04); h1->SetBarWidth(0.92);h1->SetFillColor(6);\r
+ h1->Draw("bar1");\r
+ //=========\r
+ p=cOut->cd(3);p->SetRightMargin(0.025);p->SetTopMargin(0.01);p->SetLogz();\r
+ TH2 *h2((TH2*)fContainer->At(kBC));\r
+ h1 = h2->ProjectionX("_px");\r
+ Int_t n(0); Int_t bins[3500];\r
+ for(Int_t ib(1); ib<=3500; ib++){\r
+ if(h1->GetBinContent(ib) < 1) continue;\r
+ bins[n++] = ib;\r
+ }\r
+ delete h1;\r
+\r
+ TAxis *ay(h2->GetYaxis());\r
+ TH2 *hs = new TH2I("hBC_Summary", Form("%s;%s;%s", h2->GetTitle(), h2->GetXaxis()->GetTitle(), h2->GetYaxis()->GetTitle()),\r
+ n, -0.5, n-0.5, ay->GetNbins(), ay->GetXmin(), ay->GetXmax());\r
+ hs->SetLineColor(kBlack);hs->SetLineWidth(1);\r
+ hs->SetMarkerColor(kRed);\r
+ TAxis *ax(hs->GetXaxis());\r
+ for(Int_t ib(0); ib<n; ib++){\r
+ ax->SetBinLabel(ib+1, Form("%d", bins[ib]));\r
+ for(Int_t iy(1); iy<=ay->GetNbins(); iy++){\r
+ hs->SetBinContent(ib+1, iy, h2->GetBinContent(bins[ib], iy));\r
+ }\r
+ }\r
+ hs->Draw("textbox");\r
+ cOut->SaveAs(Form("%s.gif", cOut->GetName()));\r
+}\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::SetLocalEvSelection(const AliTRDeventCuts &ec)\r
+{\r
+// Set event cuts from outside\r
+ if(!fEventCut) fEventCut = new AliTRDeventCuts(ec);\r
+ else new(fEventCut) AliTRDeventCuts(ec);\r
+ fEventCut->Print();\r
+}\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::SetLocalV0Selection(const AliTRDv0Info &v0)\r
+{\r
+// Set V0 cuts from outside\r
+\r
+ if(!fV0Cut) fV0Cut = new AliTRDv0Info(v0);\r
+ else new(fV0Cut) AliTRDv0Info(v0);\r
+ fV0Cut->Print();\r
+}\r
+\r
+//____________________________________________________________________\r
+TTreeSRedirector* AliTRDinfoGen::DebugStream()\r
+{\r
+// Manage debug stream for task\r
+ if(!fDebugStream){\r
+ TDirectory *savedir = gDirectory;\r
+ fDebugStream = new TTreeSRedirector("TRD.DebugInfoGen.root");\r
+ savedir->cd();\r
+ }\r
+ return fDebugStream;\r
+}\r
+\r
+//____________________________________________________________________\r
+void AliTRDinfoGen::Terminate(Option_t* /*option*/)\r
+{\r
+// Process run information\r
+ AliInfo("");\r
+ if(!(fContainer = dynamic_cast<TObjArray *>(GetOutputData(AliTRDpwg1Helper::kMonitor)))) return;\r
+ AliInfo(Form("fContainer(%p)", (void*)fContainer));\r
+ if(UseLocalEvSelection()){\r
+ TH1 *h1 = (TH1*)fContainer->At(kTrigger); TAxis *ax(h1->GetXaxis());\r
+ AliInfo(Form("h1(%p)", (void*)h1));\r
+ for(Int_t ix(1); ix<=ax->GetNbins(); ix++){\r
+ if(fEventCut->CheckTrigger(ax->GetBinLabel(ix))) ax->SetBinLabel(ix, Form("#color[2]{%s}", ax->GetBinLabel(ix)));\r
+ }\r
+ }\r
+}\r
+\r
git://git.uio.no / u / mrichter / AliRoot.git / blobdiff