Intermediate stage of implementing corrections in AliAnalysisHadEtReconstructed ...
authorcnattras <cnattras@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 23 Sep 2010 02:24:35 +0000 (02:24 +0000)
committercnattras <cnattras@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 23 Sep 2010 02:24:35 +0000 (02:24 +0000)
PWG4/totEt/AliAnalysisHadEt.cxx
PWG4/totEt/AliAnalysisHadEt.h
PWG4/totEt/AliAnalysisHadEtCorrections.cxx
PWG4/totEt/AliAnalysisHadEtCorrections.h
PWG4/totEt/AliAnalysisHadEtMonteCarlo.cxx
PWG4/totEt/AliAnalysisHadEtMonteCarlo.h
PWG4/totEt/AliAnalysisHadEtReconstructed.cxx
PWG4/totEt/AliAnalysisHadEtReconstructed.h
PWG4/totEt/macros/ConfigHadEtAnalysis.C [new file with mode: 0644]
PWG4/totEt/macros/CreateAlienHandlerHadEtSim.C
PWG4/totEt/macros/runHadEt.C

index d9a4420..5b79d7e 100644 (file)
@@ -8,6 +8,7 @@
 // it has daughters, AliAnalysisHadEtMonteCarlo and 
 // AliAnalysisHadEtReconstructed which loop over either Monte Carlo data or 
 // real data to get Et
+
 #include "AliAnalysisHadEt.h"
 #include "TMath.h"
 #include "TList.h"
@@ -72,6 +73,9 @@ AliAnalysisHadEt::AliAnalysisHadEt() :
        ,fEtaCode(0)
        ,fOmega0Code(0)
         ,fPionMass(0)
+        ,fKaonMass(0)
+        ,fProtonMass(0)
+        ,fElectronMass(0)
        ,fSumEt(0)
        ,fSumEtAcc(0)
        ,fTotEt(0)
@@ -110,6 +114,7 @@ void AliAnalysisHadEt::FillOutputList()
 void AliAnalysisHadEt::Init()
 {// clear variables, set up cuts and PDG info
   ResetEventValues();
+
 }
 
 void AliAnalysisHadEt::CreateHistograms()
@@ -152,6 +157,9 @@ void AliAnalysisHadEt::ResetEventValues()
 void AliAnalysisHadEt::SetParticleCodes()
 {  //the codes are defined in $ROOTSYS/etc/pdg_table.txt
     fPionMass = fPdgDB->GetParticle("pi+")->Mass();
+    fKaonMass = fPdgDB->GetParticle("K+")->Mass();
+    fProtonMass = fPdgDB->GetParticle("proton")->Mass();
+    fElectronMass = fPdgDB->GetParticle("e+")->Mass();
     fPiPlusCode = fPdgDB->GetParticle("pi+")->PdgCode();
     fPiMinusCode = fPdgDB->GetParticle("pi-")->PdgCode();
     fKPlusCode = fPdgDB->GetParticle("K+")->PdgCode();
@@ -352,3 +360,24 @@ Float_t AliAnalysisHadEt::Et(TParticle *part, float mass){//function to calculat
   }
   return 0.0;
 }
+Float_t AliAnalysisHadEt::Et(Float_t p, Float_t theta, Int_t pid, Short_t charge){//function to calculate et in the same way as it would be calculated in a calorimeter
+  if(pid==fPiPlusCode || pid==fPiMinusCode){//Nothing special for pions
+    return TMath::Sqrt(p*p + fPionMass*fPionMass) * TMath::Sin(theta);
+  }
+  if(pid==fKPlusCode || pid==fKMinusCode){//Nothing special for kaons
+    return TMath::Sqrt(p*p + fKaonMass*fKaonMass) * TMath::Sin(theta);
+  }
+  if(pid==fEPlusCode || pid==fEMinusCode){//Nothing special for electrons
+    return TMath::Sqrt(p*p + fElectronMass*fElectronMass) * TMath::Sin(theta);
+  }
+  if(pid==fProtonCode || pid==fAntiProtonCode){//But for protons we must be careful...
+    if(charge<0.0){//antiprotns: kinetic energy plus twice the rest mass
+      return (TMath::Sqrt(p*p + fProtonMass*fProtonMass) + fProtonMass) * TMath::Sin(theta);
+    }
+    if(charge>0.0){//antiprotns: kinetic energy only
+      return (TMath::Sqrt(p*p + fProtonMass*fProtonMass) - fProtonMass) * TMath::Sin(theta);
+    }
+  }
+  cerr<<"Uh-oh!  Et not set properly!"<<endl;
+  return 0.0;
+}
index 0fd33dc..3bc4c53 100644 (file)
@@ -123,6 +123,9 @@ protected:
     Int_t fEtaCode;//pdg eta code
     Int_t fOmega0Code;//pdg eta code
     Float_t fPionMass;//pdg pion mass
+    Float_t fKaonMass;//pdg kaon mass
+    Float_t fProtonMass;//pdg proton mass
+    Float_t fElectronMass;//pdg electron mass
 
     
     Double_t fSumEt;/** Sum of the total Et for all events */
@@ -149,6 +152,7 @@ protected:
     void FillHisto2D(TString histname, Float_t x, Float_t y, Float_t weight);
 
     Float_t Et(TParticle *part, float mass = -1000);
+    Float_t Et(Float_t p, Float_t theta, Int_t pid, Short_t charge);
     AliESDtrackCuts* fEsdtrackCutsITSTPC;//esd track cuts for ITS+TPC tracks
     AliESDtrackCuts* fEsdtrackCutsTPC;//esd track cuts for TPC tracks (which may also contain ITS hits)
     AliESDtrackCuts* fEsdtrackCutsITS;//esd track cuts for ITS stand alone tracks
@@ -159,6 +163,8 @@ protected:
     static Float_t fgPtAxis[117];//bins for pt axis of histograms
     static Int_t fgNumOfPtBins;//number of pt bins
     
+
+
  private:
     //Declare it private to avoid compilation warning
     AliAnalysisHadEt & operator = (const AliAnalysisHadEt & g) ;//cpy assignment
index ea4ded5..e7875b4 100644 (file)
@@ -130,6 +130,50 @@ AliAnalysisHadEtCorrections::AliAnalysisHadEtCorrections(const AliAnalysisHadEtC
   fBackgroundITS = new TH1D(*(g->fBackgroundITS));
 }
 
+
+Float_t AliAnalysisHadEtCorrections::GetConstantCorrections(Bool_t totEt, Float_t ptcut, TString type){
+  Float_t acceptance = 0.0;
+  Float_t neutral = 0.0;
+  Float_t ptcorr = 0.0;
+  float correction = 0.0;
+
+  //TString *type = new TString(mytype);
+
+  if(type.Contains("Full")) acceptance = fAcceptanceCorrectionFull;
+  if(type.Contains("EMCAL")) acceptance = fAcceptanceCorrectionEMCAL;
+  if(type.Contains("PHOS")) acceptance = fAcceptanceCorrectionPHOS;
+
+  if(type.Contains("High")){//high bound
+    if(totEt) neutral = fNotHadronicCorrectionHigh;
+    else{neutral = fNeutralCorrectionHigh;}
+    if(ptcut>0.12){ptcorr = ffpTcutCorrectionTPCHigh;}
+    else{ptcorr = ffpTcutCorrectionITSHigh;}
+    cout<<"Setting correction factor to "<<correction<<endl;
+    return correction;
+  }
+  if(type.Contains("Low")){//high bound
+    if(totEt) neutral = fNotHadronicCorrectionLow;
+    else{neutral = fNeutralCorrectionLow;}
+    if(ptcut>0.12){ptcorr = ffpTcutCorrectionTPCLow;}
+    else{ptcorr = ffpTcutCorrectionITSLow;}
+    cout<<"Setting correction factor to "<<correction<<endl;
+    return correction;
+  }
+
+  if(totEt) neutral = fNotHadronicCorrection;
+  else{neutral = fNeutralCorrection;}
+  if(ptcut>0.12){ptcorr = fpTcutCorrectionTPC;}
+  else{ptcorr = fpTcutCorrectionITS;}
+
+  correction = acceptance*neutral*ptcorr;
+  cout<<"Setting correction factor for ";
+  if(totEt) cout<<"total et";
+  else{cout<<"hadronic et";}
+  cout<<" with the pt cut off "<<ptcut<<" for "<<type<<" acceptance to "<<correction<<endl;
+  //cout<<"Acceptance "<<acceptance<<" neutral "<<neutral<<" ptcorr "<<ptcorr<<endl;
+  return correction;
+
+}
 // AliAnalysisHadEtCorrections & operator = (const AliAnalysisHadEtCorrections & g) {
 
 //   fEtaCut=g->fEtaCut;
index 94fb2af..30b8344 100644 (file)
@@ -66,8 +66,8 @@ public:
     Float_t GetITSEfficiencyCorrectionProton(const float pT){return 1.0/(fEfficiencyProtonITS->GetBinContent(fEfficiencyProtonITS->FindBin(pT)));}
     Float_t GetITSEfficiencyCorrectionHadron(const float pT){return 1.0/(fEfficiencyHadronITS->GetBinContent(fEfficiencyHadronITS->FindBin(pT)));}
     //...and these guys are too
-    Float_t GetBackgroundCorrectionTPC(const float pT){return 1.0/(fBackgroundTPC->GetBinContent(fBackgroundTPC->FindBin(pT)));}
-    Float_t GetBackgroundCorrectionITS(const float pT){return 1.0/(fBackgroundITS->GetBinContent(fBackgroundITS->FindBin(pT)));}
+    Float_t GetBackgroundCorrectionTPC(const float pT){return 1.0/(1.0-fBackgroundTPC->GetBinContent(fBackgroundTPC->FindBin(pT)));}
+    Float_t GetBackgroundCorrectionITS(const float pT){return 1.0/(1.0-fBackgroundITS->GetBinContent(fBackgroundITS->FindBin(pT)));}
 
 
     void SetEtaCut(const Float_t val){fEtaCut=val;}
@@ -100,6 +100,9 @@ public:
     void SetBackgroundCorrectionTPC(const TH1D *histo){fBackgroundTPC=(TH1D*) histo;}
     void SetBackgroundCorrectionITS(const TH1D *histo){fBackgroundITS=(TH1D*) histo;}
 
+    //Returns the factor one needs to multiply by to get the corrected et for all constant (not pt dependent) factors
+    Float_t GetConstantCorrections(Bool_t totEt, Float_t ptcut, TString type);
+
 
     AliAnalysisHadEtCorrections(const AliAnalysisHadEtCorrections *g) ; // cpy ctor
     //AliAnalysisHadEtCorrections & operator = (const AliAnalysisHadEtCorrections & g) ;//cpy assignment
index bdeab40..4747a8b 100644 (file)
@@ -915,14 +915,5 @@ void AliAnalysisHadEtMonteCarlo::CreateHistograms(){
   CreateIntHisto1D("NEvents","Number of events","number of events","Number of events",1,0,1);
 
   //CreateHisto1D("MisidentifiedPIDs","PIDs for particles misidentified that are not a #pi,K,p","PID","number of entries",3000,0.5,3000.5);
-
-
-
-//     list->Add(fHistEt);
-//     TString histname = "fHistEt" + fHistogramNameSuffix;
-
-//     fHistEt = new TH1F(histname.Data(), "Total E_{T} Distribution", 1000, 0.00, 99);
-//     fHistEt->GetXaxis()->SetTitle("E_{T} (GeV/c^{2})");
-//     fHistEt->GetYaxis()->SetTitle("dN/dE_{T} (c^{2}/GeV)");
 }
 
index b893202..4bdf9fb 100644 (file)
@@ -28,6 +28,9 @@ public:
     virtual void Init();
     
  private:
+    //Declare it private to avoid compilation warning
+    AliAnalysisHadEtMonteCarlo & operator = (const AliAnalysisHadEtMonteCarlo & g) ;//cpy assignment
+    AliAnalysisHadEtMonteCarlo(const AliAnalysisHadEtMonteCarlo & g) ; // cpy ctor
 
     ClassDef(AliAnalysisHadEtMonteCarlo, 1);
 };
index 247ae42..7532046 100644 (file)
@@ -7,14 +7,23 @@
 //Created by Christine Nattrass, Rebecca Scott, Irakli Martashvili
 //University of Tennessee at Knoxville
 //_________________________________________________________________________
+
+#include <TROOT.h>
+#include <TSystem.h>
+#include <TInterpreter.h>
 #include "AliAnalysisHadEtReconstructed.h"
 #include "AliAnalysisEtCuts.h"
 #include "AliESDtrack.h"
 #include "AliESDCaloCluster.h"
 #include "AliVEvent.h"
 #include "AliESDEvent.h"
+#include "AliESDtrackCuts.h"
+#include "AliESDpid.h"
 #include "AliVParticle.h"
 #include <iostream>
+#include "AliAnalysisHadEtCorrections.h"
+#include "TFile.h"
+#include "TString.h"
 
 using namespace std;
 
@@ -23,6 +32,26 @@ ClassImp(AliAnalysisHadEtReconstructed);
 
 AliAnalysisHadEtReconstructed::AliAnalysisHadEtReconstructed() :
         AliAnalysisHadEt()
+       ,corrections(0)
+       ,fConfigFile("ConfigHadEtAnalysis.C")
+       ,fCorrTotEtFullAcceptanceTPC(0)
+       ,fCorrTotEtFullAcceptanceITS(0)
+       ,fCorrHadEtFullAcceptanceTPC(0)
+       ,fCorrHadEtFullAcceptanceITS(0)
+       ,fCorrTotEtEMCALAcceptanceTPC(0)
+       ,fCorrTotEtEMCALAcceptanceITS(0)
+       ,fCorrHadEtEMCALAcceptanceTPC(0)
+       ,fCorrHadEtEMCALAcceptanceITS(0)
+       ,fCorrTotEtPHOSAcceptanceTPC(0)
+       ,fCorrTotEtPHOSAcceptanceITS(0)
+       ,fCorrHadEtPHOSAcceptanceTPC(0)
+       ,fCorrHadEtPHOSAcceptanceITS(0)
+       ,fRawEtFullAcceptanceTPC(0)
+       ,fRawEtFullAcceptanceITS(0)
+       ,fRawEtEMCALAcceptanceTPC(0)
+       ,fRawEtEMCALAcceptanceITS(0)
+       ,fRawEtPHOSAcceptanceTPC(0)
+       ,fRawEtPHOSAcceptanceITS(0)
 {
 
 }
@@ -34,49 +63,159 @@ AliAnalysisHadEtReconstructed::~AliAnalysisHadEtReconstructed()
 Int_t AliAnalysisHadEtReconstructed::AnalyseEvent(AliVEvent* ev)
 { // analyse ESD event
     ResetEventValues();
-    AliESDEvent *event = dynamic_cast<AliESDEvent*>(ev);
-
-    for (Int_t iTrack = 0; iTrack < event->GetNumberOfTracks(); iTrack++)
-    {
-        AliVParticle *track = event->GetTrack(iTrack);
-        if (!track)
-        {
-            Printf("ERROR: Could not get track %d", iTrack);
-            continue;
-        }
-
-        fMultiplicity++;
-
-        const Double_t *pidWeights = track->PID();
-        if (pidWeights)
-        {
-            Int_t maxpid = -1;
-            Float_t maxpidweight = 0;
-            for (Int_t p =0; p < AliPID::kSPECIES; p++)
-            {
-                if (pidWeights[p] > maxpidweight)
-                {
-                    maxpidweight = pidWeights[p];
-                    maxpid = p;
-                }
-            }
-            if (maxpid == AliPID::kProton)
-            {
-                //     massPart = -0.938*track->Charge();
-            }
-
-        }
+    fRawEtFullAcceptanceTPC=0.0;
+    fRawEtFullAcceptanceITS=0.0;
+    fRawEtEMCALAcceptanceTPC=0.0;
+    fRawEtEMCALAcceptanceITS=0.0;
+    fRawEtPHOSAcceptanceTPC=0.0;
+    fRawEtPHOSAcceptanceITS=0.0;
+
+    AliESDEvent *realEvent = dynamic_cast<AliESDEvent*>(ev);
+    //for PID
+    AliESDpid *pID = new AliESDpid();
+    pID->MakePID(realEvent);
+
+    TString *strTPC = new TString("TPC");
+    TString *strITS = new TString("ITS");
+    TString *strTPCITS = new TString("TPCITS");
+    for(Int_t cutset=0;cutset<2;cutset++){
+      TString *cutName;
+      TObjArray* list;
+      switch(cutset){
+      case 0:
+       cutName = strTPC;
+       list = fEsdtrackCutsTPC->GetAcceptedTracks(realEvent);
+       break;
+      case 1:
+       cutName = strITS;
+       list = fEsdtrackCutsITS->GetAcceptedTracks(realEvent);
+       break;
+      case 2:
+       cutName = strTPCITS;
+       list = fEsdtrackCutsITSTPC->GetAcceptedTracks(realEvent);
+       break;
+      default:
+       cerr<<"Error:  cannot fill histograms!"<<endl;
+       return -1;
+      }
+      Int_t nGoodTracks = list->GetEntries();
+      for (Int_t iTrack = 0; iTrack < nGoodTracks; iTrack++)
+       {
 
-    }
 
-    fTotNeutralEtAcc = fTotNeutralEt;
-    fTotEt = fTotChargedEt + fTotNeutralEt;
-    fTotEtAcc = fTotChargedEtAcc + fTotNeutralEtAcc;
+         AliESDtrack *track = dynamic_cast<AliESDtrack*> (list->At(iTrack));
+         if (!track)
+           {
+             Printf("ERROR: Could not get track %d", iTrack);
+             continue;
+           }
+         else{
+           Float_t nSigmaPion,nSigmaProton,nSigmaKaon,nSigmaElectron;
+           
+           if(cutset!=1){
+             nSigmaPion = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kPion));
+             nSigmaProton = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kProton));
+             nSigmaKaon = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kKaon));
+             nSigmaElectron = TMath::Abs(pID->NumberOfSigmasTPC(track,AliPID::kElectron));
+           }
+           else{
+             nSigmaPion = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kPion));
+             nSigmaProton = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kProton));
+             nSigmaKaon = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kKaon));
+             nSigmaElectron = TMath::Abs(pID->NumberOfSigmasITS(track,AliPID::kElectron));
+           }
+           bool isPion = (nSigmaPion<3.0 && nSigmaProton>2.0 && nSigmaKaon>2.0);
+           bool isElectron = (nSigmaElectron<2.0 && nSigmaPion>4.0 && nSigmaProton>3.0 && nSigmaKaon>3.0);
+           bool isKaon = (nSigmaPion>3.0 && nSigmaProton>2.0 && nSigmaKaon<2.0);
+           bool isProton = (nSigmaPion>3.0 && nSigmaProton<2.0 && nSigmaKaon>2.0);
 
-    std::cout << fTotChargedEtAcc << std::endl;
-    // Fill the histograms...
-    FillHistograms();
+           //bool IsElectron = false;
+           bool unidentified = (!isProton && !isKaon && !isElectron);
+           Float_t dEdx = track->GetTPCsignal();
+           if(cutset==1) dEdx = track->GetITSsignal();
+           FillHisto2D(Form("dEdxDataAll%s",cutName->Data()),track->P(),dEdx,1.0);
 
+           Float_t corrBkgd=0.0;
+           Float_t corrNotID=0.0;
+           Float_t corrNoID = corrections->GetNotIDCorrectionNoPID(track->Pt());
+           Float_t corrEff = 0.0;
+           Float_t corrEffNoID = 0.0;
+           if(cutset==0){//TPC
+             corrBkgd = corrections->GetBackgroundCorrectionTPC(track->Pt());
+             //corrEffNoID = corrections->GetTPCEfficiencyCorrectionHadron(track->Pt());
+             corrNotID = corrections->GetNotIDCorrectionTPC(track->Pt());
+           }
+           if(cutset==1){//ITS
+             corrBkgd = corrections->GetBackgroundCorrectionITS(track->Pt());
+             //corrEffNoID = corrections->GetITSEfficiencyCorrectionHadron(track->Pt());
+             corrNotID = corrections->GetNotIDCorrectionITS(track->Pt());
+           }
+           Float_t et = 0.0;
+           Float_t etpartialcorrected = 0.0;
+           Float_t etpartialcorrectedNoID = corrNoID*corrBkgd*corrEffNoID*Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
+           FillHisto2D(Form("EtDataRaw%sNoID",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrectedNoID);
+
+           if(isPion){
+             FillHisto2D(Form("dEdxDataPion%s",cutName->Data()),track->P(),dEdx,1.0);
+             et = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
+             //if(cutset==0){corrEff = corrections->GetTPCEfficiencyCorrectionPion(track->Pt());}
+             //else{corrEff = corrections->GetITSEfficiencyCorrectionPion(track->Pt());}
+             etpartialcorrected = et*corrBkgd*corrEff;
+             
+             if(track->Charge()>0.0){
+               FillHisto2D(Form("EtDataRaw%sPiPlus",cutName->Data()),track->Pt(),track->Eta(),et);
+               FillHisto2D(Form("EtDataCorrected%sPiPlus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+             }
+             else{
+               FillHisto2D(Form("EtDataRaw%sPiMinus",cutName->Data()),track->Pt(),track->Eta(),et);
+               FillHisto2D(Form("EtDataCorrected%sPiMinus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+             }
+           }
+           if(isKaon){
+             FillHisto2D(Form("dEdxDataKaon%s",cutName->Data()),track->P(),dEdx,1.0);
+             et = Et(track->P(),track->Theta(),fKPlusCode,track->Charge());
+             //if(cutset==0){corrEff = corrections->GetTPCEfficiencyCorrectionKaon(track->Pt());}
+             //else{corrEff = corrections->GetITSEfficiencyCorrectionKaon(track->Pt());}
+             etpartialcorrected = et*corrBkgd*corrEff;
+             
+             if(track->Charge()>0.0){
+               FillHisto2D(Form("EtDataRaw%sKPlus",cutName->Data()),track->Pt(),track->Eta(),et);
+               FillHisto2D(Form("EtDataCorrected%sKPlus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+             }
+             else{
+               FillHisto2D(Form("EtDataRaw%sKMinus",cutName->Data()),track->Pt(),track->Eta(),et);
+               FillHisto2D(Form("EtDataCorrected%sKMinus",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+             }
+           }
+           if(isProton){
+             FillHisto2D(Form("dEdxDataProton%s",cutName->Data()),track->P(),dEdx,1.0);
+             et = Et(track->P(),track->Theta(),fProtonCode,track->Charge());
+             //if(cutset==0){corrEff = corrections->GetTPCEfficiencyCorrectionProton(track->Pt());}
+             //else{corrEff = corrections->GetITSEfficiencyCorrectionProton(track->Pt());}
+             etpartialcorrected = et*corrBkgd*corrEff;
+             
+             if(track->Charge()>0.0){
+               FillHisto2D(Form("EtDataRaw%sProton",cutName->Data()),track->Pt(),track->Eta(),et);
+               FillHisto2D(Form("EtDataCorrected%sProton",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+             }
+             else{
+               FillHisto2D(Form("EtDataRaw%sAntiProton",cutName->Data()),track->Pt(),track->Eta(),et);
+               FillHisto2D(Form("EtDataCorrected%sAntiProton",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+             }
+           }
+           if(isElectron){
+             FillHisto2D(Form("dEdxDataProton%s",cutName->Data()),track->P(),dEdx,1.0);
+             //et = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
+           }
+           if(unidentified){
+             FillHisto2D(Form("dEdxDataUnidentified%s",cutName->Data()),track->P(),dEdx,1.0);
+             et = Et(track->P(),track->Theta(),fPiPlusCode,track->Charge());
+             etpartialcorrected = et*corrBkgd*corrEffNoID*corrNotID;
+             FillHisto2D(Form("EtDataCorrected%sUnidentified",cutName->Data()),track->Pt(),track->Eta(),etpartialcorrected);
+           }
+         }
+       }
+    }
     return 1;
 }
 
@@ -99,6 +238,90 @@ bool AliAnalysisHadEtReconstructed::CheckGoodVertex(AliVParticle* track)
 void AliAnalysisHadEtReconstructed::Init()
 { // Init
     AliAnalysisHadEt::Init();
+
+  if (fConfigFile.Length()) {
+    gROOT->LoadMacro(fConfigFile);
+    corrections = (AliAnalysisHadEtCorrections *) gInterpreter->ProcessLine("ConfigHadEtAnalysis()");
+    fCorrTotEtFullAcceptanceTPC = corrections->GetConstantCorrections(kTRUE,0.15,"Full");
+    fCorrTotEtFullAcceptanceITS = corrections->GetConstantCorrections(kTRUE,0.1,"Full");
+    fCorrHadEtFullAcceptanceTPC = corrections->GetConstantCorrections(kFALSE,0.15,"Full");
+    fCorrHadEtFullAcceptanceITS = corrections->GetConstantCorrections(kFALSE,0.1,"Full");
+    fCorrTotEtEMCALAcceptanceTPC = corrections->GetConstantCorrections(kTRUE,0.15,"EMCAL");
+    fCorrTotEtEMCALAcceptanceITS = corrections->GetConstantCorrections(kTRUE,0.1,"EMCAL");
+    fCorrHadEtEMCALAcceptanceTPC = corrections->GetConstantCorrections(kFALSE,0.15,"EMCAL");
+    fCorrHadEtEMCALAcceptanceITS = corrections->GetConstantCorrections(kFALSE,0.1,"EMCAL");
+    fCorrTotEtPHOSAcceptanceTPC = corrections->GetConstantCorrections(kTRUE,0.15,"PHOS");
+    fCorrTotEtPHOSAcceptanceITS = corrections->GetConstantCorrections(kTRUE,0.1,"PHOS");
+    fCorrHadEtPHOSAcceptanceTPC = corrections->GetConstantCorrections(kFALSE,0.15,"PHOS");
+    fCorrHadEtPHOSAcceptanceITS = corrections->GetConstantCorrections(kFALSE,0.1,"PHOS");
+
+  }
 }
 
 
+void AliAnalysisHadEtReconstructed::CreateHistograms(){
+
+  TString *strTPC = new TString("TPC");
+  TString *strITS = new TString("ITS");
+  TString *strTPCITS = new TString("TPCITS");
+  for(Int_t i=0;i<2;i++){
+    TString *cutName;
+    Float_t maxPtdEdx = 10;
+    Float_t mindEdx = 35;
+    Float_t maxdEdx = 150.0;
+    switch(i){
+    case 0:
+      cutName = strTPC;
+      break;
+    case 1:
+      cutName = strITS;
+      maxPtdEdx = 5;
+      maxdEdx = 500.0;
+      break;
+    case 2:
+      cutName = strTPCITS;
+      break;
+    default:
+      cerr<<"Error:  cannot make histograms!"<<endl;
+      return;
+    }
+
+    CreateEtaPtHisto2D(Form("EtDataRaw%sPiPlus",cutName->Data()),"Raw reconstructed E_{T} from identified #pi^{+}");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sPiMinus",cutName->Data()),"Raw reconstructed E_{T} from identified #pi^{-}");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sKPlus",cutName->Data()),"Raw reconstructed E_{T} from identified K^{+}");
+//     CreateEtaPtHisto2D(Form("EtDataRaw%sEMinus",cutName->Data()),"Raw reconstructed E_{T} from identified e^{-}");
+//     CreateEtaPtHisto2D(Form("EtDataRaw%sEPlus",cutName->Data()),"Raw reconstructed E_{T} from identified e^{+}");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sKMinus",cutName->Data()),"Raw reconstructed E_{T} from identified K^{-}");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sProton",cutName->Data()),"Raw reconstructed E_{T} from identified p");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sAntiProton",cutName->Data()),"Raw reconstructed E_{T} from identified #bar{p}");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sUnidentified",cutName->Data()),"Raw reconstructed E_{T} from unidentified particles using real mass");
+    CreateEtaPtHisto2D(Form("EtDataRaw%sNoID",cutName->Data()),"Raw reconstructed E_{T} from unidentified particles using real mass");
+
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sPiPlus",cutName->Data()),"Corrected reconstructed E_{T} from identified #pi^{+}");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sPiMinus",cutName->Data()),"Corrected reconstructed E_{T} from identified #pi^{-}");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sKPlus",cutName->Data()),"Corrected reconstructed E_{T} from identified K^{+}");
+//     CreateEtaPtHisto2D(Form("EtDataCorrected%sEMinus",cutName->Data()),"Corrected reconstructed E_{T} from identified e^{-}");
+//     CreateEtaPtHisto2D(Form("EtDataCorrected%sEPlus",cutName->Data()),"Corrected reconstructed E_{T} from identified e^{+}");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sKMinus",cutName->Data()),"Corrected reconstructed E_{T} from identified K^{-}");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sProton",cutName->Data()),"Corrected reconstructed E_{T} from identified p");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sAntiProton",cutName->Data()),"Corrected reconstructed E_{T} from identified #bar{p}");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sUnidentified",cutName->Data()),"Corrected reconstructed E_{T} from unidentified particles using real mass");
+    CreateEtaPtHisto2D(Form("EtDataCorrected%sNoID",cutName->Data()),"Corrected reconstructed E_{T} from unidentified particles using real mass");
+
+
+    CreateEtaPtHisto2D(Form("EtNData%sPiPlus",cutName->Data()),"Number of reconstructed #pi^{+}");
+    CreateEtaPtHisto2D(Form("EtNData%sPiMinus",cutName->Data()),"Number of reconstructed #pi^{-}");
+    CreateEtaPtHisto2D(Form("EtNData%sKPlus",cutName->Data()),"Number of reconstructed K^{+}");
+    CreateEtaPtHisto2D(Form("EtNData%sKMinus",cutName->Data()),"Number of reconstructed K^{-}");
+    CreateEtaPtHisto2D(Form("EtNData%sProton",cutName->Data()),"Number of reconstructed p");
+    CreateEtaPtHisto2D(Form("EtNData%sAntiProton",cutName->Data()),"Number of reconstructed #bar{p}");
+    CreateEtaPtHisto2D(Form("EtNData%sUnidentified",cutName->Data()),"Number of Reconstructed unidentified particles");
+
+    CreateHisto2D(Form("dEdxDataAll%s",cutName->Data()),"dE/dx for all particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxDataPion%s",cutName->Data()),"dE/dx for #pi^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxDataKaon%s",cutName->Data()),"dE/dx for K^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxDataProton%s",cutName->Data()),"dE/dx for p(#bar{p})","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxDataElectron%s",cutName->Data()),"dE/dx for e^{#pm}","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+    CreateHisto2D(Form("dEdxDataUnidentified%s",cutName->Data()),"dE/dx for unidentified particles","momentum (GeV/c)","dE/dx",400,0.0,maxPtdEdx,200,mindEdx,maxdEdx);
+  }
+}
index 2d6e759..a9dc380 100644 (file)
@@ -12,6 +12,8 @@
 #include "AliAnalysisHadEt.h"
 
 class AliVParticle;
+class AliAnalysisHadEtCorrections;
+class TString;
 
 class AliAnalysisHadEtReconstructed : public AliAnalysisHadEt
 {
@@ -21,15 +23,46 @@ public:
     AliAnalysisHadEtReconstructed();
     virtual ~AliAnalysisHadEtReconstructed();
    
+    virtual void SetConfigFile(const char *c) {fConfigFile = c;}
     virtual Int_t AnalyseEvent(AliVEvent* event);
 
+    void CreateHistograms();
     virtual void Init();
     
 protected:
 
     bool CheckGoodVertex(AliVParticle *track);
+    AliAnalysisHadEtCorrections *corrections;
+
+    TString       fConfigFile;        // the name of the ConfigFile
     //virtual bool TrackHitsCalorimeter(AliVParticle *track, Double_t magField) = 0;
 
+    Float_t fCorrTotEtFullAcceptanceTPC;
+    Float_t fCorrTotEtFullAcceptanceITS;
+    Float_t fCorrHadEtFullAcceptanceTPC;
+    Float_t fCorrHadEtFullAcceptanceITS;
+    Float_t fCorrTotEtEMCALAcceptanceTPC;
+    Float_t fCorrTotEtEMCALAcceptanceITS;
+    Float_t fCorrHadEtEMCALAcceptanceTPC;
+    Float_t fCorrHadEtEMCALAcceptanceITS;
+    Float_t fCorrTotEtPHOSAcceptanceTPC;
+    Float_t fCorrTotEtPHOSAcceptanceITS;
+    Float_t fCorrHadEtPHOSAcceptanceTPC;
+    Float_t fCorrHadEtPHOSAcceptanceITS;
+    Float_t fRawEtFullAcceptanceTPC;
+    Float_t fRawEtFullAcceptanceITS;
+    Float_t fRawEtEMCALAcceptanceTPC;
+    Float_t fRawEtEMCALAcceptanceITS;
+    Float_t fRawEtPHOSAcceptanceTPC;
+    Float_t fRawEtPHOSAcceptanceITS;
+
+
+ private:
+    //Declare it private to avoid compilation warning
+    AliAnalysisHadEtReconstructed & operator = (const AliAnalysisHadEtReconstructed & g) ;//cpy assignment
+    AliAnalysisHadEtReconstructed(const AliAnalysisHadEtReconstructed & g) ; // cpy ctor
+
+
     ClassDef(AliAnalysisHadEtReconstructed, 1);
 };
 
diff --git a/PWG4/totEt/macros/ConfigHadEtAnalysis.C b/PWG4/totEt/macros/ConfigHadEtAnalysis.C
new file mode 100644 (file)
index 0000000..1c2799e
--- /dev/null
@@ -0,0 +1,10 @@
+
+AliAnalysisHadEtCorrections * ConfigHadEtAnalysis(){
+  //cout<<"Hello I am configuring you"<<endl;
+
+  TFile *infile = new TFile("corrections.root");
+  corrections = (AliAnalysisHadEtCorrections *)infile->Get("hadCorrectionEMCAL");
+  cout<<"Setting the AliAnalysisHadEtCorrections to "<<corrections->GetName()<<endl;
+  cout<<"eta cut is "<<corrections->GetEtaCut()<<endl;
+  return corrections;
+}
index 5ac12b2..ebd8f26 100644 (file)
 // using ACLiC on the worker nodes.
    //plugin->SetAnalysisSource("AliAnalysisTaskHadEt.cxx");
    //plugin->SetAnalysisSource("AliAnalysisEt.cxx AliAnalysisEtMonteCarlo.cxx AliAnalysisEtMonteCarloPhos.cxx AliAnalysisEtReconstructed.cxx AliAnalysisEtReconstructedPhos.cxx AliAnalysisHadEt.cxx AliAnalysisHadEtMonteCarlo.cxx AliAnalysisHadEtReconstructed.cxx AliAnalysisTaskHadEt.cxx AliAnalysisTaskTotEt.cxx");
-   plugin->SetAnalysisSource("AliAnalysisEtCuts.cxx AliAnalysisHadEt.cxx AliAnalysisHadEtMonteCarlo.cxx AliAnalysisHadEtReconstructed.cxx AliAnalysisTaskHadEt.cxx");
+   plugin->SetAnalysisSource("AliAnalysisEtCuts.cxx AliAnalysisHadEtCorrections.cxx AliAnalysisHadEt.cxx AliAnalysisHadEtMonteCarlo.cxx AliAnalysisHadEtReconstructed.cxx AliAnalysisTaskHadEt.cxx");
 // Declare all libraries (other than the default ones for the framework. These will be
 // loaded by the generated analysis macro. Add all extra files (task .cxx/.h) here.
-   plugin->SetAdditionalLibs("AliAnalysisEtCuts.h AliAnalysisEtCuts.cxx AliAnalysisHadEt.cxx AliAnalysisHadEtMonteCarlo.cxx AliAnalysisHadEtReconstructed.cxx AliAnalysisTaskHadEt.cxx AliAnalysisHadEt.h AliAnalysisHadEtMonteCarlo.h AliAnalysisHadEtReconstructed.h AliAnalysisTaskHadEt.h");
+   plugin->SetAdditionalLibs("AliAnalysisEtCuts.h AliAnalysisEtCuts.cxx AliAnalysisHadEtCorrections.h AliAnalysisHadEtCorrections.cxx AliAnalysisHadEt.cxx AliAnalysisHadEtMonteCarlo.cxx AliAnalysisHadEtReconstructed.cxx AliAnalysisTaskHadEt.cxx AliAnalysisHadEt.h AliAnalysisHadEtMonteCarlo.h AliAnalysisHadEtReconstructed.h AliAnalysisTaskHadEt.h corrections.root ConfigHadEtAnalysis.C");
 // No need for output file names. Procedure is automatic. <-- not true
    plugin->SetDefaultOutputs(kFALSE);
    plugin->SetOutputFiles("Et.ESD.new.sim.root");
index 8b5bfa9..f68385d 100644 (file)
@@ -21,6 +21,7 @@ void runHadEt(bool submit = false) {
 
     gSystem->AddIncludePath("-I$ALICE_ROOT/include");
    gROOT->ProcessLine(".L AliAnalysisEtCuts.cxx+g");
+   gROOT->ProcessLine(".L AliAnalysisHadEtCorrections.cxx+g");
    gROOT->ProcessLine(".L AliAnalysisHadEt.cxx+g");
    gROOT->ProcessLine(".L AliAnalysisHadEtMonteCarlo.cxx+g");
    gROOT->ProcessLine(".L AliAnalysisHadEtReconstructed.cxx+g");