jet embedding response + average E

diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.cxx b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.cxx
index f92344f5e2cc338d93f294a29c9391900a61ff4e..78a1915e8cfe97d40c304ad3c3eec800bb6ab445 100644 (file)
--- a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.cxx
+++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.cxx
@@ -54,6 +54,10 @@ AliAnalysisTaskEmcalJetMassBkg::AliAnalysisTaskEmcalJetMassBkg() :
   fh2PtVsMassRCExLJDPhi(0),
   fh2PtVsMassPerpConeLJ(0),
   fh2PtVsMassPerpConeTJ(0),
+  fh2PtVsERC(0),
+  fh2PtVsERCExLJDPhi(0),
+  fh2PtVsEPerpConeLJ(0),
+  fh2PtVsEPerpConeTJ(0),
   fpPtVsMassRC(0),
   fpPtVsMassRCExLJ(0),
   fpPtVsMassPerpConeLJ(0),
@@ -90,6 +94,11 @@ AliAnalysisTaskEmcalJetMassBkg::AliAnalysisTaskEmcalJetMassBkg() :
   fh2PtVsMassPerpConeLJ    = new TH2F*[fNcentBins];
   fh2PtVsMassPerpConeTJ    = new TH2F*[fNcentBins];
 
+  fh2PtVsERC               = new TH2F*[fNcentBins];
+  fh2PtVsERCExLJDPhi       = new TH3F*[fNcentBins];
+  fh2PtVsEPerpConeLJ       = new TH2F*[fNcentBins];
+  fh2PtVsEPerpConeTJ       = new TH2F*[fNcentBins];
+
   fpPtVsMassRC             = new TProfile*[fNcentBins];
   fpPtVsMassRCExLJ         = new TProfile*[fNcentBins];
   fpPtVsMassPerpConeLJ     = new TProfile*[fNcentBins];
@@ -106,6 +115,11 @@ AliAnalysisTaskEmcalJetMassBkg::AliAnalysisTaskEmcalJetMassBkg() :
     fh2PtVsMassPerpConeLJ[i]     = 0;
     fh2PtVsMassPerpConeTJ[i]     = 0;
 
+    fh2PtVsERC[i]                = 0;
+    fh2PtVsERCExLJDPhi[i]        = 0;
+    fh2PtVsEPerpConeLJ[i]        = 0;
+    fh2PtVsEPerpConeTJ[i]        = 0;
+
     fpPtVsMassRC[i]              = 0;
     fpPtVsMassRCExLJ[i]          = 0;
     fpPtVsMassPerpConeLJ[i]      = 0;
@@ -139,6 +153,10 @@ AliAnalysisTaskEmcalJetMassBkg::AliAnalysisTaskEmcalJetMassBkg(const char *name)
   fh2PtVsMassRCExLJDPhi(0),
   fh2PtVsMassPerpConeLJ(0),
   fh2PtVsMassPerpConeTJ(0),
+  fh2PtVsERC(0),
+  fh2PtVsERCExLJDPhi(0),
+  fh2PtVsEPerpConeLJ(0),
+  fh2PtVsEPerpConeTJ(0),
   fpPtVsMassRC(0),
   fpPtVsMassRCExLJ(0),
   fpPtVsMassPerpConeLJ(0),
@@ -175,6 +193,11 @@ AliAnalysisTaskEmcalJetMassBkg::AliAnalysisTaskEmcalJetMassBkg(const char *name)
   fh2PtVsMassPerpConeLJ    = new TH2F*[fNcentBins];
   fh2PtVsMassPerpConeTJ    = new TH2F*[fNcentBins];
 
+  fh2PtVsERC               = new TH2F*[fNcentBins];
+  fh2PtVsERCExLJDPhi       = new TH3F*[fNcentBins];
+  fh2PtVsEPerpConeLJ       = new TH2F*[fNcentBins];
+  fh2PtVsEPerpConeTJ       = new TH2F*[fNcentBins];
+
   fpPtVsMassRC             = new TProfile*[fNcentBins];
   fpPtVsMassRCExLJ         = new TProfile*[fNcentBins];
   fpPtVsMassPerpConeLJ     = new TProfile*[fNcentBins];
@@ -191,6 +214,11 @@ AliAnalysisTaskEmcalJetMassBkg::AliAnalysisTaskEmcalJetMassBkg(const char *name)
     fh2PtVsMassPerpConeLJ[i]     = 0;
     fh2PtVsMassPerpConeTJ[i]     = 0;
 
+    fh2PtVsERC[i]                = 0;
+    fh2PtVsERCExLJDPhi[i]        = 0;
+    fh2PtVsEPerpConeLJ[i]        = 0;
+    fh2PtVsEPerpConeTJ[i]        = 0;
+
     fpPtVsMassRC[i]              = 0;
     fpPtVsMassRCExLJ[i]          = 0;
     fpPtVsMassPerpConeLJ[i]      = 0;
@@ -229,6 +257,10 @@ void AliAnalysisTaskEmcalJetMassBkg::UserCreateOutputObjects()
   const Double_t minPt = -50.;
   const Double_t maxPt = 200.;
 
+  const Int_t nBinsE  = 250;
+  const Double_t minE = -50.;
+  const Double_t maxE = 200.;
+
   const Int_t nBinsM  = 150;
   const Double_t minM = -50.;
   const Double_t maxM = 100.;
@@ -327,7 +359,27 @@ void AliAnalysisTaskEmcalJetMassBkg::UserCreateOutputObjects()
     histTitle = TString::Format("%s;#it{p}_{T,PerpConeTJ};#it{M}_{PerpConeTJ}",histName.Data());
     fh2PtVsMassPerpConeTJ[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
     fOutput->Add(fh2PtVsMassPerpConeTJ[i]);
+    //
+    histName = TString::Format("fh2PtVsERC_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,RC};#it{M}_{RC}",histName.Data());
+    fh2PtVsERC[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsE,minE,maxE);
+    fOutput->Add(fh2PtVsERC[i]);
 
+    histName = TString::Format("fh2PtVsERCExLJDPhi_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,RC};#it{M}_{RC}",histName.Data());
+    fh2PtVsERCExLJDPhi[i] = new TH3F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,72,-0.5*TMath::Pi(),1.5*TMath::Pi());
+    fOutput->Add(fh2PtVsERCExLJDPhi[i]);
+
+    histName = TString::Format("fh2PtVsEPerpConeLJ_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,PerpConeLJ};#it{M}_{PerpConeLJ}",histName.Data());
+    fh2PtVsEPerpConeLJ[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsE,minE,maxE);
+    fOutput->Add(fh2PtVsEPerpConeLJ[i]);
+
+    histName = TString::Format("fh2PtVsEPerpConeTJ_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,PerpConeTJ};#it{M}_{PerpConeTJ}",histName.Data());
+    fh2PtVsEPerpConeTJ[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsE,minE,maxE);
+    fOutput->Add(fh2PtVsEPerpConeTJ[i]);
+    //
     histName = TString::Format("fh2EtaVsMassRC_%d",i);
     histTitle = TString::Format("%s;#eta_{RC};#it{M}_{RC}",histName.Data());
     fh2EtaVsMassRC[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsEta,minEta,maxEta,nBinsM,minM,maxM);
@@ -414,6 +466,7 @@ Bool_t AliAnalysisTaskEmcalJetMassBkg::FillHistograms()
   Float_t RCeta = 0;
   Float_t RCphi = 0;
   Float_t RCmass = 0.;  
+  Float_t RCE = 0.;  
 
   static Double_t massvecRC[999];
   static Double_t massPerAreavecRC[999];
@@ -432,8 +485,10 @@ Bool_t AliAnalysisTaskEmcalJetMassBkg::FillHistograms()
     RCphi = 0;
     GetRandomCone(lvRC,RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, 0);
     RCmass = lvRC.M();
+    RCE = lvRC.E();
     if (RCpt > 0) {
       fh2PtVsMassRC[fCentBin]->Fill(RCpt - rho*rcArea,RCmass);
+      fh2PtVsERC[fCentBin]->Fill(RCpt - rho*rcArea,RCE);
       fpPtVsMassRC[fCentBin]->Fill(RCpt - rho*rcArea,RCmass);
       fh2EtaVsMassRC[fCentBin]->Fill(RCeta,RCmass);
       fh2CentVsMassRC->Fill(fCent,RCmass);
@@ -452,11 +507,13 @@ Bool_t AliAnalysisTaskEmcalJetMassBkg::FillHistograms()
       RCphi = 0;
       GetRandomCone(lvRC,RCpt, RCeta, RCphi, fTracksCont, fCaloClustersCont, jet);
       RCmass = lvRC.M();
+      RCE = lvRC.E();
       if (RCpt > 0 && jet) {
        Float_t dphi = RCphi - jet->Phi();
        if (dphi > 1.5*TMath::Pi()) dphi -= TMath::Pi() * 2;
        if (dphi < -0.5*TMath::Pi()) dphi += TMath::Pi() * 2;
        fh2PtVsMassRCExLJDPhi[fCentBin]->Fill(RCpt - rho*rcArea,RCmass,dphi);
+       fh2PtVsERCExLJDPhi[fCentBin]->Fill(RCpt - rho*rcArea,RCE,dphi);
        fpPtVsMassRCExLJ[fCentBin]->Fill(RCpt - rho*rcArea,RCmass);
        fh2EtaVsMassRCExLJ[fCentBin]->Fill(RCeta,RCmass);
        fh2CentVsMassRCExLJ->Fill(fCent,RCmass);
@@ -506,12 +563,15 @@ Bool_t AliAnalysisTaskEmcalJetMassBkg::FillHistograms()
     Float_t PCpt = 0;
     Float_t PCeta = 0;
     Float_t PCphi = 0;
-    Float_t PCmass = 0.;  
+    Float_t PCmass = 0.;
+    Float_t PCE = 0.;
     if(jet) {
       GetPerpCone(lvPC,PCpt, PCeta, PCphi, fTracksCont, fCaloClustersCont, jet);
       PCmass = lvPC.M();
+      PCE = lvPC.E();
       if(PCpt>0.) {
        fh2PtVsMassPerpConeLJ[fCentBin]->Fill(PCpt-rho*rcArea,PCmass);
+       fh2PtVsEPerpConeLJ[fCentBin]->Fill(PCpt-rho*rcArea,PCE);
        fpPtVsMassPerpConeLJ[fCentBin]->Fill(PCpt-rho*rcArea,PCmass);
        fh2EtaVsMassPerpConeLJ[fCentBin]->Fill(PCeta,PCmass);
        fh2CentVsMassPerpConeLJ->Fill(fCent,PCmass);
@@ -532,8 +592,10 @@ Bool_t AliAnalysisTaskEmcalJetMassBkg::FillHistograms()
       PCphi = 0;
       GetPerpCone(lvPC,PCpt, PCeta, PCphi, fTracksCont, fCaloClustersCont, jet);
       PCmass = lvPC.M();
+      PCE = lvPC.E();
       if(PCpt>0.) {
        fh2PtVsMassPerpConeTJ[fCentBin]->Fill(PCpt-rho*rcArea,PCmass);
+       fh2PtVsEPerpConeTJ[fCentBin]->Fill(PCpt-rho*rcArea,PCE);
        fpPtVsMassPerpConeTJ[fCentBin]->Fill(PCpt-rho*rcArea,PCmass);
        fh2EtaVsMassPerpConeTJ[fCentBin]->Fill(PCeta,PCmass);
        fh2CentVsMassPerpConeTJ->Fill(fCent,PCmass);

diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.h b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.h
index a3d90938587df5d2485310632c466e510685a946..23bad2b1232b2c769f4387da1384cb18010014b3 100644 (file)
--- a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.h
+++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassBkg.h
@@ -64,6 +64,11 @@ class AliAnalysisTaskEmcalJetMassBkg : public AliAnalysisTaskEmcalJet {
   TH2F            **fh2PtVsMassPerpConeLJ;        //!pT vs mass of cone perpendicular to leading jet
   TH2F            **fh2PtVsMassPerpConeTJ;        //!pT vs mass of cone perpendicular to all tagged jets
 
+  TH2F            **fh2PtVsERC;                   //!E vs mass of RC
+  TH3F            **fh2PtVsERCExLJDPhi;           //!E vs mass of RC excluding area around leading jet
+  TH2F            **fh2PtVsEPerpConeLJ;           //!E vs mass of cone perpendicular to leading jet
+  TH2F            **fh2PtVsEPerpConeTJ;           //!E vs mass of cone perpendicular to all tagged jets
+
   TProfile        **fpPtVsMassRC;                 //!pT vs Avg mass of RC
   TProfile        **fpPtVsMassRCExLJ;             //!pT vs Avg mass of RC excluding area around leading jet
   TProfile        **fpPtVsMassPerpConeLJ;         //!pT vs Avg mass of cone perpendicular to leading jet
@@ -105,7 +110,7 @@ class AliAnalysisTaskEmcalJetMassBkg : public AliAnalysisTaskEmcalJet {
   AliAnalysisTaskEmcalJetMassBkg(const AliAnalysisTaskEmcalJetMassBkg&);            // not implemented
   AliAnalysisTaskEmcalJetMassBkg &operator=(const AliAnalysisTaskEmcalJetMassBkg&); // not implemented
 
-  ClassDef(AliAnalysisTaskEmcalJetMassBkg, 3)
+  ClassDef(AliAnalysisTaskEmcalJetMassBkg, 4)
 };
 #endif
 

diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.cxx b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.cxx
index 2a30d827e67c75124776f0aab511a67f9fb455d7..724185a180faa7744ee02b369f9bf66b49b773d9 100644 (file)
--- a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.cxx
+++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.cxx
@@ -8,6 +8,7 @@
 #include <TH2F.h>
 #include <TH3F.h>
 #include <THnSparse.h>
+#include <TF1.h>
 #include <TList.h>
 #include <TLorentzVector.h>
 #include <TProfile.h>
@@ -40,24 +41,39 @@ AliAnalysisTaskEmcalJetMassResponse::AliAnalysisTaskEmcalJetMassResponse() :
   AliAnalysisTaskEmcalJet("AliAnalysisTaskEmcalJetMassResponse", kTRUE),
   fContainerBase(0),
   fMinFractionShared(0),
-  fJetMassAvg(0),
+  f1JetMassAvg(0),
+  fh2PtJet1DeltaMNoSub(0),
+  fh2PtJet2DeltaMNoSub(0),
+  fh3PtJet1DeltaPtDeltaMCheat(0),
+  fh3PtJet2DeltaPtDeltaMCheat(0),
   fh3PtJet1DeltaPtDeltaM(0),
   fh3PtJet2DeltaPtDeltaM(0),
   fh3PtJet1MJet1MJet2(0),
-  fh3PtJet2MJet1MJet2(0)
+  fh3PtJet2MJet1MJet2(0),
+  fh2PtJet1DeltaPtVecSub(0)
 {
   // Default constructor.
 
+  fh2PtJet1DeltaMNoSub         = new TH2F*[fNcentBins];
+  fh2PtJet2DeltaMNoSub         = new TH2F*[fNcentBins];
+  fh3PtJet1DeltaPtDeltaMCheat  = new TH3F*[fNcentBins];
+  fh3PtJet2DeltaPtDeltaMCheat  = new TH3F*[fNcentBins];
   fh3PtJet1DeltaPtDeltaM       = new TH3F*[fNcentBins];
   fh3PtJet2DeltaPtDeltaM       = new TH3F*[fNcentBins];
   fh3PtJet1MJet1MJet2          = new TH3F*[fNcentBins];
   fh3PtJet2MJet1MJet2          = new TH3F*[fNcentBins];
+  fh2PtJet1DeltaPtVecSub       = new TH2F*[fNcentBins];
  
   for (Int_t i = 0; i < fNcentBins; i++) {
-    fh3PtJet1DeltaPtDeltaM[i]     = 0; 
-    fh3PtJet2DeltaPtDeltaM[i]     = 0;
-    fh3PtJet1MJet1MJet2[i]        = 0;
-    fh3PtJet2MJet1MJet2[i]        = 0;
+    fh2PtJet1DeltaMNoSub[i]        = 0;
+    fh2PtJet2DeltaMNoSub[i]        = 0;
+    fh3PtJet1DeltaPtDeltaMCheat[i] = 0;
+    fh3PtJet2DeltaPtDeltaMCheat[i] = 0;
+    fh3PtJet1DeltaPtDeltaM[i]      = 0; 
+    fh3PtJet2DeltaPtDeltaM[i]      = 0;
+    fh3PtJet1MJet1MJet2[i]         = 0;
+    fh3PtJet2MJet1MJet2[i]         = 0;
+    fh2PtJet1DeltaPtVecSub[i]      = 0;
   }
   SetMakeGeneralHistograms(kTRUE);
 }
@@ -67,24 +83,39 @@ AliAnalysisTaskEmcalJetMassResponse::AliAnalysisTaskEmcalJetMassResponse(const c
   AliAnalysisTaskEmcalJet(name, kTRUE),  
   fContainerBase(0),
   fMinFractionShared(0),
-  fJetMassAvg(0),
+  f1JetMassAvg(0),
+  fh2PtJet1DeltaMNoSub(0),
+  fh2PtJet2DeltaMNoSub(0),
+  fh3PtJet1DeltaPtDeltaMCheat(0),
+  fh3PtJet2DeltaPtDeltaMCheat(0),
   fh3PtJet1DeltaPtDeltaM(0),
   fh3PtJet2DeltaPtDeltaM(0),
   fh3PtJet1MJet1MJet2(0),
-  fh3PtJet2MJet1MJet2(0)
+  fh3PtJet2MJet1MJet2(0),
+  fh2PtJet1DeltaPtVecSub(0)
 {
   // Standard constructor.
 
+  fh2PtJet1DeltaMNoSub         = new TH2F*[fNcentBins];
+  fh2PtJet2DeltaMNoSub         = new TH2F*[fNcentBins];
+  fh3PtJet1DeltaPtDeltaMCheat  = new TH3F*[fNcentBins];
+  fh3PtJet2DeltaPtDeltaMCheat  = new TH3F*[fNcentBins];
   fh3PtJet1DeltaPtDeltaM       = new TH3F*[fNcentBins];
   fh3PtJet2DeltaPtDeltaM       = new TH3F*[fNcentBins];
   fh3PtJet1MJet1MJet2          = new TH3F*[fNcentBins];
   fh3PtJet2MJet1MJet2          = new TH3F*[fNcentBins];
+  fh2PtJet1DeltaPtVecSub       = new TH2F*[fNcentBins];
  
   for (Int_t i = 0; i < fNcentBins; i++) {
-    fh3PtJet1DeltaPtDeltaM[i]     = 0; 
-    fh3PtJet2DeltaPtDeltaM[i]     = 0; 
-    fh3PtJet1MJet1MJet2[i]        = 0;
-    fh3PtJet2MJet1MJet2[i]        = 0;
+    fh2PtJet1DeltaMNoSub[i]        = 0;
+    fh2PtJet2DeltaMNoSub[i]        = 0;
+    fh3PtJet1DeltaPtDeltaMCheat[i] = 0;
+    fh3PtJet2DeltaPtDeltaMCheat[i] = 0;
+    fh3PtJet1DeltaPtDeltaM[i]      = 0; 
+    fh3PtJet2DeltaPtDeltaM[i]      = 0;
+    fh3PtJet1MJet1MJet2[i]         = 0;
+    fh3PtJet2MJet1MJet2[i]         = 0;
+    fh2PtJet1DeltaPtVecSub[i]      = 0;
   }
 
   SetMakeGeneralHistograms(kTRUE);
@@ -106,25 +137,37 @@ void AliAnalysisTaskEmcalJetMassResponse::UserCreateOutputObjects()
   Bool_t oldStatus = TH1::AddDirectoryStatus();
   TH1::AddDirectory(kFALSE);
 
-  const Int_t nBinsPt  = 250;
+  const Int_t nBinsPt  = 200;
   const Double_t minPt = -50.;
-  const Double_t maxPt = 200.;
+  const Double_t maxPt = 150.;
 
   const Int_t nBinsM  = 150;
   const Double_t minM = -50.;
   const Double_t maxM = 100.;
 
-  // const Int_t nBinsArea = 100;
-  // const Double_t minArea = 0.;
-  // const Double_t maxArea = 1.;
-
-  // const Int_t nBinsNConst = 100;
-  // const Double_t minNConst = 0.;
-  // const Double_t maxNConst = 500.;
-
   TString histName = "";
   TString histTitle = "";
   for (Int_t i = 0; i < fNcentBins; i++) {
+    histName = TString::Format("fh2PtJet1DeltaMNoSub_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,jet1};#delta#it{M}_{jet}",histName.Data());
+    fh2PtJet1DeltaMNoSub[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
+    fOutput->Add(fh2PtJet1DeltaMNoSub[i]);
+
+    histName = TString::Format("fh2PtJet2DeltaMNoSub_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,jet1};#delta#it{M}_{jet}",histName.Data());
+    fh2PtJet2DeltaMNoSub[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
+    fOutput->Add(fh2PtJet2DeltaMNoSub[i]);
+
+    histName = TString::Format("fh3PtJet1DeltaPtDeltaMCheat_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,jet1};#delta#it{p}_{T};#delta#it{M}_{jet}",histName.Data());
+    fh3PtJet1DeltaPtDeltaMCheat[i] = new TH3F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
+    fOutput->Add(fh3PtJet1DeltaPtDeltaMCheat[i]);
+
+    histName = TString::Format("fh3PtJet2DeltaPtDeltaMCheat_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,jet1};#delta#it{p}_{T};#delta#it{M}_{jet}",histName.Data());
+    fh3PtJet2DeltaPtDeltaMCheat[i] = new TH3F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
+    fOutput->Add(fh3PtJet2DeltaPtDeltaMCheat[i]);
+
     histName = TString::Format("fh3PtJet1DeltaPtDeltaM_%d",i);
     histTitle = TString::Format("%s;#it{p}_{T,jet1};#delta#it{p}_{T};#delta#it{M}_{jet}",histName.Data());
     fh3PtJet1DeltaPtDeltaM[i] = new TH3F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt,nBinsM,minM,maxM);
@@ -144,6 +187,11 @@ void AliAnalysisTaskEmcalJetMassResponse::UserCreateOutputObjects()
     histTitle = TString::Format("%s;#it{p}_{T,jet2};#it{M}_{jet1};#it{M}_{jet2}",histName.Data());
     fh3PtJet2MJet1MJet2[i] = new TH3F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsM,minM,maxM,nBinsM,minM,maxM);
     fOutput->Add(fh3PtJet2MJet1MJet2[i]);
+
+    histName = TString::Format("fh2PtJet1DeltaPtVecSub_%d",i);
+    histTitle = TString::Format("%s;#it{p}_{T,jet1};#delta#it{p}_{T}",histName.Data());
+    fh2PtJet1DeltaPtVecSub[i] = new TH2F(histName.Data(),histTitle.Data(),nBinsPt,minPt,maxPt,nBinsPt,minPt,maxPt);
+    fOutput->Add(fh2PtJet1DeltaPtVecSub[i]);
   }
 
   // =========== Switch on Sumw2 for all histos ===========
@@ -197,11 +245,18 @@ Bool_t AliAnalysisTaskEmcalJetMassResponse::FillHistograms()
       Double_t deltaPt = ptJet1 - ptJet2;
       Double_t deltaM  = massJet1 - massJet2;
 
+      fh2PtJet1DeltaMNoSub[fCentBin]->Fill(ptJet1,jet1->M()-jet2->M());
+      fh2PtJet2DeltaMNoSub[fCentBin]->Fill(ptJet2,jet1->M()-jet2->M());
+
       fh3PtJet1DeltaPtDeltaM[fCentBin]->Fill(ptJet1,deltaPt,deltaM);
       fh3PtJet2DeltaPtDeltaM[fCentBin]->Fill(ptJet2,deltaPt,deltaM);
 
       fh3PtJet1MJet1MJet2[fCentBin]->Fill(ptJet1,massJet1,massJet2);
       fh3PtJet2MJet1MJet2[fCentBin]->Fill(ptJet2,massJet1,massJet2);
+
+      TLorentzVector vpC = GetSubtractedVectorCheat(jet1);
+      fh3PtJet1DeltaPtDeltaMCheat[fCentBin]->Fill(vpC.Perp(),vpC.Perp()-jet2->Pt(),vpC.M()-jet2->M());
+      fh3PtJet2DeltaPtDeltaMCheat[fCentBin]->Fill(ptJet2,vpC.Perp()-jet2->Pt(),vpC.M()-jet2->M());
     }
   }
 
@@ -209,22 +264,53 @@ Bool_t AliAnalysisTaskEmcalJetMassResponse::FillHistograms()
 }
 
 //________________________________________________________________________
-Double_t AliAnalysisTaskEmcalJetMassResponse::GetJetMass(AliEmcalJet *jet) {
+TLorentzVector AliAnalysisTaskEmcalJetMassResponse::GetSubtractedVector(AliEmcalJet *jet) {
+  //get subtracted vector
+  TLorentzVector vpS;
+  if(f1JetMassAvg) {
+    Double_t pt = jet->Pt() - GetRhoVal(fContainerBase)*jet->Area();
+    TLorentzVector vpB; vpB.SetPtEtaPhiE(GetRhoVal(fContainerBase)*jet->Area(),0.,0.,f1JetMassAvg->Eval(pt));
+    TLorentzVector vpAAboost; vpAAboost.SetPtEtaPhiM(jet->Pt(),0.,0.,jet->M());
+    TLorentzVector vpSboost = vpAAboost - vpB;
+    vpS.SetPtEtaPhiM(vpSboost.Perp(),jet->Eta(),jet->Phi(),vpSboost.M());
+  }
+  return vpS;
+}
 
-  Double_t subM = jet->M() - fJetMassAvg;
-  // Double_t scale = deltaM / jet->M();
+//________________________________________________________________________
+TLorentzVector AliAnalysisTaskEmcalJetMassResponse::GetSubtractedVectorCheat(AliEmcalJet *jet) {
+  //get subtracted vector taking pT and mass difference from MC match
+  TLorentzVector vpS;
+  AliEmcalJet *jet2 = jet->ClosestJet();
+  if(jet2) {
+    TLorentzVector vpAAboost; vpAAboost.SetPtEtaPhiM(jet->Pt(),0.,0.,jet->M());
+    TLorentzVector vpPPboost; vpPPboost.SetPtEtaPhiM(jet2->Pt(),0.,0.,jet2->M());
+    Double_t dpt = vpAAboost.Perp()-vpPPboost.Perp();
+    /*
+      Double_t dm = vpAAboost.M() - vpPPboost.M();
+      Double_t dE = TMath::Sqrt(dpt*dpt + dm*dm);
+    */
+    Double_t dE = vpAAboost.E()-vpPPboost.E();
+    TLorentzVector vpB; vpB.SetPtEtaPhiE(dpt,0.,0.,dE);
+    TLorentzVector vpSboost = vpAAboost - vpB;
+    vpS.SetPtEtaPhiM(vpSboost.Perp(),jet->Eta(),jet->Phi(),vpSboost.M());
+  }
+  return vpS;
+}
 
-  // Double_t mt2 = jet->E()*jet->E() - jet->Pt()*jet->Pt();
-  // Double_t et2 = jet->M()*jet->M() + jet->Pt()*jet->Pt();
+//________________________________________________________________________
+Double_t AliAnalysisTaskEmcalJetMassResponse::GetJetMass(AliEmcalJet *jet) {
 
-  // Double_t pxScale = jet->Px()*scale;
-  // Double_t pyScale = jet->Py()*scale;
-  // Double_t pzScale = jet->Pz()*scale;
+  if(f1JetMassAvg) {
+    TLorentzVector vpS = GetSubtractedVector(jet);
   
-  // Printf("scaled jet 4-vector: %f-%f-%f-%f",pxScale,pyScale,pzScale);
-
-  return subM;
+    AliEmcalJet *jet2 = jet->ClosestJet();
+    if(jet2) fh2PtJet1DeltaPtVecSub[fCentBin]->Fill(vpS.Perp(),vpS.Perp()-jet2->Pt());
 
+    return vpS.M();
+  }
+  else
+    return jet->M();
 }
 
 //________________________________________________________________________

diff --git a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.h b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.h
index cf1e0ffad234655339f7cf81e36a083226ae1d04..f6a45601be3c37b121e548d07ece813d34b39a1b 100644 (file)
--- a/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.h
+++ b/PWGJE/EMCALJetTasks/UserTasks/AliAnalysisTaskEmcalJetMassResponse.h
@@ -6,6 +6,8 @@ class TH2;
 class TH3;
 class TH3F;
 class THnSparse;
+class TF1;
+class TLorentzVector;
 class TClonesArray;
 class TArrayI;
 class AliAnalysisManager;
@@ -26,24 +28,34 @@ class AliAnalysisTaskEmcalJetMassResponse : public AliAnalysisTaskEmcalJet {
   //Setters
   void SetJetContainerBase(Int_t c)                             { fContainerBase     = c   ; }
   void SetMinFractionShared(Double_t f)                         { fMinFractionShared = f   ; }
-  void SetJetMassAverage(Double_t m)                            { fJetMassAvg        = m   ; }
+  void SetJetMassAverageFunc(TF1 *f)                            { f1JetMassAvg       = f   ; }
 
  protected:
   Bool_t                              RetrieveEventObjects();
   Bool_t                              Run();
   Bool_t                              FillHistograms();
 
+  TLorentzVector                      GetSubtractedVector(AliEmcalJet *jet);
+  TLorentzVector                      GetSubtractedVectorCheat(AliEmcalJet *jet);
   Double_t                            GetJetMass(AliEmcalJet *jet);
 
   Int_t                               fContainerBase;              // jets to be analyzed
   Double_t                            fMinFractionShared;          // only fill histos for jets if shared fraction larger than X
-  Double_t                            fJetMassAvg;                 // average jet mass
+  TF1                                *f1JetMassAvg;                // parametrization of average jet mass
+
+  TH2F            **fh2PtJet1DeltaMNoSub;                          //!pt jet1 vs delta-pt vs delta-M
+  TH2F            **fh2PtJet2DeltaMNoSub;                          //!pt jet2 vs delta-pt vs delta-M
+
+  TH3F            **fh3PtJet1DeltaPtDeltaMCheat;                   //!pt jet1 vs delta-pt vs delta-M
+  TH3F            **fh3PtJet2DeltaPtDeltaMCheat;                   //!pt jet2 vs delta-pt vs delta-M
 
   TH3F            **fh3PtJet1DeltaPtDeltaM;                        //!pt jet1 vs delta-pt vs delta-M
   TH3F            **fh3PtJet2DeltaPtDeltaM;                        //!pt jet2 vs delta-pt vs delta-M
   TH3F            **fh3PtJet1MJet1MJet2;                           //!pt jet1 vs jet mass jet1 vs jet mass jet2
   TH3F            **fh3PtJet2MJet1MJet2;                           //!pt jet2 vs jet mass jet1 vs jet mass jet2
 
+  TH2F            **fh2PtJet1DeltaPtVecSub;                        //!pt jet1 (AA) vs delta pT while using vector subtraction
+
  private:
   AliAnalysisTaskEmcalJetMassResponse(const AliAnalysisTaskEmcalJetMassResponse&);            // not implemented
   AliAnalysisTaskEmcalJetMassResponse &operator=(const AliAnalysisTaskEmcalJetMassResponse&); // not implemented