update to reflect latest changes
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 11 Apr 2011 22:10:57 +0000 (22:10 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Mon, 11 Apr 2011 22:10:57 +0000 (22:10 +0000)
PWG4/CaloCalib/AliAnalysisTaskEMCALPi0PbPb.cxx
PWG4/CaloCalib/AliAnalysisTaskEMCALPi0PbPb.h
PWG4/PWG4CaloCalibLinkDef.h

index 9cc61c4..908bb95 100644 (file)
@@ -73,8 +73,13 @@ AliAnalysisTaskEMCALPi0PbPb::AliAnalysisTaskEMCALPi0PbPb(const char *name)
     fAodClusters(0),
     fAodCells(0),
     fPtRanges(0),
-    fNtuple(0),
     fSelTracks(0),
+    fNtuple(0),
+    fHeader(0),
+    fPrimVert(0),
+    fSpdVert(0),
+    fTpcVert(0),
+    fClusters(0),
     fHCuts(0x0),
     fHVertexZ(0x0),
     fHVertexZ2(0x0),
@@ -153,11 +158,22 @@ void AliAnalysisTaskEMCALPi0PbPb::UserCreateOutputObjects()
     TFile *f = OpenFile(1);
     if (f) {
       f->SetCompressionLevel(2);
-      fNtuple = new TNtuple(Form("nt%.0fto%.0f",fCentFrom,fCentTo),"nt",
-                            "run:evt:l0:tcls:cent:pt:eta:phi:e:emax:n:n1:idmax:nsm:db:disp:mn:ms:ecc:sig:tkdz:tkdr:tkep:tkiso:ceiso");
+      fNtuple = new TTree(Form("tree%.0fto%.0f",fCentFrom,fCentTo), "StandaloneTree");
       fNtuple->SetDirectory(f);
       fNtuple->SetAutoFlush(-1024*1024*1024);
       fNtuple->SetAutoSave(-1024*1024*1024);
+      fHeader = new AliStaHeader;
+      fNtuple->Branch("header", &fHeader, 16*1024, 99);
+      fPrimVert = new AliStaVertex;
+      fNtuple->Branch("primv", &fPrimVert, 16*1024, 99);
+      fSpdVert = new AliStaVertex;
+      fNtuple->Branch("spdv", &fSpdVert, 16*1024, 99);
+      fTpcVert = new AliStaVertex;
+      fNtuple->Branch("tpcv", &fTpcVert, 16*1024, 99);
+      if (TClass::GetClass("AliStaCluster"))
+        TClass::GetClass("AliStaCluster")->IgnoreTObjectStreamer();
+      fClusters = new TClonesArray("AliStaCluster",1000);
+      fNtuple->Branch("clusters", &fClusters, 8*16*1024, 99);
     }  
   }
 
@@ -554,6 +570,7 @@ void AliAnalysisTaskEMCALPi0PbPb::UserExec(Option_t *)
   FillClusHists();
   FillPionHists();
   FillOtherHists();
+  FillNtuple();
 
   PostData(1, fOutput);
 }      
@@ -966,6 +983,7 @@ void AliAnalysisTaskEMCALPi0PbPb::FillClusHists()
     clus->GetMomentum(clusterVec,vertex);
     Double_t maxAxis = 0, minAxis = 0;
     GetSigma(clus,maxAxis,minAxis);
+    clus->SetTOF(maxAxis);     // store sigma in TOF
     Double_t clusterEcc = 0;
     if (maxAxis > 0)
       clusterEcc = TMath::Sqrt(1.0 - minAxis*minAxis/(maxAxis*maxAxis));
@@ -976,59 +994,121 @@ void AliAnalysisTaskEMCALPi0PbPb::FillClusHists()
     fHClustEnergySigma->Fill(clus->E()*maxAxis,clus->E());
     fHClustSigmaSigma->Fill(max(clus->GetM02(),clus->GetM20()),clus->E()*maxAxis);
     fHClustNCellEnergyRatio->Fill(clus->GetNCells(),GetMaxCellEnergy(clus)/clus->E());
-    if (fNtuple) {
-      if (clus->E()<fMinE)
-        continue;
-      Float_t vals[25];
-      TString trgclasses;
-      vals[0]  = InputEvent()->GetRunNumber();
-      vals[1]  = (((UInt_t)InputEvent()->GetOrbitNumber()  << 12) | (UInt_t)InputEvent()->GetBunchCrossNumber()); 
-      if (vals[1]<=0) 
-        vals[1] = fNEvs;
-      AliESDHeader *h = dynamic_cast<AliESDHeader*>(InputEvent()->GetHeader());
-      if (h) {
-        vals[2]  = h->GetL0TriggerInputs();
-        trgclasses = fEsdEv->GetFiredTriggerClasses();
-      }
-      else {
-        AliAODHeader *h2 = dynamic_cast<AliAODHeader*>(InputEvent()->GetHeader());
-        if (h2) {
-          vals[2]  = h2->GetL0TriggerInputs();
-          trgclasses = h2->GetFiredTriggerClasses();
-        } else 
-          vals[2] = 0;
-      }
-      vals[3]  = 0;
+  }
+}
 
-      for (Int_t j = 0; j<fTrClassNamesArr->GetEntries(); ++j) {
-        const char *name = fTrClassNamesArr->At(j)->GetName();
-        if (trgclasses.Contains(name))
-          vals[3] += TMath::Power(2,j);
-      }
-      vals[4]  = InputEvent()->GetCentrality()->GetCentralityPercentileUnchecked(fCentVar);
-      vals[5]  = clusterVec.Pt();
-      vals[6]  = clusterVec.Eta();
-      vals[7]  = clusterVec.Phi();
-      vals[8]  = clusterVec.E();
-      vals[9]  = GetMaxCellEnergy(clus);
-      vals[10] = clus->GetNCells();
-      vals[11] = GetNCells(clus,0.100);
-      vals[12] = clus->GetCellAbsId(0);
-      vals[13] = fGeom->GetSuperModuleNumber(clus->GetCellAbsId(0));
-      vals[14] = clus->GetDistanceToBadChannel();
-      vals[15] = clus->GetDispersion();
-      vals[16] = clus->GetM20();
-      vals[17] = clus->GetM02();
-      vals[18] = clusterEcc;
-      vals[19] = maxAxis;
-      vals[20] = fClusProps[i].fTrDz; 
-      vals[21] = fClusProps[i].fTrDr;
-      vals[22] = fClusProps[i].fTrEp;
-      vals[23] = fClusProps[i].fTrIso;
-      vals[24] = fClusProps[i].fCellIso;
-      fNtuple->Fill(vals);
-    }
+//________________________________________________________________________
+void AliAnalysisTaskEMCALPi0PbPb::FillNtuple()
+{
+  // Fill ntuple.
+
+  if (!fNtuple)
+    return;
+
+  AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager();
+  if (fAodEv) {
+    fHeader->fRun            = fAodEv->GetRunNumber();
+    fHeader->fOrbit          = fAodEv->GetHeader()->GetOrbitNumber(); 
+    fHeader->fPeriod         = fAodEv->GetHeader()->GetPeriodNumber();
+    fHeader->fBx             = fAodEv->GetHeader()->GetBunchCrossNumber();
+    fHeader->fL0             = fAodEv->GetHeader()->GetL0TriggerInputs();
+    fHeader->fL1             = fAodEv->GetHeader()->GetL1TriggerInputs();
+    fHeader->fL2             = fAodEv->GetHeader()->GetL2TriggerInputs();
+    fHeader->fTrClassMask    = fAodEv->GetHeader()->GetTriggerMask();
+    fHeader->fTrCluster      = fAodEv->GetHeader()->GetTriggerCluster();
+    fHeader->fOffTriggers    = fAodEv->GetHeader()->GetOfflineTrigger();
+    fHeader->fFiredTriggers  = fAodEv->GetHeader()->GetFiredTriggerClasses();
+  } else {
+    fHeader->fRun            = fEsdEv->GetRunNumber();
+    fHeader->fOrbit          = fEsdEv->GetHeader()->GetOrbitNumber(); 
+    fHeader->fPeriod         = fEsdEv->GetESDRun()->GetPeriodNumber();
+    fHeader->fBx             = fEsdEv->GetHeader()->GetBunchCrossNumber();
+    fHeader->fL0             = fEsdEv->GetHeader()->GetL0TriggerInputs();
+    fHeader->fL1             = fEsdEv->GetHeader()->GetL1TriggerInputs();
+    fHeader->fL2             = fEsdEv->GetHeader()->GetL2TriggerInputs();
+    fHeader->fTrClassMask    = fEsdEv->GetHeader()->GetTriggerMask();
+    fHeader->fTrCluster      = fEsdEv->GetHeader()->GetTriggerCluster();
+    fHeader->fOffTriggers    = ((AliInputEventHandler*)(am->GetInputEventHandler()))->IsEventSelected();
+    fHeader->fFiredTriggers  = fEsdEv->GetFiredTriggerClasses();
+  }
+  AliCentrality *cent = InputEvent()->GetCentrality();
+  fHeader->fV0Cent    = cent->GetCentralityPercentileUnchecked("V0M");
+  fHeader->fCl1Cent   = cent->GetCentralityPercentileUnchecked("CL1");
+  fHeader->fTrCent    = cent->GetCentralityPercentileUnchecked("TRK");
+  fHeader->fCqual     = cent->GetQuality();
+  Double_t val = 0;
+  TString trgclasses(fHeader->fFiredTriggers);
+  for (Int_t j = 0; j<fTrClassNamesArr->GetEntries(); ++j) {
+    const char *name = fTrClassNamesArr->At(j)->GetName();
+    if (trgclasses.Contains(name))
+      val += TMath::Power(2,j);
+  }
+  fHeader->fTcls = (UInt_t)val;
+
+  if (fAodEv) { 
+    am->LoadBranch("vertices");
+    AliAODVertex *pv = fAodEv->GetPrimaryVertex();
+    FillVertex(fPrimVert, pv);
+    AliAODVertex *sv = fAodEv->GetPrimaryVertexSPD();
+    FillVertex(fSpdVert, sv);
+  } else {
+    am->LoadBranch("PrimaryVertex.");
+    const AliESDVertex *pv = fEsdEv->GetPrimaryVertexTracks();
+    FillVertex(fPrimVert, pv);
+    am->LoadBranch("SPDVertex.");
+    const AliESDVertex *sv = fEsdEv->GetPrimaryVertexSPD();
+    FillVertex(fSpdVert, sv);
+    am->LoadBranch("TPCVertex.");
+    const AliESDVertex *tv = fEsdEv->GetPrimaryVertexTPC();
+    FillVertex(fTpcVert, tv);
   }
+
+  TObjArray *clusters = fEsdClusters;
+  if (!clusters)
+    clusters = fAodClusters;
+  if (!clusters)
+    return;
+
+  fClusters->Clear();
+  Int_t nclus = clusters->GetEntries();
+  for(Int_t i=0,ncl=0; i<nclus; ++i) {
+    AliVCluster *clus = static_cast<AliVCluster*>(clusters->At(i));
+    if (!clus)
+      continue;
+    if (!clus->IsEMCAL()) 
+      continue;
+    if (clus->E()<fMinE)
+      continue;
+
+    AliStaCluster *cl = static_cast<AliStaCluster*>(fClusters->New(ncl++));
+    cl->fE        = clus->E();
+    Float_t pos[3];
+    clus->GetPosition(pos);
+    TVector3 vpos(pos);
+    cl->fR        = vpos.Perp();
+    cl->fEta      = vpos.Eta();
+    cl->fPhi      = vpos.Phi();
+    cl->fN        =  clus->GetNCells();
+    cl->fN1       = GetNCells(clus,0.100);
+    cl->fN3       = GetNCells(clus,0.300);
+    Short_t id    = -1;
+    Double_t emax = GetMaxCellEnergy(clus, id);
+    cl->fIdMax    = id;
+    cl->fEmax     = emax;
+    cl->fDbc      = clus->GetDistanceToBadChannel();;
+    cl->fDisp     = clus->GetDispersion();
+    cl->fM20      = clus->GetM20();
+    cl->fM02      = clus->GetM02();
+    cl->fEcc      = clus->Chi2();   //eccentricity
+    cl->fSig      = clus->GetTOF(); //sigma
+    cl->fTrDz     = fClusProps[i].fTrDz;
+    cl->fTrDr     = fClusProps[i].fTrDr;;
+    cl->fTrEp     = fClusProps[i].fTrEp;;
+    cl->fTrIso    = fClusProps[i].fTrIso;
+    cl->fCeIso    = fClusProps[i].fCellIso;
+  }
+  fNtuple->Fill();
 }
 
 //________________________________________________________________________
@@ -1098,7 +1178,36 @@ void AliAnalysisTaskEMCALPi0PbPb::FillPionHists()
 //________________________________________________________________________
 void AliAnalysisTaskEMCALPi0PbPb::FillOtherHists()
 {
-  // Fill histograms related to cell properties.
+  // Fill other histograms.
+}
+
+//__________________________________________________________________________________________________
+void AliAnalysisTaskEMCALPi0PbPb::FillVertex(AliStaVertex *v, const AliESDVertex *esdv)
+{
+  // Fill vertex from ESD vertex info.
+
+  v->fVx   = esdv->GetX();
+  v->fVy   = esdv->GetY();
+  v->fVz   = esdv->GetZ();
+  v->fVc   = esdv->GetNContributors();
+  v->fDisp = esdv->GetDispersion();
+  v->fZres = esdv->GetZRes();
+  v->fChi2 = esdv->GetChi2();
+  v->fSt   = esdv->GetStatus();
+  v->fIs3D = esdv->IsFromVertexer3D();
+  v->fIsZ  = esdv->IsFromVertexerZ();
+}
+
+//__________________________________________________________________________________________________
+void AliAnalysisTaskEMCALPi0PbPb::FillVertex(AliStaVertex *v, const AliAODVertex *aodv)
+{
+  // Fill vertex from AOD vertex info.
+
+  v->fVx   = aodv->GetX();
+  v->fVy   = aodv->GetY();
+  v->fVz   = aodv->GetZ();
+  v->fVc   = aodv->GetNContributors();
+  v->fChi2 = aodv->GetChi2();
 }
 
 //________________________________________________________________________
@@ -1130,10 +1239,11 @@ Double_t AliAnalysisTaskEMCALPi0PbPb::GetCellIsolation(Double_t cEta, Double_t c
 }
 
 //________________________________________________________________________
-Double_t AliAnalysisTaskEMCALPi0PbPb::GetMaxCellEnergy(AliVCluster *cluster) const
+Double_t AliAnalysisTaskEMCALPi0PbPb::GetMaxCellEnergy(AliVCluster *cluster, Short_t &id) const
 {
   // Get maximum energy of attached cell.
 
+  id = -1;
   Double_t maxe = 0;
   Int_t ncells = cluster->GetNCells();
   if (fEsdCells) {
@@ -1141,6 +1251,7 @@ Double_t AliAnalysisTaskEMCALPi0PbPb::GetMaxCellEnergy(AliVCluster *cluster) con
       Double_t e = fEsdCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
       if (e>maxe) {
         maxe = e;
+        id   = cluster->GetCellAbsId(i);
       }
     }
   } else {
@@ -1148,6 +1259,7 @@ Double_t AliAnalysisTaskEMCALPi0PbPb::GetMaxCellEnergy(AliVCluster *cluster) con
       Double_t e = fAodCells->GetCellAmplitude(TMath::Abs(cluster->GetCellAbsId(i)));
       if (e>maxe)
         maxe = e;
+        id   = cluster->GetCellAbsId(i);
     }
   }
   return maxe;
index 83ceaff..1d68077 100644 (file)
@@ -13,13 +13,17 @@ class AliAODCaloCells;
 class AliAODCaloCluster;
 class AliAODEvent;
 class AliAODTrack;
+class AliAODVertex;
 class AliEMCALGeoUtils;
 class AliEMCALRecoUtils;
 class AliESDCaloCells;
 class AliESDCaloCluster;
 class AliESDEvent;
 class AliESDTrack;
+class AliESDVertex;
 class AliESDtrackCuts;
+class AliStaHeader;
+class AliStaVertex;
 
 #include "AliAnalysisTaskSE.h"
 
@@ -59,10 +63,15 @@ class AliAnalysisTaskEMCALPi0PbPb : public AliAnalysisTaskSE {
   virtual void ClusterAfterburner();
   virtual void FillCellHists();
   virtual void FillClusHists();
-  virtual void FillPionHists();
+  virtual void FillNtuple();
   virtual void FillOtherHists();
+  virtual void FillPionHists();
+  void         FillVertex(AliStaVertex *v, const AliESDVertex *esdv);
+  void         FillVertex(AliStaVertex *v, const AliAODVertex *aodv);
+
   Double_t     GetCellIsolation(Double_t cEta, Double_t cPhi, Double_t radius=0.2)                const;
-  Double_t     GetMaxCellEnergy(AliVCluster *c)                                                   const;
+  Double_t     GetMaxCellEnergy(AliVCluster *c) const { Short_t id=-1; return GetMaxCellEnergy(c,id); }
+  Double_t     GetMaxCellEnergy(AliVCluster *c, Short_t &id)                                      const;
   Int_t        GetNCells(AliVCluster *c, Double_t emin=0.)                                        const;
   void         GetSigma(AliVCluster *c, Double_t &sigmaMax, Double_t &sigmaMin)                   const;
   Double_t     GetTrackIsolation(Double_t cEta, Double_t cPhi, Double_t radius=0.2)               const;
@@ -119,9 +128,15 @@ class AliAnalysisTaskEMCALPi0PbPb : public AliAnalysisTaskSE {
   TObjArray             *fAodClusters;            //!pointer to aod clusters
   AliAODCaloCells       *fAodCells;               //!pointer to aod cells
   TAxis                 *fPtRanges;               //!pointer to pt ranges
-  TNtuple               *fNtuple;                 //!pointer to ntuple
   TObjArray             *fSelTracks;              //!pointer to selected tracks
   ClusProps              fClusProps[1000];        //!array of cluster properties
+    // ntuple
+  TTree                 *fNtuple;                 //!pointer to ntuple
+  AliStaHeader          *fHeader;                 //!pointer to header
+  AliStaVertex          *fPrimVert;               //!pointer to primary vertex
+  AliStaVertex          *fSpdVert;                //!pointer to SPD vertex
+  AliStaVertex          *fTpcVert;                //!pointer to TPC vertex
+  TClonesArray          *fClusters;               //!pointer to clusters
     // histograms
   TH1                   *fHCuts;                  //!histo for cuts
   TH1                   *fHVertexZ;               //!histo for vtxz
@@ -165,3 +180,96 @@ class AliAnalysisTaskEMCALPi0PbPb : public AliAnalysisTaskSE {
   ClassDef(AliAnalysisTaskEMCALPi0PbPb, 5); // Analysis task for neutral pions in Pb+Pb
 };
 #endif
+
+#ifndef AliStaObjs_h
+#define AliStaObjs_h
+class AliStaHeader
+{
+ public:
+  AliStaHeader() : fRun(0), fOrbit(0), fPeriod(0), fBx(0), fL0(0), fL1(0), fL2(0),
+                   fTrClassMask(0), fTrCluster(0), fOffTriggers(0), fFiredTriggers(),
+                   fTcls(0), fV0Cent(0), fCl1Cent(0), fTrCent(0), fCqual(-1) {;}
+  virtual ~AliStaHeader() {;}
+  ULong64_t     GetEventId() const {
+                  return (((ULong64_t)fPeriod << 36) |
+                          ((ULong64_t)fOrbit  << 12) |
+                          (ULong64_t)fBx); 
+                }
+
+ public:
+  Int_t         fRun;            //         run number
+  UInt_t        fOrbit;          //         orbit number
+  UInt_t        fPeriod;         //         period number
+  UShort_t      fBx;             //         bunch crossing id
+  UInt_t        fL0;             //         l0 trigger bits
+  UInt_t        fL1;             //         l1 trigger bits
+  UShort_t      fL2;             //         l2 trigger bits
+  ULong64_t     fTrClassMask;    //         trigger class mask
+  UChar_t       fTrCluster;      //         trigger cluster mask
+  UInt_t        fOffTriggers;    //         fired offline triggers for this event
+  TString       fFiredTriggers;  //         string with fired triggers
+  UInt_t        fTcls;           //         custom trigger definition
+  Double32_t    fV0Cent;         //[0,0,16] v0 cent
+  Double32_t    fCl1Cent;        //[0,0,16] cl1 cent
+  Double32_t    fTrCent;         //[0,0,16] tr cent
+  Int_t         fCqual;          //         centrality quality
+
+  ClassDef(AliStaHeader,1) // Header class
+};
+
+class AliStaVertex
+{
+ public:
+  AliStaVertex(Double_t x=0, Double_t y=0, Double_t z=0) : fVx(x), fVy(y), fVz(z), fVc(-1), fDisp(0), fZres(0),
+                                                           fChi2(0), fSt(0), fIs3D(0), fIsZ(0) {;}
+  virtual ~AliStaVertex() {;}
+
+ public:
+  Double_t      fVx;          //[0,0,16] vertex x
+  Double_t      fVy;          //[0,0,16] vertex y
+  Double_t      fVz;          //[0,0,16] vertex z
+  Double_t      fVc;          //[0,0,16] number of contributors to vertex
+  Double_t      fDisp;        //[0,0,16] dispersion
+  Double_t      fZres;        //[0,0,16] z-resolution
+  Double_t      fChi2;        //[0,0,16] chi2 of fit
+  Bool_t        fSt;          //         status bit
+  Bool_t        fIs3D;        //         is vertex from 3D
+  Bool_t        fIsZ;         //         is vertex from Z only
+
+  ClassDef(AliStaVertex,1) // Vertex class
+};
+
+class AliStaCluster : public TObject
+{
+ public:
+    AliStaCluster() : TObject(), fE(0), fR(0), fEta(0), fPhi(0), fN(0), fN1(0), fN3(0), fIdMax(0), fEmax(0),  
+                      fDbc(0), fDisp(0), fM20(0), fM02(0), fEcc(0), fSig(0), fTrDz(0), fTrDr(-1), fTrEp(0), 
+                      fTrIso(0), fCeIso(0) {;}
+
+//  void          GetMom(TLorentzVector& p, Double_t *vertex=0);
+
+ public:
+  Double32_t    fE;                //[0,0,16] energy
+  Double32_t    fR;                //[0,0,16] radius
+  Double32_t    fEta;              //[0,0,16] eta
+  Double32_t    fPhi;              //[0,0,16] phi
+  UChar_t       fN;                //         number of cells
+  UChar_t       fN1;               //         number of cells > 100 MeV
+  UChar_t       fN3;               //         number of cells > 300 MeV
+  UShort_t      fIdMax;            //         id maximum cell
+  Double32_t    fEmax;             //[0,0,16] energy of maximum cell
+  Double32_t    fDbc;              //[0,0,16] distance to nearest bad channel
+  Double32_t    fDisp;             //[0,0,16] cluster dispersion, for shape analysis
+  Double32_t    fM20;              //[0,0,16] 2-nd moment along the main eigen axis
+  Double32_t    fM02;              //[0,0,16] 2-nd moment along the second eigen axis
+  Double32_t    fEcc;              //[0,0,16] eccentricity
+  Double32_t    fSig;              //[0,0,16] sigma
+  Double32_t    fTrDz;             //[0,0,16] dZ to nearest track
+  Double32_t    fTrDr;             //[0,0,16] dR to nearest track (in x,y; if neg then no match)
+  Double32_t    fTrEp;             //[0,0,16] E/P to nearest track 
+  Double32_t    fTrIso;            //[0,0,16] track isolation
+  Double32_t    fCeIso;            //[0,0,16] cell isolation
+
+  ClassDef(AliStaCluster,1) // Cluster class
+};
+#endif
index d69d5aa..8fec4c5 100644 (file)
@@ -10,5 +10,8 @@
 #pragma link C++ class AliAnalysisTaskEMCALClusterize+;
 #pragma link C++ class AliAnalysisTaskEMCALClusterizeFast+;
 #pragma link C++ class AliAnalysisTaskEMCALPi0PbPb+;
+#pragma link C++ class AliStaHeader+;
+#pragma link C++ class AliStaCluster+;
+#pragma link C++ class AliStaVertex+;
 
 #endif