updates from Salvatore
authorloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 11 Aug 2012 18:26:24 +0000 (18:26 +0000)
committerloizides <loizides@f7af4fe6-9843-0410-8265-dc069ae4e863>
Sat, 11 Aug 2012 18:26:24 +0000 (18:26 +0000)
PWGJE/EMCALJetTasks/AliJetRandomizerTask.h
PWGJE/EMCALJetTasks/AliJetResponseMaker.cxx
PWGJE/EMCALJetTasks/AliJetResponseMaker.h
PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskSAJF.cxx
PWGJE/EMCALJetTasks/macros/AddTaskJetRespPtHard.C [new file with mode: 0644]
PWGJE/EMCALJetTasks/macros/AddTaskJetResponseMaker.C

index efb5e63..cda88ca 100644 (file)
@@ -27,6 +27,6 @@ class AliJetRandomizerTask : public AliJetModelBaseTask {
   AliJetRandomizerTask(const AliJetRandomizerTask&);            // not implemented
   AliJetRandomizerTask &operator=(const AliJetRandomizerTask&); // not implemented
 
-  ClassDef(AliJetRandomizerTask, 1) // Jet randomizer task
+  ClassDef(AliJetRandomizerTask, 2) // Jet randomizer task
 };
 #endif
index 9fb3d82..b7619ac 100644 (file)
@@ -31,20 +31,21 @@ AliJetResponseMaker::AliJetResponseMaker() :
   fMCJetsName("MCJets"),
   fMaxDistance(0.25),
   fDoWeighting(kFALSE),
-  fVertexCut(10),
+  fEventWeightHist(kFALSE),
   fMCFiducial(kFALSE),
   fMCminEta(0),
   fMCmaxEta(0),
   fMCminPhi(0),
   fMCmaxPhi(0),
+  fSelectPtHardBin(-999),
   fPythiaHeader(0),
   fEventWeight(0),
   fPtHardBin(0),
   fNTrials(0),
   fMCTracks(0),
   fMCJets(0),
+  fHistZVertex(0),
   fHistNTrials(0),
-  fHistAcceptedEvents(0),
   fHistEvents(0),
   fHistMCJetPhiEta(0),
   fHistMCJetsPt(0),
@@ -66,6 +67,10 @@ AliJetResponseMaker::AliJetResponseMaker() :
   fHistMissedMCJets(0)
 {
   // Default constructor.
+
+  for (Int_t i = 0; i < 11; i++) {
+    fHistEventWeight[i] = 0;
+  }
 }
 
 //________________________________________________________________________
@@ -75,20 +80,21 @@ AliJetResponseMaker::AliJetResponseMaker(const char *name) :
   fMCJetsName("MCJets"),
   fMaxDistance(0.25),
   fDoWeighting(kFALSE),
-  fVertexCut(10),
+  fEventWeightHist(kFALSE),
   fMCFiducial(kFALSE),
   fMCminEta(0),
   fMCmaxEta(0),
   fMCminPhi(0),
   fMCmaxPhi(0),
+  fSelectPtHardBin(-999),
   fPythiaHeader(0),
   fEventWeight(0),
   fPtHardBin(0),
   fNTrials(0),
   fMCTracks(0),
   fMCJets(0),
+  fHistZVertex(0),
   fHistNTrials(0),
-  fHistAcceptedEvents(0),
   fHistEvents(0),
   fHistMCJetPhiEta(0),
   fHistMCJetsPt(0),
@@ -110,6 +116,10 @@ AliJetResponseMaker::AliJetResponseMaker(const char *name) :
   fHistMissedMCJets(0)
 {
   // Standard constructor.
+
+  for (Int_t i = 0; i < 11; i++) {
+    fHistEventWeight[i] = 0;
+  }
 }
 
 //________________________________________________________________________
@@ -130,24 +140,31 @@ void AliJetResponseMaker::UserCreateOutputObjects()
   const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
   const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
 
+  fHistZVertex = new TH1F("fHistZVertex","Z vertex position", 60, -30, 30);
+  fHistZVertex->GetXaxis()->SetTitle("z");
+  fHistZVertex->GetYaxis()->SetTitle("counts");
+  fOutput->Add(fHistZVertex);
+
   fHistNTrials = new TH1F("fHistNTrials", "fHistNTrials", 11, 0, 11);
   fHistNTrials->GetXaxis()->SetTitle("p_{T} hard bin");
   fHistNTrials->GetYaxis()->SetTitle("trials");
   fOutput->Add(fHistNTrials);
 
-  fHistAcceptedEvents = new TH1F("fHistAcceptedEvents", "fHistAcceptedEvents", 11, 0, 11);
-  fHistAcceptedEvents->GetXaxis()->SetTitle("p_{T} hard bin");
-  fHistAcceptedEvents->GetYaxis()->SetTitle("accepted events");
-  fOutput->Add(fHistAcceptedEvents);
-
   fHistEvents = new TH1F("fHistEvents", "fHistEvents", 11, 0, 11);
   fHistEvents->GetXaxis()->SetTitle("p_{T} hard bin");
   fHistEvents->GetYaxis()->SetTitle("total events");
   fOutput->Add(fHistEvents);
 
+  if (fEventWeightHist) {
+    for (Int_t i = 0; i < 11; i++) {
+      TString name(Form("fHistEventWeight_%d", i+1));
+      fHistEventWeight[i] = new TH1F(name, name, 10, 0, 10);
+      fOutput->Add(fHistEventWeight[i]);
+    }
+  }
+
   for (Int_t i = 1; i < 12; i++) {
     fHistNTrials->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
-    fHistAcceptedEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
     fHistEvents->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
   }
 
@@ -312,6 +329,17 @@ void AliJetResponseMaker::ExecOnce()
 }
 
 //________________________________________________________________________
+Bool_t AliJetResponseMaker::IsEventSelected()
+{
+  // Check if event is selected
+
+  if (fSelectPtHardBin != -999 && fSelectPtHardBin != fPtHardBin) 
+    return kFALSE;
+
+  return AliAnalysisTaskEmcalJet::IsEventSelected();
+}
+
+//________________________________________________________________________
 Bool_t AliJetResponseMaker::RetrieveEventObjects()
 {
   // Retrieve event objects.
@@ -349,14 +377,6 @@ Bool_t AliJetResponseMaker::Run()
 {
   // Find the closest jets
 
-  fHistEvents->SetBinContent(fPtHardBin + 1, fHistEvents->GetBinContent(fPtHardBin + 1) + 1);
-
-  if (TMath::Abs(fVertex[2]) > fVertexCut)
-    return kFALSE;
-  
-  fHistAcceptedEvents->SetBinContent(fPtHardBin + 1, fHistEvents->GetBinContent(fPtHardBin + 1) + 1);
-  fHistNTrials->SetBinContent(fPtHardBin + 1, fHistNTrials->GetBinContent(fPtHardBin + 1) + fNTrials);
-
   DoJetLoop(fJets, fMCJets, kFALSE);
   DoJetLoop(fMCJets, fJets, kTRUE);
 
@@ -442,6 +462,12 @@ Bool_t AliJetResponseMaker::FillHistograms()
 {
   // Fill histograms.
 
+  fHistEvents->SetBinContent(fPtHardBin + 1, fHistEvents->GetBinContent(fPtHardBin + 1) + 1);
+  fHistNTrials->SetBinContent(fPtHardBin + 1, fHistNTrials->GetBinContent(fPtHardBin + 1) + fNTrials);
+  if (fEventWeightHist)
+    fHistEventWeight[fPtHardBin]->Fill(fPythiaHeader->EventWeight());
+  fHistZVertex->Fill(fVertex[2]);
+
   const Int_t nMCJets = fMCJets->GetEntriesFast();
 
   for (Int_t i = 0; i < nMCJets; i++) {
index cd015fa..92f9a16 100644 (file)
@@ -18,14 +18,16 @@ class AliJetResponseMaker : public AliAnalysisTaskEmcalJet {
 
   void                        UserCreateOutputObjects();
 
-  void                        SetMCJetsName(const char *n)       { fMCJetsName    = n; }
-  void                        SetMCTracksName(const char *n)     { fMCTracksName  = n; }
-  void                        SetMaxDistance(Double_t d)         { fMaxDistance   = d; }
-  void                        SetDoWeighting(Bool_t d = kTRUE)   { fDoWeighting   = d; }
-  void                        SetVertexCut(Double_t v)           { fVertexCut     = v; }
-  void                        SetMCFiducial(Bool_t f = kTRUE)    { fMCFiducial    = f; }
+  void                        SetMCJetsName(const char *n)       { fMCJetsName      = n; }
+  void                        SetMCTracksName(const char *n)     { fMCTracksName    = n; }
+  void                        SetMaxDistance(Double_t d)         { fMaxDistance     = d; }
+  void                        SetDoWeighting(Bool_t d = kTRUE)   { fDoWeighting     = d; }
+  void                        SetMCFiducial(Bool_t f = kTRUE)    { fMCFiducial      = f; }
+  void                        SetEventWeightHist(Bool_t h)       { fEventWeightHist = h; }
+  void                        SetPtHardBin(Int_t b)              { fSelectPtHardBin = b; }
 
  protected:
+  Bool_t                      IsEventSelected();
   Bool_t                      AcceptJet(AliEmcalJet* jet, Bool_t /*bias*/ = kFALSE, Bool_t /*upCut*/ = kFALSE)   const;
   void                        ExecOnce();
   void                        DoJetLoop(TClonesArray *jets1, TClonesArray *jets2, Bool_t mc);
@@ -37,12 +39,13 @@ class AliJetResponseMaker : public AliAnalysisTaskEmcalJet {
   TString                     fMCJetsName;                // name of MC jet collection
   Double_t                    fMaxDistance;               // maximum distance between matched particle and detector level jets
   Bool_t                      fDoWeighting;               // = true, weight using trials and given x section
-  Double_t                    fVertexCut;                 // z vertex cut
+  Bool_t                      fEventWeightHist;           // = true create event weight histogram
   Bool_t                      fMCFiducial;                // = true MC jets in fiducial acceptance
   Double_t                    fMCminEta;                  // MC jets minimum eta
   Double_t                    fMCmaxEta;                  // MC jets maximum eta
   Double_t                    fMCminPhi;                  // MC jets minimum phi
   Double_t                    fMCmaxPhi;                  // MC jets maximum phi
+  Int_t                       fSelectPtHardBin;           // Select one pt hard bin for analysis
 
   AliGenPythiaEventHeader    *fPythiaHeader;              //!event Pythia header
   Double_t                    fEventWeight;               //!event weight
@@ -51,9 +54,10 @@ class AliJetResponseMaker : public AliAnalysisTaskEmcalJet {
   TClonesArray               *fMCTracks;                  //!MC particles
   TClonesArray               *fMCJets;                    //!MC jets
   // General histograms
+  TH1F                       *fHistZVertex;               //!Z vertex position
   TH1F                       *fHistNTrials;               //!total number of trials per pt hard bin
-  TH1F                       *fHistAcceptedEvents;        //!number of accepted events per pt hard bin
   TH1F                       *fHistEvents;                //!total number of events per pt hard bin
+  TH1F                       *fHistEventWeight[11];       //!event weight
   // Particle level jets
   TH2F                       *fHistMCJetPhiEta;           //!phi-eta distribution of jets
   TH1F                       *fHistMCJetsPt;              //!inclusive jet pt spectrum
@@ -80,6 +84,6 @@ class AliJetResponseMaker : public AliAnalysisTaskEmcalJet {
   AliJetResponseMaker(const AliJetResponseMaker&);            // not implemented
   AliJetResponseMaker &operator=(const AliJetResponseMaker&); // not implemented
 
-  ClassDef(AliJetResponseMaker, 5) // Jet response matrix producing task
+  ClassDef(AliJetResponseMaker, 6) // Jet response matrix producing task
 };
 #endif
index 725ad5d..9504be3 100644 (file)
@@ -812,8 +812,6 @@ void AliAnalysisTaskSAJF::DoJetLoop()
   if (!fJets)
     return;
 
-  const Double_t rho = fRho->GetVal();
-
   const Int_t njets = fJets->GetEntriesFast();
 
   TH1F constituents("constituents", "constituents", 
@@ -831,7 +829,7 @@ void AliAnalysisTaskSAJF::DoJetLoop()
     if (!AcceptJet(jet))
       continue;
 
-    Float_t corrPt = jet->Pt() - rho * jet->Area();
+    Float_t corrPt = jet->Pt() - fRhoVal * jet->Area();
 
     fHistJetsPt[fCentBin]->Fill(jet->Pt());
     fHistJetsPtArea[fCentBin]->Fill(jet->Pt(), jet->Area());
diff --git a/PWGJE/EMCALJetTasks/macros/AddTaskJetRespPtHard.C b/PWGJE/EMCALJetTasks/macros/AddTaskJetRespPtHard.C
new file mode 100644 (file)
index 0000000..657cc6b
--- /dev/null
@@ -0,0 +1,32 @@
+// $Id$
+
+AliJetResponseMaker* AddTaskJetRespPtHard(const char *ntracks            = "Tracks",
+                                         const char *nclusters          = "CaloClusters",
+                                         const char *njets              = "Jets",
+                                         const char *nmcjets            = "MCJets",
+                                         const char *nmctracks          = "MCParticles",
+                                         Double_t    jetradius          = 0.4,
+                                         Double_t    jetptcut           = 1,
+                                         Double_t    jetareacut         = 0.8,
+                                         Double_t    ptcut              = 0.15,
+                                         Double_t    jetBiasTrack       = 10,
+                                         Double_t    jetBiasClus        = 10,
+                                         Double_t    maxDistance        = 0.25,
+                                         UInt_t      type               = AliAnalysisTaskEmcal::kTPC,
+                                         Int_t minPtBin                 = 1, 
+                                         Int_t maxPtBin                 = 11,
+                                         const char *taskname           = "AliJetResponseMaker"
+)
+{  
+  gROOT->LoadMacro("$ALICE_ROOT/PWGJE/macros/AddTaskJetResponseMaker.C");
+  
+  AliJetResponseMaker *jetTask = 0;
+
+  for (Int_t i = minPtBin; i <= maxPtBin; i++) {
+    jetTask = AddTaskJetResponseMaker(ntracks, nclusters, njets, nmcjets, nmctracks, 
+                                     jetradius, jetptcut, jetareacut, ptcut, jetBiasTrack, 
+                                     jetBiasClus, maxDistance, type, i, taskname);
+  }
+  
+  return jetTask;
+}
index d0e7338..1421270 100644 (file)
@@ -1,7 +1,6 @@
 // $Id$
 
 AliJetResponseMaker* AddTaskJetResponseMaker(
-  const char *taskname           = "AliJetResponseMaker",
   const char *ntracks            = "Tracks",
   const char *nclusters          = "CaloClusters",
   const char *njets              = "Jets",
@@ -14,7 +13,9 @@ AliJetResponseMaker* AddTaskJetResponseMaker(
   Double_t    jetBiasTrack       = 10,
   Double_t    jetBiasClus        = 10,
   Double_t    maxDistance        = 0.25,
-  UInt_t      type               = AliAnalysisTaskEmcal::kTPC
+  UInt_t      type               = AliAnalysisTaskEmcal::kTPC,
+  Int_t       ptHardBin          = -999,
+  const char *taskname           = "AliJetResponseMaker"
 )
 {  
   // Get the pointer to the existing analysis manager via the static access method.
@@ -38,7 +39,32 @@ AliJetResponseMaker* AddTaskJetResponseMaker(
   // Init the task and do settings
   //-------------------------------------------------------
 
-  AliJetResponseMaker* jetTask = new AliJetResponseMaker(taskname);
+  TString name(taskname);
+  name += "_";
+  name += njets;
+  name += "_Track";
+  name += jetBiasTrack;
+  name += "_Clus";
+  name += jetBiasClus;
+  name += "_R0";
+  name += floor(jetradius*10+0.5);
+  name += "_PtCut";
+  name += floor(ptcut*1000+0.5);
+  name += "_";
+  if (type == AliAnalysisTaskEmcal::kTPC) 
+    name += "TPC";
+  else if (type == AliAnalysisTaskEmcal::kEMCAL) 
+    name += "EMCAL";
+  else if (type == AliAnalysisTaskEmcal::kTPCSmall) 
+    name += "TPCSmall";
+  else if (type == AliAnalysisTaskEmcal::kEMCALOnly) 
+    name += "EMCALOnly";
+  if (ptHardBin != -999) {
+    name += "_PtHard";
+    name += ptHardBin;
+  }
+
+  AliJetResponseMaker* jetTask = new AliJetResponseMaker(name);
   jetTask->SetAnaType(type);
   jetTask->SetTracksName(ntracks);
   jetTask->SetClusName(nclusters);
@@ -52,6 +78,8 @@ AliJetResponseMaker* AddTaskJetResponseMaker(
   jetTask->SetPtBiasJetTrack(jetBiasTrack);
   jetTask->SetPtBiasJetClus(jetBiasClus);
   jetTask->SetMaxDistance(maxDistance);
+  jetTask->SetVzRange(-10,10);
+  jetTask->SetPtHardBin(ptHardBin);
   
   //-------------------------------------------------------
   // Final settings, pass to manager and set the containers
@@ -61,7 +89,7 @@ AliJetResponseMaker* AddTaskJetResponseMaker(
   
   // Create containers for input/output
   AliAnalysisDataContainer *cinput1  = mgr->GetCommonInputContainer()  ;
-  TString contname(taskname);
+  TString contname(name);
   contname += "_histos";
   AliAnalysisDataContainer *coutput1 = mgr->CreateContainer(contname.Data(), 
                                                            TList::Class(),AliAnalysisManager::kOutputContainer,