]> git.uio.no Git - u/mrichter/AliRoot.git/commitdiff
updates jet mass detector response
authormverweij <marta.verweij@cern.ch>
Mon, 1 Sep 2014 09:05:02 +0000 (11:05 +0200)
committermverweij <marta.verweij@cern.ch>
Tue, 2 Sep 2014 13:25:26 +0000 (15:25 +0200)
add track fast simulation class

PWGJE/CMakelibPWGJEEMCALJetTasks.pkg
PWGJE/EMCALJetTasks/AliJetEmbeddingFromGenTask.cxx
PWGJE/EMCALJetTasks/AliJetFastSimulation.cxx [new file with mode: 0644]
PWGJE/EMCALJetTasks/AliJetFastSimulation.h [new file with mode: 0644]
PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskJetMassResponseDet.cxx
PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskJetMassResponseDet.h
PWGJE/EMCALJetTasks/macros/AddTaskFastSimulation.C [new file with mode: 0644]
PWGJE/EMCALJetTasks/macros/AddTaskJetMassResponseDet.C
PWGJE/PWGJEEMCALJetTasksLinkDef.h

index 001d5fef758063b8ed079eef2e2da52fe6359b25..e58fe3111bdab39519bc9d106d1309cea877fda4 100644 (file)
@@ -42,6 +42,7 @@ set ( SRCS
  EMCALJetTasks/AliJetEmbeddingFromGenTask.cxx
  EMCALJetTasks/AliJetEmbeddingFromPYTHIATask.cxx
  EMCALJetTasks/AliJetEmbeddingTask.cxx
+ EMCALJetTasks/AliJetFastSimulation.cxx
  EMCALJetTasks/AliJetModelBaseTask.cxx
  EMCALJetTasks/AliJetModelMergeBranches.cxx
  EMCALJetTasks/AliJetRandomizerTask.cxx
index 6fe2b3ddf77825818ab158ca542ca2311d5eea1e..beb9614fdf7470cd6f8d61d0f4af8be8bb5349d2 100644 (file)
@@ -82,7 +82,7 @@ void AliJetEmbeddingFromGenTask::UserCreateOutputObjects()
   fHistPt = new TH1F("fHistpt","fHistPt;#it{p}_{T};N",100,0.,100.);
   fOutput->Add(fHistPt);
 
-  fHistEtaPhi = new TH2F("fHistEtapHI","fHistEtaPhi;#eta;#varphi",100,-3.,3.,100.,0.,TMath::TwoPi());
+  fHistEtaPhi = new TH2F("fHistEtaPhi","fHistEtaPhi;#eta;#varphi",100,-3.,3.,100.,0.,TMath::TwoPi());
   fOutput->Add(fHistEtaPhi);
 
   fHistTrials = new TH1F("fHistTrials", "fHistTrials", 1, 0, 1);
diff --git a/PWGJE/EMCALJetTasks/AliJetFastSimulation.cxx b/PWGJE/EMCALJetTasks/AliJetFastSimulation.cxx
new file mode 100644 (file)
index 0000000..0ac3f1d
--- /dev/null
@@ -0,0 +1,411 @@
+// $Id$
+//
+// Jet model task to merge to existing branches
+// only implemented for track branches
+//
+// Author: M. Verweij
+
+#include "AliJetFastSimulation.h"
+
+#include <TClonesArray.h>
+#include <TFolder.h>
+#include <TLorentzVector.h>
+#include <TParticle.h>
+#include <TParticlePDG.h>
+#include <TRandom3.h>
+#include <TProfile.h>
+#include <TGrid.h>
+#include <TFile.h>
+#include <TF1.h>
+#include "AliAnalysisManager.h"
+#include "AliEMCALDigit.h"
+#include "AliEMCALGeometry.h"
+#include "AliEMCALRecPoint.h"
+#include "AliGenerator.h"
+#include "AliHeader.h"
+#include "AliLog.h"
+#include "AliPicoTrack.h"
+#include "AliRun.h"
+#include "AliRunLoader.h"
+#include "AliStack.h"
+#include "AliStack.h"
+#include "AliVCluster.h"
+#include "AliVEvent.h"
+
+ClassImp(AliJetFastSimulation)
+
+//________________________________________________________________________
+AliJetFastSimulation::AliJetFastSimulation() : 
+AliAnalysisTaskEmcal("AliJetFastSimulation",kTRUE),
+  fTracksOutName(""),
+  fTracksOut(0x0),
+  fNTrackClasses(2),
+  fRandom(0),
+  fEfficiencyFixed(1.),
+  fMomResH1(0x0),
+  fMomResH2(0x0),
+  fMomResH3(0x0),
+  fMomResH1Fit(0x0),
+  fMomResH2Fit(0x0),
+  fMomResH3Fit(0x0),
+  fhEffH1(0x0),
+  fhEffH2(0x0),
+  fhEffH3(0x0),
+  fUseTrPtResolutionSmearing(kFALSE),
+  fUseDiceEfficiency(kFALSE),
+  fDiceEfficiencyMinPt(-1.),
+  fUseTrPtResolutionFromOADB(kFALSE),
+  fUseTrEfficiencyFromOADB(kFALSE),
+  fPathTrPtResolution(""),
+  fPathTrEfficiency(""),
+  fHistPtDet(0),
+  fh2PtGenPtSmeared(0),
+  fp1Efficiency(0),
+  fp1PtResolution(0)
+{
+  // Default constructor.
+  SetMakeGeneralHistograms(kTRUE);
+}
+
+//________________________________________________________________________
+AliJetFastSimulation::AliJetFastSimulation(const char *name) : 
+  AliAnalysisTaskEmcal(name,kTRUE),
+  fTracksOutName(""),
+  fTracksOut(0x0),
+  fNTrackClasses(2),
+  fRandom(0),
+  fEfficiencyFixed(1.),
+  fMomResH1(0x0),
+  fMomResH2(0x0),
+  fMomResH3(0x0),
+  fMomResH1Fit(0x0),
+  fMomResH2Fit(0x0),
+  fMomResH3Fit(0x0),
+  fhEffH1(0x0),
+  fhEffH2(0x0),
+  fhEffH3(0x0),
+  fUseTrPtResolutionSmearing(kFALSE),
+  fUseDiceEfficiency(kFALSE),
+  fDiceEfficiencyMinPt(-1.),
+  fUseTrPtResolutionFromOADB(kFALSE),
+  fUseTrEfficiencyFromOADB(kFALSE),
+  fPathTrPtResolution(""),
+  fPathTrEfficiency(""),
+  fHistPtDet(0),
+  fh2PtGenPtSmeared(0),
+  fp1Efficiency(0),
+  fp1PtResolution(0)
+{
+  // Standard constructor.
+  SetMakeGeneralHistograms(kTRUE);
+}
+
+//________________________________________________________________________
+AliJetFastSimulation::~AliJetFastSimulation()
+{
+  // Destructor
+
+  delete fRandom;
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::ExecOnce() 
+{
+  // Exec only once.
+
+  AliAnalysisTaskEmcal::ExecOnce();
+
+  if(!fRandom) fRandom = new TRandom3(0);
+
+  if (!fTracksOutName.IsNull()) {
+    fTracksOut = new TClonesArray("AliPicoTrack");
+    fTracksOut->SetName(fTracksOutName);
+    if (InputEvent()->FindListObject(fTracksOutName)) {
+      AliFatal(Form("%s: Collection %s is already present in the event!", GetName(), fTracksOutName.Data()));
+      return;
+    }
+    else {
+      InputEvent()->AddObject(fTracksOut);
+    }
+  }
+}
+//________________________________________________________________________
+void AliJetFastSimulation::LocalInit() {
+  //initialize track response
+  if(fUseTrPtResolutionFromOADB) LoadTrPtResolutionRootFileFromOADB();
+  if(fUseTrEfficiencyFromOADB)   LoadTrEfficiencyRootFileFromOADB();
+  FitMomentumResolution();
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::UserCreateOutputObjects() 
+{
+  AliAnalysisTaskEmcal::UserCreateOutputObjects();
+
+  const Int_t nBinPt = 100;
+  Double_t binLimitsPt[nBinPt+1];
+  for(Int_t iPt = 0;iPt <= nBinPt;iPt++){
+    if(iPt == 0){
+      binLimitsPt[iPt] = 0.0;
+    } else {// 1.0
+      binLimitsPt[iPt] =  binLimitsPt[iPt-1] + 1.0;
+    }
+  }
+
+  fHistPtDet = new TH1F("fHistpt","fHistPtDet;#it{p}_{T};N",nBinPt,binLimitsPt);
+  fOutput->Add(fHistPtDet);
+
+  fh2PtGenPtSmeared = new TH2F("fh2PtGenPtSmeared","fh2PtGenPtSmeared",nBinPt,binLimitsPt,nBinPt,binLimitsPt);
+  fOutput->Add(fh2PtGenPtSmeared);
+
+  fp1Efficiency = new TProfile("fp1Efficiency","fp1Efficiency",nBinPt,binLimitsPt);
+  fOutput->Add(fp1Efficiency);
+
+  fp1PtResolution = new TProfile("fp1PtResolution","fp1PtResolution",nBinPt,binLimitsPt);
+  fOutput->Add(fp1PtResolution);
+
+  PostData(1, fOutput); // Post data for ALL output slots > 0 here.
+}
+
+
+//________________________________________________________________________
+Bool_t AliJetFastSimulation::Run() 
+{
+  //Check if information is provided detector level effects
+  if(!fMomResH1 && fNTrackClasses>0 ) 
+    fUseTrPtResolutionSmearing = 0;
+  if(!fMomResH2 && fNTrackClasses>1 ) 
+    fUseTrPtResolutionSmearing = 0;
+  if(!fMomResH3 && fNTrackClasses>2 ) 
+    fUseTrPtResolutionSmearing = 0;
+
+  if(fEfficiencyFixed < 1. && !fUseDiceEfficiency)
+    fUseDiceEfficiency = 1; // 1 is the default; 2 can be set by user but not implemented
+  else {
+    if(!fhEffH1 && fNTrackClasses>0 ) 
+      fUseDiceEfficiency = 0;
+    if(!fhEffH2 && fNTrackClasses>1 ) 
+      fUseDiceEfficiency = 0;
+    if(!fhEffH3 && fNTrackClasses>2 ) 
+      fUseDiceEfficiency = 0;
+  }
+
+  SimulateTracks();
+  return kTRUE;
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::SimulateTracks()
+{
+  //Apply toy detector simulation to tracks
+  const Int_t nTracks = fTracks->GetEntriesFast();
+   for (Int_t i = 0; i < nTracks; ++i) {
+    AliPicoTrack *picotrack = static_cast<AliPicoTrack*>(fTracks->At(i));
+    if (!picotrack)
+      continue;
+
+    Bool_t accept = kTRUE;
+    Double_t eff[3] = {0};
+    Double_t rnd = fRandom->Uniform(1.);  
+    if(fUseDiceEfficiency) accept = DiceEfficiency(picotrack,eff,rnd);
+    if(!accept) continue;
+
+    AliPicoTrack *track = NULL;
+    if(fUseTrPtResolutionSmearing) {
+      track = SmearPt(picotrack,eff,rnd);
+      (*fTracksOut)[i] = track;
+    } else
+      track = new ((*fTracksOut)[i]) AliPicoTrack(*picotrack);
+
+    track->SetBit(TObject::kBitMask,1);
+    fHistPtDet->Fill(track->Pt());
+   }
+}
+
+//________________________________________________________________________
+Bool_t AliJetFastSimulation::DiceEfficiency(AliPicoTrack *vp, Double_t eff[3], Double_t rnd)
+{
+  // Dice to decide if particle is kept or not - toy  model for efficiency
+  //
+  Double_t sumEff = 0.;
+  Double_t pT = 0.;
+  if(fEfficiencyFixed<1.)
+    sumEff = fEfficiencyFixed;
+  else {
+    pT = vp->Pt();
+    Double_t pTtmp = pT;
+    if(pT>10.) pTtmp = 10.;
+    if(fhEffH1) eff[0] = fhEffH1->GetBinContent(fhEffH1->FindBin(pTtmp));
+    if(fhEffH2) eff[1] = fhEffH2->GetBinContent(fhEffH2->FindBin(pTtmp));
+    if(fhEffH3) eff[2] = fhEffH3->GetBinContent(fhEffH3->FindBin(pTtmp));
+
+    sumEff = eff[0]+eff[1]+eff[2];
+  }
+  fp1Efficiency->Fill(vp->Pt(),sumEff);
+  if(rnd>sumEff && pT > fDiceEfficiencyMinPt) return kFALSE;
+  return kTRUE;
+}
+
+//________________________________________________________________________
+AliPicoTrack* AliJetFastSimulation::SmearPt(AliPicoTrack *vp, Double_t eff[3], Double_t rnd)
+{
+  //Smear momentum of generated particle
+  Double_t smear = 1.;
+  Double_t  pT = vp->Pt();
+  //Select hybrid track category
+
+  //Sort efficiencies from large to small
+  Int_t cat[3] = {0};
+  TMath::Sort(3,eff,cat);
+  if(rnd<=eff[cat[2]])
+    smear = GetMomentumSmearing(cat[2],pT);
+  else if(rnd<=(eff[cat[2]]+eff[cat[1]]))
+    smear = GetMomentumSmearing(cat[1],pT);
+  else
+    smear = GetMomentumSmearing(cat[0],pT);
+  
+  fp1PtResolution->Fill(vp->Pt(),smear);
+
+  Double_t sigma = vp->Pt()*smear;
+  Double_t pTrec = fRandom->Gaus(vp->Pt(),sigma);
+  fh2PtGenPtSmeared->Fill(vp->Pt(),pTrec);
+
+  AliPicoTrack *picotrack = new AliPicoTrack(pTrec,
+                                            vp->Eta(),
+                                            vp->Phi(),
+                                            vp->Charge(),
+                                            vp->GetLabel(),
+                                            AliPicoTrack::GetTrackType(vp),
+                                            vp->GetTrackEtaOnEMCal(),
+                                            vp->GetTrackPhiOnEMCal(),
+                                            vp->GetTrackPtOnEMCal(),
+                                            vp->IsEMCAL(),
+                                            0.13957); //assume pion mass
+    return picotrack;
+}
+
+//________________________________________________________________________
+Double_t AliJetFastSimulation::GetMomentumSmearing(Int_t cat, Double_t pt) {
+
+  //
+  // Get smearing on generated momentum
+  //
+
+  TProfile *fMomRes = 0x0;
+  if(cat==1 && fMomResH1) fMomRes = (TProfile*)fMomResH1->Clone("fMomRes");
+  if(cat==2 && fMomResH2) fMomRes = (TProfile*)fMomResH2->Clone("fMomRes");
+  if(cat==3 && fMomResH3) fMomRes = (TProfile*)fMomResH3->Clone("fMomRes");
+
+  if(!fMomRes)
+    return 0.;
+
+  Double_t smear = 0.;
+  if(pt>20.) {
+    if(cat==1 && fMomResH1Fit) smear = fMomResH1Fit->Eval(pt);
+    if(cat==2 && fMomResH2Fit) smear = fMomResH2Fit->Eval(pt);
+    if(cat==3 && fMomResH3Fit) smear = fMomResH3Fit->Eval(pt);
+  }
+  else {
+    Int_t bin = fMomRes->FindBin(pt);
+    smear = fRandom->Gaus(fMomRes->GetBinContent(bin),fMomRes->GetBinError(bin));
+  }
+  
+  if(fMomRes) delete fMomRes;
+
+  return smear;
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::LoadTrPtResolutionRootFileFromOADB() {
+
+  if (!gGrid && fPathTrPtResolution.Contains("alien://")) {
+     AliInfo("Trying to connect to AliEn ...");
+     TGrid::Connect("alien://");
+   }
+
+  TFile *f = TFile::Open(fPathTrPtResolution.Data());
+  if(!f)return;
+  TProfile *fProfPtPtSigma1PtGlobSt     = NULL;
+  TProfile *fProfPtPtSigma1PtGlobCnoSPD = NULL;
+  TProfile *fProfPtPtSigma1PtGlobCnoITS = NULL;
+  if(fNTrackClasses>0) fProfPtPtSigma1PtGlobSt     = (TProfile*)f->Get("fProfPtPtSigma1PtGlobSt");
+  if(fNTrackClasses>1) fProfPtPtSigma1PtGlobCnoSPD = (TProfile*)f->Get("fProfPtPtSigma1PtGlobCnoSPD");
+  if(fNTrackClasses>2) fProfPtPtSigma1PtGlobCnoITS = (TProfile*)f->Get("fProfPtPtSigma1PtGlobCnoITS");
+
+  SetSmearResolution(kTRUE);
+  SetMomentumResolutionHybrid(fProfPtPtSigma1PtGlobSt,fProfPtPtSigma1PtGlobCnoSPD,fProfPtPtSigma1PtGlobCnoITS);
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::LoadTrEfficiencyRootFileFromOADB() {
+
+   if (!gGrid && fPathTrPtResolution.Contains("alien://")) {
+     AliInfo("Trying to connect to AliEn ...");
+     TGrid::Connect("alien://");
+   }
+
+  TFile *f = TFile::Open(fPathTrEfficiency.Data());
+  if(!f)return;
+  TH1D *hEffPosGlobSt     = NULL;
+  TH1D *hEffPosGlobCnoSPD = NULL;
+  TH1D *hEffPosGlobCnoITS = NULL;
+  if(fNTrackClasses>0) hEffPosGlobSt = (TH1D*)f->Get("hEffPosGlobSt");
+  if(fNTrackClasses>1) hEffPosGlobCnoSPD = (TH1D*)f->Get("hEffPosGlobCnoSPD");
+  if(fNTrackClasses>2) hEffPosGlobCnoITS = (TH1D*)f->Get("hEffPosGlobCnoITS");
+
+  SetDiceEfficiency(kTRUE);
+  SetEfficiencyHybrid(hEffPosGlobSt,hEffPosGlobCnoSPD,hEffPosGlobCnoITS);
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::SetMomentumResolutionHybrid(TProfile *p1, TProfile *p2, TProfile *p3) {
+  //
+  // set mom res profiles
+  //
+  if(fMomResH1) delete fMomResH1;
+  if(fMomResH2) delete fMomResH2;
+  if(fMomResH3) delete fMomResH3;
+  
+  if(p1) fMomResH1 = new TProfile(*p1);//(TProfile*)p1->Clone("fMomResH1");
+  if(p2) fMomResH2 = new TProfile(*p2);//(TProfile*)p2->Clone("fMomResH2");
+  if(p3) fMomResH3 = new TProfile(*p3);//(TProfile*)p3->Clone("fMomResH3");
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation:: SetEfficiencyHybrid(TH1 *h1, TH1 *h2, TH1 *h3) {
+  //
+  // set tracking efficiency histos
+  //
+  if(h1) fhEffH1 = (TH1*)h1->Clone("fhEffH1");
+  if(h2) fhEffH2 = (TH1*)h2->Clone("fhEffH2");
+  if(h3) fhEffH3 = (TH1*)h3->Clone("fhEffH3");
+}
+
+//________________________________________________________________________
+void AliJetFastSimulation::FitMomentumResolution() {
+  //
+  // Fit linear function on momentum resolution at high pT
+  //
+
+  if(!fMomResH1Fit && fMomResH1) {
+    fMomResH1Fit = new TF1("fMomResH1Fit","[0]+[1]*x",0.,200.);
+    fMomResH1->Fit(fMomResH1Fit,"LL V0","",5.,30.);
+    fMomResH1Fit ->SetRange(5.,100.);
+  }
+
+  if(!fMomResH2Fit && fMomResH2) {
+    fMomResH2Fit = new TF1("fMomResH2Fit","[0]+[1]*x",0.,200.);
+    fMomResH2->Fit(fMomResH2Fit,"LL V0","",5.,30.);
+    fMomResH2Fit ->SetRange(5.,100.);
+  }
+
+  if(!fMomResH3Fit && fMomResH3) {
+    fMomResH3Fit = new TF1("fMomResH3Fit","[0]+[1]*x",0.,200.);
+    fMomResH3->Fit(fMomResH3Fit,"LL V0","",5.,30.);
+    fMomResH3Fit ->SetRange(5.,100.);
+  }
+
+}
+
+
diff --git a/PWGJE/EMCALJetTasks/AliJetFastSimulation.h b/PWGJE/EMCALJetTasks/AliJetFastSimulation.h
new file mode 100644 (file)
index 0000000..be2b8df
--- /dev/null
@@ -0,0 +1,81 @@
+#ifndef ALIJETFASTSIMULATION_H
+#define ALIJETFASTSIMULATION_H
+
+// $Id$
+
+class TClonesArray;
+class TRandom3;
+class AliVParticle;
+class AliPicoTrack;
+
+#include "AliAnalysisTaskEmcal.h"
+
+class AliJetFastSimulation : public AliAnalysisTaskEmcal {
+ public:
+  AliJetFastSimulation();
+  AliJetFastSimulation(const char *name); 
+  virtual ~AliJetFastSimulation();
+
+  virtual void           LocalInit();
+  virtual void           UserCreateOutputObjects();
+
+  void                   SetTracksOutName(const char *n)          { fTracksOutName   = n;    }
+  void                   SetNTrackClasses(Int_t i)                { fNTrackClasses   = i;    }
+  void                   SetFixedTrackEfficiency(Double_t eff)    { fEfficiencyFixed = eff ; }
+
+  void                   SetUseTrResolutionFromOADB(Bool_t b=kTRUE, TString path="$ALICE_ROOT/OADB/PWGJE/Resolution/PtResol_LHCh_Cent0-10_v1.root") {fUseTrPtResolutionFromOADB = b; fPathTrPtResolution=path;}
+  void                   SetUseTrEfficiencyFromOADB(Bool_t b=kTRUE, TString path="$ALICE_ROOT/OADB/PWGJE/Efficiency/Efficiency_LHC11a2aj_Cent0_v1.root") {fUseTrEfficiencyFromOADB = b; fPathTrEfficiency=path;}
+  void                   SetSmearResolution(Bool_t b)                               { fUseTrPtResolutionSmearing = b ;}
+  void                   SetDiceEfficiency(Int_t b)                                 { fUseDiceEfficiency         = b ;}
+  void                   SetDiceEfficiencyMinPt(Double_t pt)                        { fDiceEfficiencyMinPt       = pt;}
+
+ protected:
+  void                   ExecOnce();
+  Bool_t                 Run();
+
+  void                   SimulateTracks();
+  Bool_t                 DiceEfficiency(AliPicoTrack *vp, Double_t eff[3], Double_t rnd);
+  AliPicoTrack          *SmearPt(AliPicoTrack *vp, Double_t eff[3], Double_t rnd);
+  Double_t               GetMomentumSmearing(Int_t cat, Double_t pt);
+  void                   FitMomentumResolution();
+  void                   LoadTrEfficiencyRootFileFromOADB();
+  void                   LoadTrPtResolutionRootFileFromOADB();
+  void                   SetMomentumResolutionHybrid(TProfile *p1, TProfile *p2, TProfile *p3);
+  void                   SetEfficiencyHybrid(TH1 *h1, TH1 *h2, TH1 *h3);
+
+  TString                fTracksOutName;       // name of output track collection
+  TClonesArray          *fTracksOut;           //!output track collection
+  Int_t                  fNTrackClasses;       // number of track classes
+  TRandom3 *fRandom;                           //! random number generator
+  Double_t  fEfficiencyFixed;                  // fixed efficiency for all pT and all types of tracks
+  TProfile *fMomResH1;                         // Momentum resolution from TrackQA Hybrid Category 1
+  TProfile *fMomResH2;                         // Momentum resolution from TrackQA Hybrid Category 2
+  TProfile *fMomResH3;                         // Momentum resolution from TrackQA Hybrid Category 3
+  TF1      *fMomResH1Fit;                      // fit to momentum resolution
+  TF1      *fMomResH2Fit;                      // fit to momentum resolution
+  TF1      *fMomResH3Fit;                      // fit to momentum resolution
+  TH1      *fhEffH1;                           // Efficiency for Spectra Hybrid Category 1
+  TH1      *fhEffH2;                           // Efficiency for Spectra Hybrid Category 2
+  TH1      *fhEffH3;                           // Efficiency for Spectra Hybrid Category 3
+  Bool_t    fUseTrPtResolutionSmearing;        // Apply momentum smearing on track level
+  Int_t     fUseDiceEfficiency;                // Flag to apply efficiency on track level by dicing 0: no dicing; 1: dicing wrt to input branch;
+  Double_t  fDiceEfficiencyMinPt;              // Only do efficiency dicing for tracks above this pt
+  Bool_t    fUseTrPtResolutionFromOADB;        // Load track pt resolution root file from OADB path
+  Bool_t    fUseTrEfficiencyFromOADB;          // Load tracking efficiency root file from OADB path
+  TString   fPathTrPtResolution;               // OADB path to root file
+  TString   fPathTrEfficiency;                 // OADB path to root file
+
+  //Output objects
+  TH1F     *fHistPtDet;                        //!pT spectrum of detector level particles
+  TH2F     *fh2PtGenPtSmeared;                 //! Control histo smeared momentum
+  TProfile *fp1Efficiency;                     //! Control profile efficiency
+  TProfile *fp1PtResolution;                   //! Control profile for pT resolution
+
+  
+ private:
+  AliJetFastSimulation(const AliJetFastSimulation&);            // not implemented
+  AliJetFastSimulation &operator=(const AliJetFastSimulation&); // not implemented
+
+  ClassDef(AliJetFastSimulation, 1) // Jet fast simulation task
+};
+#endif
index b0697dc2c2403497e4fc674f2e9b2d918fd5e2d1..6785cedd48aab7efcd952468cea658fe2b0bb956 100644 (file)
@@ -47,6 +47,8 @@ AliAnalysisTaskJetMassResponseDet::AliAnalysisTaskJetMassResponseDet() :
   fh2PtVsMassJetPartTaggedMatch(0),
   fh2PtVsMassJetDetAll(0),
   fh2PtVsMassJetDetTagged(0),
+  fh2EtaPhiMatchedDet(0),
+  fh2EtaPhiMatchedPart(0),
   fhnMassResponse(0)
 {
   // Default constructor.
@@ -66,6 +68,8 @@ AliAnalysisTaskJetMassResponseDet::AliAnalysisTaskJetMassResponseDet(const char
   fh2PtVsMassJetPartTaggedMatch(0),
   fh2PtVsMassJetDetAll(0),
   fh2PtVsMassJetDetTagged(0),
+  fh2EtaPhiMatchedDet(0),
+  fh2EtaPhiMatchedPart(0),
   fhnMassResponse(0)
 {
   // Standard constructor.
@@ -97,19 +101,19 @@ void AliAnalysisTaskJetMassResponseDet::UserCreateOutputObjects()
   const Double_t minM = 0.;
   const Double_t maxM = 50.;
 
-  const Int_t nBinsMT  = 50;
-  const Double_t minMT = 0.;
-  const Double_t maxMT = 50.;
+  const Int_t nBinsConstEff  = 24;
+  const Double_t minConstEff = 0.;
+  const Double_t maxConstEff = 1.2;
 
-  const Int_t nBinsConst = 26;
-  const Double_t minConst = -5.5;
-  const Double_t maxConst = 20.5;
+  // const Int_t nBinsConst = 26;
+  // const Double_t minConst = -5.5;
+  // const Double_t maxConst = 20.5;
 
   //Binning for THnSparse
   const Int_t nBinsSparse0 = 5;
-  const Int_t nBins0[nBinsSparse0] = {nBinsM,nBinsM,nBinsPt,nBinsPt,nBinsMT};
-  const Double_t xmin0[nBinsSparse0]  = { minM, minM, minPt, minPt, minMT};
-  const Double_t xmax0[nBinsSparse0]  = { maxM, maxM, maxPt, maxPt, maxMT};
+  const Int_t nBins0[nBinsSparse0] = {nBinsM,nBinsM,nBinsPt,nBinsPt,nBinsConstEff};
+  const Double_t xmin0[nBinsSparse0]  = { minM, minM, minPt, minPt, minConstEff};
+  const Double_t xmax0[nBinsSparse0]  = { maxM, maxM, maxPt, maxPt, maxConstEff};
 
   //Create histograms
   TString histName = "";
@@ -145,6 +149,16 @@ void AliAnalysisTaskJetMassResponseDet::UserCreateOutputObjects()
   fh2PtVsMassJetDetTagged = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
   fOutput->Add(fh2PtVsMassJetDetTagged);
 
+  histName = "fh2EtaPhiMatchedDet";
+  histTitle = TString::Format("%s;#eta;#varphi",histName.Data());
+  fh2EtaPhiMatchedDet = new TH2F(histName.Data(),histTitle.Data(),100,-1.,1.,72,0.,TMath::TwoPi());
+  fOutput->Add(fh2EtaPhiMatchedDet);
+
+  histName = "fh2EtaPhiMatchedPart";
+  histTitle = TString::Format("%s;#eta;#varphi",histName.Data());
+  fh2EtaPhiMatchedPart = new TH2F(histName.Data(),histTitle.Data(),100,-1.,1.,72,0.,TMath::TwoPi());
+  fOutput->Add(fh2EtaPhiMatchedPart);
+
   histName = "fhnMassResponse";
   histTitle = Form("%s;#it{M}_{det};#it{M}_{part};#it{p}_{T,det};#it{p}_{T,part};#it{M}_{det}^{tagged}",histName.Data());
   fhnMassResponse = new THnSparseF(histName.Data(),histTitle.Data(),nBinsSparse0,nBins0,xmin0,xmax0);
@@ -171,7 +185,6 @@ void AliAnalysisTaskJetMassResponseDet::UserCreateOutputObjects()
 Bool_t AliAnalysisTaskJetMassResponseDet::Run()
 {
   // Run analysis code here, if needed. It will be executed before FillHistograms().
-
   return kTRUE;
 }
 
@@ -211,10 +224,16 @@ Bool_t AliAnalysisTaskJetMassResponseDet::FillHistograms()
        //fill detector response
        jPart = jDet->ClosestJet();
        if(jPart) {
-        AliEmcalJet *jDetT = jDet->GetTaggedJet();
-        Double_t mdetT = 0.;
-        if(jDetT) mdetT = jDetT->M();
-        Double_t var[5] = {GetJetMass(jDet),jPart->M(),jDet->Pt(),jPart->Pt(),mdetT};
+
+        fh2EtaPhiMatchedDet->Fill(jDet->Eta(),jDet->Phi());
+        fh2EtaPhiMatchedPart->Fill(jPart->Eta(),jPart->Phi());
+
+        Int_t nConstPart = jPart->GetNumberOfConstituents();
+        Int_t nConstDet = jDet->GetNumberOfConstituents();
+        Int_t diff = nConstPart-nConstDet;
+        Double_t eff = -1.;
+        if(nConstPart>0) eff = (Double_t)nConstDet/((Double_t)nConstPart);
+        Double_t var[5] = {GetJetMass(jDet),jPart->M(),jDet->Pt(),jPart->Pt(),eff};
         fhnMassResponse->Fill(var);
 
         if(jPart->Pt()>40. && jPart->Pt()<50.) {
@@ -223,9 +242,7 @@ Bool_t AliAnalysisTaskJetMassResponseDet::FillHistograms()
           else Printf("correct");
           Printf("pT Part: %f Det: %f",jPart->Pt(),jDet->Pt());
           Printf("mass Part: %f Det: %f",jPart->M(),jDet->M());
-          Int_t nConstPart = jPart->GetNumberOfConstituents();
-          Int_t nConstDet = jDet->GetNumberOfConstituents();
-          Int_t diff = nConstPart-nConstDet;
+
           Printf("nConst Part: %d  Det: %d  diff: %d",nConstPart,nConstDet,diff);
         }
        }
index d9d5edf94200c6e5499c7c8c415c084af79fb67f..b91f65f46a76e641dba77bcda0256ba706e2d8e9 100644 (file)
@@ -49,13 +49,15 @@ class AliAnalysisTaskJetMassResponseDet : public AliAnalysisTaskEmcalJet {
   TH2F            *fh2PtVsMassJetPartTaggedMatch;    //!pT vs mass of tagged particle level jets matched to a detector level jet
   TH2F            *fh2PtVsMassJetDetAll;             //!pT vs mass of all detector level jets
   TH2F            *fh2PtVsMassJetDetTagged;          //!pT vs mass of tagged detector level jets
+  TH2F            *fh2EtaPhiMatchedDet;              //!eta,phi of matched detector level jets
+  TH2F            *fh2EtaPhiMatchedPart;             //!eta,phi of matched particle level jets
   THnSparse       *fhnMassResponse;                  //!response matrix
 
  private:
   AliAnalysisTaskJetMassResponseDet(const AliAnalysisTaskJetMassResponseDet&);            // not implemented
   AliAnalysisTaskJetMassResponseDet &operator=(const AliAnalysisTaskJetMassResponseDet&); // not implemented
 
-  ClassDef(AliAnalysisTaskJetMassResponseDet, 1)
+  ClassDef(AliAnalysisTaskJetMassResponseDet, 2)
 };
 #endif
 
diff --git a/PWGJE/EMCALJetTasks/macros/AddTaskFastSimulation.C b/PWGJE/EMCALJetTasks/macros/AddTaskFastSimulation.C
new file mode 100644 (file)
index 0000000..5979353
--- /dev/null
@@ -0,0 +1,52 @@
+// $Id$
+
+AliJetFastSimulation* AddTaskFastSimulation(
+                                           const char     *tracksName1   = "Tracks",
+                                           const char     *tracksName2   = "Tracks2",
+                                           const char     *taskName      = "JetFastSimulation"
+                                           )
+{  
+  // Get the pointer to the existing analysis manager via the static access method.
+  //==============================================================================
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr)
+  {
+    ::Error("AddTaskMergeBranches", "No analysis manager to connect to.");
+    return NULL;
+  }  
+  
+  // Check the analysis type using the event handlers connected to the analysis manager.
+  //==============================================================================
+  if (!mgr->GetInputEventHandler())
+  {
+    ::Error("AddTaskMergeBranches", "This task requires an input event handler");
+    return NULL;
+  }
+  
+  //-------------------------------------------------------
+  // Init the task and do settings
+  //-------------------------------------------------------
+
+  AliJetFastSimulation *jetFastSim = new AliJetFastSimulation(taskName);
+  jetFastSim->SetTracksName(tracksName1);
+  jetFastSim->SetTracksOutName(tracksName2);
+
+  //-------------------------------------------------------
+  // Final settings, pass to manager and set the containers
+  //-------------------------------------------------------
+  mgr->AddTask(jetFastSim);
+    
+  // Create containers for input/output
+  mgr->ConnectInput (jetFastSim, 0, mgr->GetCommonInputContainer() );
+
+  TString contName = taskName;
+  contName += "_histos";
+  TString outputfile = Form("%s",AliAnalysisManager::GetCommonFileName());
+  AliAnalysisDataContainer *outc = mgr->CreateContainer(contName.Data(),
+                                                       TList::Class(),
+                                                       AliAnalysisManager::kOutputContainer,
+                                                       outputfile);
+  mgr->ConnectOutput(jetFastSim, 1, outc);
+
+  return jetFastSim;
+}
index 92ea7bf94ced7011a3deb6b2c2f9e35d580d9bbc..64310e59a42e1e2ba72db7551664edb951b1e50a 100644 (file)
@@ -30,7 +30,7 @@ AliAnalysisTaskJetMassResponseDet* AddTaskJetMassResponseDet(const char * njetsP
   AliAnalysisTaskJetMassResponseDet *task = new AliAnalysisTaskJetMassResponseDet(wagonName.Data());
 
   task->SetNCentBins(1);
-  task->SetVzRange(-10.,10.);
+  //  task->SetVzRange(-10.,10.);
 
   task->SetJetContainerPart(0);
   task->SetJetContainerDet(1);
index 0397dea971540c362350802cacccfd990e6c2a42..ce9534c9eea1d01abb1623e1fd81b2dd424b7ddc 100644 (file)
@@ -20,6 +20,7 @@
 #pragma link C++ class AliJetEmbeddingTask+;
 #pragma link C++ class AliJetEmbeddingFromGenTask+;
 #pragma link C++ class AliJetEmbeddingFromPYTHIATask+;
+#pragma link C++ class AliJetFastSimulation+;
 #pragma link C++ class AliJetModelBaseTask+;
 #pragma link C++ class AliJetModelMergeBranches+;
 #pragma link C++ class AliJetRandomizerTask+;