First version of the mean opt fluctuations code (Stefan Heckel)
authorpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 11 May 2011 08:06:26 +0000 (08:06 +0000)
committerpchrist <pchrist@f7af4fe6-9843-0410-8265-dc069ae4e863>
Wed, 11 May 2011 08:06:26 +0000 (08:06 +0000)
PWG2/CMakelibPWG2ebye.pkg
PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.cxx [new file with mode: 0755]
PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.h [new file with mode: 0755]
PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.cxx [new file with mode: 0644]
PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.h [new file with mode: 0644]
PWG2/EBYE/macros/AddTaskPtFluc.C [new file with mode: 0644]
PWG2/EBYE/macros/AddTaskPtFlucPbPb.C [new file with mode: 0644]
PWG2/PWG2ebyeLinkDef.h

index 749b56b..9602738 100644 (file)
@@ -25,7 +25,7 @@
 # SHLIBS - Shared Libraries and objects for linking (Executables only)           #
 #--------------------------------------------------------------------------------#
 
-set ( SRCS   EBYE/AliAnalysisTaskBF.cxx EBYE/AliBalance.cxx EBYE/LRC/AliAnalysisTaskLRC.cxx EBYE/LRC/AliLRCAnalysis.cxx EBYE/LRC/AliLRCFit.cxx EBYE/LRC/AliLRCNN.cxx EBYE/LRC/AliLRCPtN.cxx EBYE/LRC/AliLRCPtPt.cxx EBYE/LRC/AliLRCProcess.cxx EBYE/Fluctuations/AliEbyEFluctuationAnalysisTask.cxx)
+set ( SRCS   EBYE/AliAnalysisTaskBF.cxx EBYE/AliBalance.cxx EBYE/LRC/AliAnalysisTaskLRC.cxx EBYE/LRC/AliLRCAnalysis.cxx EBYE/LRC/AliLRCFit.cxx EBYE/LRC/AliLRCNN.cxx EBYE/LRC/AliLRCPtN.cxx EBYE/LRC/AliLRCPtPt.cxx EBYE/LRC/AliLRCProcess.cxx EBYE/Fluctuations/AliEbyEFluctuationAnalysisTask.cxx EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.cxx EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.cxx)
 
 string ( REPLACE ".cxx" ".h" HDRS "${SRCS}" )
 
diff --git a/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.cxx b/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.cxx
new file mode 100755 (executable)
index 0000000..8f0f47e
--- /dev/null
@@ -0,0 +1,431 @@
+#include "TChain.h"
+#include "TTree.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TCanvas.h"
+#include "TF1.h"
+#include "TProfile.h"
+#include "TList.h"
+#include "iostream"
+
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisManager.h"
+
+#include "AliESD.h"
+#include "AliESDEvent.h"
+#include "AliESDInputHandler.h"
+#include "AliESDtrackCuts.h"
+
+#include "AliLog.h"
+
+#include "AliAnalysisTaskPtFluc.h"
+
+
+using namespace std;
+
+// Analysis of Pt FLuctuations (pp)
+// Author: Stefan Heckel
+// Version of pp task:   8.0, 19.04.2011
+
+
+ClassImp(AliAnalysisTaskPtFluc)
+
+//________________________________________________________________________
+AliAnalysisTaskPtFluc::AliAnalysisTaskPtFluc(const char *name)
+  :AliAnalysisTaskSE(name),
+  fESD(0),
+  fOutputList(0),
+  fPtSpec(0),
+  fMult(0),
+  fEta(0),
+  fEtaPhiPlus(0),
+  fEtaPhiMinus(0),
+  fVtxZ(0),
+  fVtxZCut(0),
+  fVtxZCont(0),
+  fVtxZTrackCuts(0),
+  fEventMeanPt(0),
+  fEventMeanPtSq(0),
+  fMultEventMeanPt(0),
+  fMultEventMeanPtSq(0),
+  fTwoPartCorrEv(0),
+  fTwoPartCorrEvSq(0),
+  fTwoPartCorrEvSample(0),
+  fTwoPartCorrEvSampleSq(0),
+  fESDTrackCuts(0),
+  fMaxVertexZ(0),
+  fNContributors(0),
+  fMC(0)
+{
+  DefineOutput(1, TList::Class());
+}
+
+//________________________________________________________________________
+AliAnalysisTaskPtFluc::~AliAnalysisTaskPtFluc()
+{
+  if(fOutputList) delete fOutputList;  fOutputList =0;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskPtFluc::UserCreateOutputObjects()
+{
+  // Create histograms
+  // Called once
+
+  OpenFile(1, "RECREATE");
+  fOutputList = new TList();
+  fOutputList->SetOwner();
+
+
+  fPtSpec = new TH1F("fPtSpec","Pt spectrum",50,0,2.5);
+
+  fMult = new TH1F("fMult","Multiplicity distribution",120,-0.5,119.5);
+
+  fEta = new TH1F("fEta","Eta distribution",80,-2,2);
+
+  fEtaPhiPlus = new TH1F("fEtaPhiPlus","Phi distribution for positive eta",62,0,6.2);
+
+  fEtaPhiMinus = new TH1F("fEtaPhiMinus","Phi distribution for negative eta",62,0,6.2);
+
+  fVtxZ = new TH1F("fVtxZ","Vertex Z distribution before cuts",100,-20,20);
+
+  fVtxZCut = new TH1F("fVtxZCut","Vertex Z distribution after vtxZ cut",110,-11,11);
+
+  fVtxZCont = new TH1F("fVtxZCont","Vertex Z distribution after nCont cut",110,-11,11);
+
+  fVtxZTrackCuts = new TH1F("fVtxZTrackCuts","Vertex Z distribution after track cuts",110,-11,11);
+
+
+  fEventMeanPt = new TH1F("fEventMeanPt","Mean-Pt event by event",50,0,2.5);
+
+  fEventMeanPtSq = new TH1F("fEventMeanPtSq","Mean-Pt event by event squared",100,0,5);
+
+  fMultEventMeanPt = new TH1F("fMultEventMeanPt","Mean-Pt event by event vs. Multiplicity",100,-0.5,99.5);
+
+  fMultEventMeanPtSq = new TH1F("fMultEventMeanPtSq","Mean-Pt event by event squared vs. Multiplicity",100,-0.5,99.5);
+
+
+  fTwoPartCorrEv = new TH1F("fTwoPartCorrEv","Two-particle correlator vs. Multiplicity",100,-0.5,99.5);
+
+  fTwoPartCorrEvSq = new TH1F("fTwoPartCorrEvSq","Two-particle correlator squared vs. Multiplicity",100,-0.5,99.5);
+
+  fTwoPartCorrEvSample = new TH1F("fTwoPartCorrEvSample","Two-particle correlator vs. Multiplicity",100,-0.5,99.5);
+
+  fTwoPartCorrEvSampleSq = new TH1F("fTwoPartCorrEvSampleSq","Two-particle correlator squared vs. Multiplicity",100,-0.5,99.5);
+
+
+  // Add histograms to the output list
+  fOutputList->Add(fPtSpec);
+  fOutputList->Add(fMult);
+  fOutputList->Add(fEta);
+  fOutputList->Add(fEtaPhiPlus);
+  fOutputList->Add(fEtaPhiMinus);
+  fOutputList->Add(fVtxZ);
+  fOutputList->Add(fVtxZCut);
+  fOutputList->Add(fVtxZCont);
+  fOutputList->Add(fVtxZTrackCuts);
+  fOutputList->Add(fEventMeanPt);
+  fOutputList->Add(fEventMeanPtSq);
+  fOutputList->Add(fMultEventMeanPt);
+  fOutputList->Add(fMultEventMeanPtSq);
+  fOutputList->Add(fTwoPartCorrEv);
+  fOutputList->Add(fTwoPartCorrEvSq);
+  fOutputList->Add(fTwoPartCorrEvSample);
+  fOutputList->Add(fTwoPartCorrEvSampleSq);
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskPtFluc::UserExec(Option_t *)
+{
+  // Main loop
+  // Called for each event
+
+
+  // --- ESD event handler ---
+
+  AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+
+  if (!esdH) {
+    Printf("ERROR: Could not get ESDInputHandler");
+  }
+  else {
+    fESD = esdH->GetEvent();
+  }
+
+  if (!fESD) {
+    Printf("ERROR: fESD not available");
+    return;
+  }
+
+  if(!fESDTrackCuts) Printf("ERROR: No esd track cut");
+
+  // --- End event handler ---
+
+
+  // --- Vertex cuts ---
+
+  // TPC+ITS vertex
+  const AliESDVertex* vtxESD = fESD->GetPrimaryVertexTracks();
+  Double_t vtxZ = vtxESD->GetZv();
+  Double_t vtxNCont = vtxESD->GetNContributors();
+
+//   // TPConly vertex
+//   const AliESDVertex* vtxESDTPC = fESD->GetPrimaryVertexTPC();
+//   Double_t vtxZ = vtxESDTPC->GetZv();
+//   Double_t vtxNCont = vtxESDTPC->GetNContributors();
+
+  fVtxZ->Fill(vtxZ); // VtxZ before cuts
+
+  // Event cut on the z-Position of the vertex
+  if(vtxZ > fMaxVertexZ || vtxZ < (-1.*fMaxVertexZ)) {
+//     Printf("VertexZ out of range, Zv = %f",vtxZ);
+    return;
+  }
+
+  fVtxZCut->Fill(vtxZ); // VtxZ after cut on vtxZ
+
+  // Event cut on the number of contributors
+  if(vtxNCont < fNContributors) {
+//     Printf("Vertex has no contributors");
+    return;
+  }
+//   else {
+//     Printf("Vertex nContributors = %i",vtxNCont);
+//   }
+
+  fVtxZCont->Fill(vtxZ); // VtxZ after cut on nContributors
+
+  // --- End vertex cuts ---
+
+
+  Int_t nrTracks = 0;          // count for number of tracks which pass the track cuts
+  Int_t nrESDTracks = 0;       // number of tracks in the ESD file
+  Double_t sumPt = 0.;         // sum of track Pts
+  Double_t twoPartCorrEv = 0.; // Two-particle correlator of one event for multiplicity bin analysis
+  Double_t twoPartCorrEvSample = 0.; // Two-part. corr. of one event for whole sample analysis
+
+  Double_t tracks[120] = {0.}; // array of track Pts, needed for the two-particle correlator
+
+  Double_t trackPt, trackEta, trackPhi;
+  Double_t eventMeanPt, eventMeanPtSq, evMptMult;
+  Double_t twoPartCorrPair, twoPartCorrEvSq;
+  Double_t twoPartCorrPairSample, twoPartCorrEvSampleSq;
+
+  Double_t *nbins = 0x0; // Mean pT values for multiplicity bin analysis
+  Double_t evMptSample;  // Mean pT value for whole sample analysis
+
+  Double_t energy; // beam energy
+
+  energy = fESD->GetBeamEnergy();
+
+
+//   Printf("MC: %d,  %f",fMC,energy);
+
+
+// --- Mean pT values ---
+
+// !!!!! Have to be calculated in a first run - for each sample, which should be analysed, separately! !!!!!
+
+
+// -- Mean pT values for the multiplicity bins (starting with N = 0, therefore the first value has to be 0.000) --
+
+// MC 900 GeV 10e13 (Pythia6 - Perugia0) (Pt-Range: 0.15 - 2)
+Double_t nbinsMC900[100] = {0.000, 0.456, 0.475, 0.490, 0.507, 0.523, 0.535, 0.546, 0.554, 0.561, 0.568, 0.572, 0.578, 0.582, 0.586, 0.590, 0.593, 0.596, 0.602, 0.601, 0.606, 0.605, 0.610, 0.618, 0.610, 0.620, 0.633, 0.615, 0.609, 0.620, 0.632, 0.634, 0.661, 0.647, 0.560, 0.671, 0.660, 0.678, 0.000, 0.612, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+// MC 7 TeV 10f6a (Pythia6 - Perugia0) (Pt-Range: 0.15 - 2)
+Double_t nbinsMC7[100] = {0.000, 0.453, 0.475, 0.498, 0.523, 0.546, 0.563, 0.577, 0.587, 0.595, 0.603, 0.610, 0.616, 0.622, 0.627, 0.632, 0.636, 0.640, 0.644, 0.648, 0.651, 0.655, 0.658, 0.662, 0.664, 0.667, 0.669, 0.671, 0.674, 0.677, 0.679, 0.682, 0.683, 0.686, 0.688, 0.690, 0.692, 0.695, 0.696, 0.696, 0.699, 0.700, 0.701, 0.705, 0.709, 0.705, 0.701, 0.705, 0.714, 0.718, 0.727, 0.702, 0.725, 0.738, 0.724, 0.728, 0.710, 0.729, 0.711, 0.748, 0.702, 0.708, 0.732, 0.000, 0.683, 0.000, 0.000, 0.675, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+
+// // MC 900 GeV 10e12 (Phojet) (Pt-Range: 0.15 - 2)
+// Double_t nbinsMC900[100] = {0.000, 0.478, 0.483, 0.488, 0.494, 0.498, 0.502, 0.505, 0.507, 0.509, 0.513, 0.516, 0.520, 0.522, 0.526, 0.528, 0.531, 0.535, 0.539, 0.540, 0.544, 0.540, 0.548, 0.548, 0.549, 0.556, 0.562, 0.549, 0.557, 0.550, 0.561, 0.559, 0.554, 0.547, 0.539, 0.564, 0.546, 0.603, 0.532, 0.583, 0.656, 0.531, 0.458, 0.547, 0.000, 0.633, 0.000, 0.603, 0.588, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+// // MC 7 TeV 10f6 (Phojet) (Pt-Range: 0.15 - 2)
+// Double_t nbinsMC7[100] = {0.000, 0.488, 0.499, 0.512, 0.526, 0.536, 0.544, 0.549, 0.553, 0.557, 0.560, 0.563, 0.566, 0.569, 0.572, 0.574, 0.577, 0.579, 0.580, 0.582, 0.584, 0.585, 0.586, 0.587, 0.588, 0.589, 0.589, 0.590, 0.591, 0.591, 0.590, 0.591, 0.592, 0.591, 0.589, 0.591, 0.590, 0.588, 0.590, 0.589, 0.590, 0.588, 0.589, 0.589, 0.591, 0.588, 0.582, 0.580, 0.582, 0.586, 0.575, 0.575, 0.560, 0.570, 0.569, 0.604, 0.616, 0.549, 0.614, 0.611, 0.480, 0.628, 0.568, 0.000, 0.677, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+
+// Data 900 GeV 10c.pass3 (Pt-Range: 0.15 - 2)
+Double_t nbinsData900[100] = {0.000, 0.490, 0.491, 0.495, 0.503, 0.512, 0.521, 0.529, 0.536, 0.542, 0.548, 0.551, 0.556, 0.561, 0.564, 0.567, 0.570, 0.573, 0.575, 0.578, 0.580, 0.582, 0.583, 0.590, 0.591, 0.591, 0.590, 0.598, 0.593, 0.608, 0.609, 0.609, 0.587, 0.615, 0.603, 0.613, 0.560, 0.583, 0.614, 0.649, 0.584, 0.584, 0.594, 0.000, 0.622, 0.000, 0.000, 0.703, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+// Data 2.76 TeV 11a.pass1 (Pt-Range: 0.15 - 2)
+Double_t nbinsData276[100] = {0.000, 0.491, 0.500, 0.511, 0.524, 0.536, 0.548, 0.558, 0.567, 0.574, 0.581, 0.586, 0.592, 0.596, 0.601, 0.605, 0.609, 0.613, 0.616, 0.619, 0.622, 0.625, 0.628, 0.629, 0.632, 0.634, 0.637, 0.639, 0.641, 0.645, 0.645, 0.646, 0.649, 0.652, 0.655, 0.658, 0.656, 0.655, 0.659, 0.668, 0.660, 0.672, 0.672, 0.665, 0.669, 0.662, 0.666, 0.669, 0.677, 0.673, 0.595, 0.642, 0.700, 0.623, 0.670, 0.704, 0.694, 0.000, 0.000, 0.667, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+// Data 7 TeV 10d.pass2 + 10e.pass2 (Pt-Range: 0.15 - 2)
+Double_t nbinsData7[100] = {0.000, 0.495, 0.498, 0.505, 0.516, 0.528, 0.539, 0.550, 0.559, 0.566, 0.573, 0.579, 0.584, 0.589, 0.593, 0.597, 0.601, 0.605, 0.608, 0.611, 0.614, 0.617, 0.619, 0.622, 0.624, 0.626, 0.629, 0.631, 0.633, 0.634, 0.636, 0.638, 0.640, 0.642, 0.643, 0.645, 0.646, 0.648, 0.649, 0.651, 0.652, 0.654, 0.654, 0.656, 0.657, 0.658, 0.659, 0.660, 0.662, 0.662, 0.663, 0.665, 0.664, 0.665, 0.668, 0.667, 0.670, 0.666, 0.670, 0.670, 0.677, 0.673, 0.672, 0.679, 0.680, 0.667, 0.682, 0.687, 0.675, 0.658, 0.677, 0.663, 0.686, 0.671, 0.682, 0.683, 0.667, 0.000, 0.000, 0.662, 0.720, 0.683, 0.720, 0.622, 0.000, 0.729, 0.000, 0.603, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000, 0.000};
+
+// -- End mean pT values for the multiplicity bins --
+
+
+// -- Selection of MC/Data and energy; whole sample values --
+
+if (fMC) { // - MC -
+
+  if (energy < 1000.) {
+
+//     Printf(" -- MC, 900 GeV -- ");
+
+    evMptSample = 0.5307; // MC 900 GeV 10e13 (Pythia6 - Perugia0) (Pt-Range: 0.15 - 2)
+//     evMptSample = 0.5044; // MC 900 GeV 10e12 (Phojet) (Pt-Range: 0.15 - 2)
+
+    nbins = nbinsMC900;
+  }
+  else {
+
+//     Printf(" -- MC, 7 TeV -- ");
+
+    evMptSample = 0.5847; // MC 7 TeV 10f6a (Pythia6 - Perugia0) (Pt-Range: 0.15 - 2)
+//     evMptSample = 0.5517; // MC 7 TeV 10f6 (Phojet) (Pt-Range: 0.15 - 2)
+
+    nbins = nbinsMC7;
+  }
+
+} // - End MC -
+else { // - Data -
+
+  if (energy < 1000.) {
+
+//     Printf(" -- Data, 900 GeV -- ");
+
+    evMptSample = 0.5267; // Data 900 GeV 10c.pass3 (Pt-Range: 0.15 - 2)
+    nbins = nbinsData900;
+  }
+  else if (energy > 1000. && energy < 3000.) {
+
+//     Printf(" -- Data, 2.76 TeV -- ");
+
+    evMptSample = 0.5622; // Data 2.76 TeV 11a.pass1 (Pt-Range: 0.15 - 2)
+    nbins = nbinsData276;
+  }
+  else {
+
+//     Printf(" -- Data, 7 TeV -- ");
+
+    evMptSample = 0.5743; // Data 7 TeV 10d.pass2 + 10e.pass2 (Pt-Range: 0.15 - 2)
+    nbins = nbinsData7;
+  }
+
+} // - End data -
+
+// -- End selection of MC/Data and energy; whole sample values --
+
+// --- End mean pT values ---
+
+
+  nrESDTracks = fESD->GetNumberOfTracks();
+//   if( nrESDTracks ) { printf("Found event with %i tracks.\n",nrESDTracks); }
+
+
+  // --- Loop over all tracks of one event ---
+
+  for (Int_t iTracks = 0; iTracks < nrESDTracks; iTracks++) {
+    AliESDtrack* track = fESD->GetTrack(iTracks);
+    if (!track) {
+      Printf("ERROR: Could not receive track %d\n", iTracks);
+      continue;
+    }
+
+      if(!fESDTrackCuts->AcceptTrack(track))continue;
+
+       trackPt = track->Pt();
+       fPtSpec->Fill(trackPt);
+       tracks[nrTracks] = trackPt;
+
+       trackEta = track->Eta();
+       fEta->Fill(trackEta);
+
+       trackPhi = track->Phi();
+
+       if (trackEta > 0.) {
+         fEtaPhiPlus->Fill(trackPhi);
+       }
+       else if (trackEta < 0.) {
+         fEtaPhiMinus->Fill(trackPhi);
+       }
+
+       sumPt = sumPt + trackPt;
+       nrTracks++;
+
+//     printf("Track Pt = %.3f   Track Eta = %.3f   Track Phi = %3f\n",trackPt,trackEta,trackPhi);
+
+  } // --- End track loop ---
+
+
+  // --- Calculation of various values and filling of histograms (for remaining events with N_acc > 0) ---
+
+  if(nrTracks != 0) {
+
+    // VtxZ after all track cuts
+    fVtxZTrackCuts->Fill(vtxZ);
+
+    // Multiplicity distribution
+    fMult->Fill(nrTracks);
+
+    // Calculation of mean Pt and mean Pt Squared
+    eventMeanPt = sumPt / nrTracks;
+    eventMeanPtSq = eventMeanPt * eventMeanPt;
+
+    // Mean-Pt and Mean-Pt squared
+    fEventMeanPt->Fill(eventMeanPt);
+    fEventMeanPtSq->Fill(eventMeanPtSq);
+
+    // Mean-Pt and Mean-Pt squared depending on multiplicity
+    fMultEventMeanPt->Fill(nrTracks,eventMeanPt);
+    fMultEventMeanPtSq->Fill(nrTracks,eventMeanPtSq);
+
+//     printf("nrTracks: %i   sumPt: %.8f   meanPt: %.8f   meanPtSq: %.8f\n",nrTracks,sumPt,eventMeanPt,eventMeanPtSq);
+
+
+    // --- Two-particle correlator ---
+
+    evMptMult = nbins[nrTracks];
+//     printf("nrTracks = %3i   evMptMult = %.3f   evMptSample = %.3f \n\n",nrTracks,evMptMult,evMptSample);
+
+    for (int i=0; i<nrTracks; i++) {
+
+//       printf("--- TrackPt = %.3f ",tracks[i]);
+
+      for (int j=i+1; j<nrTracks; j++) {
+
+       twoPartCorrPair = (tracks[i] - evMptMult) * (tracks[j] - evMptMult); // Multiplicity bins
+       twoPartCorrEv = twoPartCorrEv + twoPartCorrPair;
+
+       twoPartCorrPairSample = (tracks[i] - evMptSample) * (tracks[j] - evMptSample); // Whole sample
+       twoPartCorrEvSample = twoPartCorrEvSample + twoPartCorrPairSample;
+
+//     printf("twoPartCorrPair = %.3f   twoPartCorrPairSample = %.3f \n",twoPartCorrPair,twoPartCorrPairSample);
+      }
+//       printf("twoPartCorrEv = %.3f   twoPartCorrEvSample = %.3f \n",twoPartCorrEv,twoPartCorrEvSamples);
+    }
+//     printf("Total Event: twoPartCorrEv = %.3f   twoPartCorrEvSample = %.3f \n\n",twoPartCorrEv,twoPartCorrEvSample);
+
+    // Multiplicity bins
+    twoPartCorrEvSq = twoPartCorrEv * twoPartCorrEv;
+    fTwoPartCorrEv->Fill(nrTracks,twoPartCorrEv);
+    fTwoPartCorrEvSq->Fill(nrTracks,twoPartCorrEvSq);
+
+    // Whole sample
+    twoPartCorrEvSampleSq = twoPartCorrEvSample * twoPartCorrEvSample;
+    fTwoPartCorrEvSample->Fill(nrTracks,twoPartCorrEvSample);
+    fTwoPartCorrEvSampleSq->Fill(nrTracks,twoPartCorrEvSampleSq);
+
+    // --- End two-particle correlator ---
+
+
+  } // --- End calculation of various values and filling of histograms ---
+
+
+  // Post output data
+  PostData(1, fOutputList);
+}
+
+void AliAnalysisTaskPtFluc::Terminate(Option_t *)
+{
+  // Draw result to the screen
+  // Called once at the end of the query
+
+  fOutputList = dynamic_cast<TList*> (GetOutputData(1));
+  if (!fOutputList) {
+    Printf("ERROR: fOutputList not available\n");
+    return;
+  }
+
+}
diff --git a/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.h b/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFluc.h
new file mode 100755 (executable)
index 0000000..02ad9d4
--- /dev/null
@@ -0,0 +1,69 @@
+#ifndef AliAnalysisTaskPtFluc_cxx\r
+#define AliAnalysisTaskPtFluc_cxx\r
+\r
+// Analysis of Pt FLuctuations (pp)\r
+// Author: Stefan Heckel\r
+// Version of pp task:   8.0, 19.04.2011\r
+\r
+\r
+class TList;\r
+class TH1F;\r
+class TProfile;\r
+\r
+class AliESDEvent;\r
+class AliESDtrack;\r
+class AliESDtrackCuts;\r
+\r
+\r
+#include "AliAnalysisTaskSE.h"\r
+\r
+\r
+class AliAnalysisTaskPtFluc : public AliAnalysisTaskSE {\r
+  public:\r
+    AliAnalysisTaskPtFluc(const char *name = "AliAnalysisTaskPtFluc");\r
+    virtual ~AliAnalysisTaskPtFluc();\r
+\r
+    virtual void   UserCreateOutputObjects();\r
+    virtual void   UserExec(Option_t *option);\r
+    virtual void   Terminate(Option_t *);\r
+\r
+    void SetAliESDtrackCuts(AliESDtrackCuts* esdTrackCuts) {fESDTrackCuts = esdTrackCuts;}\r
+    void SetMaxVertexZ(Float_t vZ) {fMaxVertexZ = vZ;}\r
+    void SetNContributors(Int_t nCont) {fNContributors = nCont;}\r
+    void SetMC(Bool_t bMC) {fMC = bMC;}\r
+\r
+\r
+  private:\r
+    AliESDEvent             *fESD;             // ESD object\r
+    TList*           fOutputList;      // List where all the output files are stored\r
+       TH1F*           fPtSpec;        // Pt spectrum\r
+       TH1F*           fMult;          // Multiplicity distribution\r
+       TH1F*           fEta;           // Eta distribution\r
+       TH1F*           fEtaPhiPlus;    // Phi distribution for positive eta\r
+       TH1F*           fEtaPhiMinus;   // Phi distribution for negative eta\r
+       TH1F*           fVtxZ;          // Vertex Z distribution after physics selection before any further cuts\r
+       TH1F*           fVtxZCut;       // Vertex Z dist. after vertex Z cut\r
+       TH1F*           fVtxZCont;      // Vertex Z dist. after vertex cut on nContributors\r
+       TH1F*           fVtxZTrackCuts; // Vertex Z dist. after all event and track cuts\r
+       TH1F*           fEventMeanPt;           // Event mean pT distribution\r
+       TH1F*           fEventMeanPtSq;         // Event mean pT squared dist.\r
+       TH1F*           fMultEventMeanPt;       // Event mean pT for multiplicity bins\r
+       TH1F*           fMultEventMeanPtSq;     // Event mean pT squared for mult. bins\r
+       TH1F*           fTwoPartCorrEv;         // Two-particle correlator for multiplicity bins\r
+       TH1F*           fTwoPartCorrEvSq;       // Two-part. corr. squared for mult. bins\r
+       TH1F*           fTwoPartCorrEvSample;   // Two-part. corr. for the whole sample\r
+       TH1F*           fTwoPartCorrEvSampleSq; // Two-part. corr. squared for the whole sample\r
+\r
+    AliESDtrackCuts* fESDTrackCuts;    // Esd track cuts\r
+    Float_t          fMaxVertexZ;      // Maximum value for Vertex Z position\r
+    Int_t            fNContributors;   // Minimum contributors to the vertex\r
+    Bool_t           fMC;              // Check for MC\r
+\r
+    AliAnalysisTaskPtFluc(const AliAnalysisTaskPtFluc&); // not implemented\r
+    AliAnalysisTaskPtFluc& operator=(const AliAnalysisTaskPtFluc&); // not implemented\r
+\r
+    ClassDef(AliAnalysisTaskPtFluc, 1);\r
+\r
+};\r
+\r
+#endif\r
diff --git a/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.cxx b/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.cxx
new file mode 100644 (file)
index 0000000..ad6e5ef
--- /dev/null
@@ -0,0 +1,525 @@
+#include "TChain.h"
+#include "TTree.h"
+#include "TH1F.h"
+#include "TH2F.h"
+#include "TCanvas.h"
+#include "TF1.h"
+#include "TProfile.h"
+#include "TList.h"
+#include "iostream"
+
+#include "AliAnalysisTaskSE.h"
+#include "AliAnalysisManager.h"
+
+#include "AliESD.h"
+#include "AliESDEvent.h"
+#include "AliESDInputHandler.h"
+#include "AliESDtrackCuts.h"
+#include "AliCentrality.h"
+
+#include "AliLog.h"
+
+#include "AliAnalysisTaskPtFlucPbPb.h"
+
+
+using namespace std;
+
+// Analysis of Pt FLuctuations (PbPb)
+// Author: Stefan Heckel
+// Version of PbPb task: 5.0, 18.04.2011
+
+
+ClassImp(AliAnalysisTaskPtFlucPbPb)
+
+//________________________________________________________________________
+AliAnalysisTaskPtFlucPbPb::AliAnalysisTaskPtFlucPbPb(const char *name)
+  :AliAnalysisTaskSE(name),
+  fESD(0),
+  fOutputList(0),
+  fPtSpec(0),
+  fMult(0),
+  fMultSum(0),
+  fMultSumPt(0),
+  fMultNrPairs(0),
+  fCent(0),
+  fCentSum(0),
+  fCentSumPt(0),
+  fCentNrPairs(0),
+  fEta(0),
+  fEtaPhiPlus(0),
+  fEtaPhiMinus(0),
+  fVtxZ(0),
+  fVtxZCut(0),
+  fVtxZCont(0),
+  fVtxZTrackCuts(0),
+  fEventMeanPt(0),
+  fEventMeanPtSq(0),
+  fMultEventMeanPt(0),
+  fMultEventMeanPtSq(0),
+  fCentEventMeanPt(0),
+  fCentEventMeanPtSq(0),
+  fTwoPartCorrEv(0),
+  fTwoPartCorrEvSq(0),
+  fTwoPartCorrEvCent(0),
+  fTwoPartCorrEvCentSq(0),
+  fESDTrackCuts(0),
+  fMaxVertexZ(0),
+  fNContributors(0),
+  fUseCentrality(0),
+  fMC(0)
+{
+  DefineOutput(1, TList::Class());
+}
+
+//________________________________________________________________________
+AliAnalysisTaskPtFlucPbPb::~AliAnalysisTaskPtFlucPbPb()
+{
+  if(fOutputList) delete fOutputList;  fOutputList =0;
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskPtFlucPbPb::UserCreateOutputObjects()
+{
+  // Create histograms
+  // Called once
+
+  OpenFile(1, "RECREATE");
+  fOutputList = new TList();
+  fOutputList->SetOwner();
+
+
+  fPtSpec = new TH1F("fPtSpec","Pt spectrum",50,0,2.5);
+
+  fMult = new TH1F("fMult","Multiplicity distribution",30,0,3000);
+
+  fMultSum = new TH1F("fMultSum","Sum of nrtracks of all events in multiplicity bins",30,0,3000);
+
+  fMultSumPt = new TH1F("fMultSumPt","Sum of pTs of all events in multiplicity bins",30,0,3000);
+
+  fMultNrPairs = new TH1F("fMultNrPairs","Sum of number of pairs in multiplicity bins",30,0,3000);
+
+  fCent = new TH1F("fCent","Centrality distribution",21,0,105);
+
+  fCentSum = new TH1F("fCentSum","Sum of nrtracks of all events in centrality bins",21,0,105);
+
+  fCentSumPt = new TH1F("fCentSumPt","Sum of pTs of all events in centrality bins",21,0,105);
+
+  fCentNrPairs = new TH1F("fCentNrPairs","Sum of number of pairs in centrality bins",21,0,105);
+
+
+  fEta = new TH1F("fEta","Eta distribution",80,-2,2);
+
+  fEtaPhiPlus = new TH1F("fEtaPhiPlus","Phi distribution for positive eta",62,0,6.2);
+
+  fEtaPhiMinus = new TH1F("fEtaPhiMinus","Phi distribution for negative eta",62,0,6.2);
+
+  fVtxZ = new TH1F("fVtxZ","Vertex Z distribution before cuts",100,-20,20);
+
+  fVtxZCut = new TH1F("fVtxZCut","Vertex Z distribution after vtxZ cut",110,-11,11);
+
+  fVtxZCont = new TH1F("fVtxZCont","Vertex Z distribution after nCont cut",110,-11,11);
+
+  fVtxZTrackCuts = new TH1F("fVtxZTrackCuts","Vertex Z distribution after track cuts",110,-11,11);
+
+
+  fEventMeanPt = new TH1F("fEventMeanPt","Mean-Pt event by event",50,0,2.5);
+
+  fEventMeanPtSq = new TH1F("fEventMeanPtSq","Mean-Pt event by event squared",100,0,5);
+
+  fMultEventMeanPt = new TH1F("fMultEventMeanPt","Mean-Pt event by event vs. multiplicity",30,0,3000);
+
+  fMultEventMeanPtSq = new TH1F("fMultEventMeanPtSq","Mean-Pt event by event squared vs. multiplicity",30,0,3000);
+
+  fCentEventMeanPt = new TH1F("fCentEventMeanPt","Mean-Pt event by event vs. centrality",21,0,105);
+
+  fCentEventMeanPtSq = new TH1F("fCentEventMeanPtSq","Mean-Pt event by event squared vs. centrality",21,0,105);
+
+
+  fTwoPartCorrEv = new TH1F("fTwoPartCorrEv","Two-particle correlator vs. multiplicity",30,0,3000);
+
+  fTwoPartCorrEvSq = new TH1F("fTwoPartCorrEvSq","Two-particle correlator squared vs. multiplicity",30,0,3000);
+
+  fTwoPartCorrEvCent = new TH1F("fTwoPartCorrEvCent","Two-particle correlator vs. centrality",21,0,105);
+
+  fTwoPartCorrEvCentSq = new TH1F("fTwoPartCorrEvCentSq","Two-particle correlator squared vs. centrality",21,0,105);
+
+
+  // Add histograms to the output list
+  fOutputList->Add(fPtSpec);
+  fOutputList->Add(fMult);
+  fOutputList->Add(fMultSum);
+  fOutputList->Add(fMultSumPt);
+  fOutputList->Add(fMultNrPairs);
+  fOutputList->Add(fCent);
+  fOutputList->Add(fCentSum);
+  fOutputList->Add(fCentSumPt);
+  fOutputList->Add(fCentNrPairs);
+  fOutputList->Add(fEta);
+  fOutputList->Add(fEtaPhiPlus);
+  fOutputList->Add(fEtaPhiMinus);
+  fOutputList->Add(fVtxZ);
+  fOutputList->Add(fVtxZCut);
+  fOutputList->Add(fVtxZCont);
+  fOutputList->Add(fVtxZTrackCuts);
+  fOutputList->Add(fEventMeanPt);
+  fOutputList->Add(fEventMeanPtSq);
+  fOutputList->Add(fMultEventMeanPt);
+  fOutputList->Add(fMultEventMeanPtSq);
+  fOutputList->Add(fCentEventMeanPt);
+  fOutputList->Add(fCentEventMeanPtSq);
+  fOutputList->Add(fTwoPartCorrEv);
+  fOutputList->Add(fTwoPartCorrEvSq);
+  fOutputList->Add(fTwoPartCorrEvCent);
+  fOutputList->Add(fTwoPartCorrEvCentSq);
+}
+
+//________________________________________________________________________
+void AliAnalysisTaskPtFlucPbPb::UserExec(Option_t *)
+{
+  // Main loop
+  // Called for each event
+
+
+  // --- ESD event handler ---
+
+  AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
+
+  if (!esdH) {
+    Printf("ERROR: Could not get ESDInputHandler");
+  }
+  else {
+    fESD = esdH->GetEvent();
+  }
+
+  if (!fESD) {
+    Printf("ERROR: fESD not available");
+    return;
+  }
+
+  if(!fESDTrackCuts) Printf("ERROR: No esd track cut");
+
+  // --- End event handler ---
+
+
+  // --- Vertex cuts ---
+
+//   // TPC+ITS vertex
+//   const AliESDVertex* vtxESD = fESD->GetPrimaryVertexTracks();
+//   Double_t vtxZ = vtxESD->GetZv();
+//   Double_t vtxNCont = vtxESD->GetNContributors();
+
+  // TPConly vertex
+  const AliESDVertex* vtxESDTPC = fESD->GetPrimaryVertexTPC();
+  Double_t vtxZ = vtxESDTPC->GetZv();
+  Double_t vtxNCont = vtxESDTPC->GetNContributors();
+
+  fVtxZ->Fill(vtxZ); // VtxZ before cuts
+
+  // Event cut on the z-Position of the vertex
+  if(vtxZ > fMaxVertexZ || vtxZ < (-1.*fMaxVertexZ)) {
+//     Printf("VertexZ out of range, Zv = %f",vtxZ);
+    return;
+  }
+
+  fVtxZCut->Fill(vtxZ); // VtxZ after cut on vtxZ
+
+  // Event cut on the number of contributors
+  if(vtxNCont < fNContributors) {
+//     Printf("Vertex has no contributors");
+    return;
+  }
+//   else {
+//     Printf("Vertex nContributors = %i",vtxNCont);
+//   }
+
+  fVtxZCont->Fill(vtxZ); // VtxZ after cut on nContributors
+
+  // --- End vertex cuts ---
+
+
+  Int_t nrTracks = 0;          // count for number of tracks which pass the track cuts
+  Int_t nrESDTracks = 0;       // number of tracks in the ESD file
+  Double_t sumPt = 0.;         // sum of track Pts
+  Double_t twoPartCorrEv = 0.; // Two-particle correlator of one event for multiplicity bin analysis
+  Double_t twoPartCorrEvCent = 0.; // Two-particle correlator of one event for centrality bin analysis
+
+  Double_t tracks[3000] = {0.};        // array of track Pts, needed for the two-particle correlator
+
+  Double_t trackPt = 0., trackEta = 0., trackPhi = 0.;
+  Double_t eventMeanPt = 0., eventMeanPtSq = 0., evMptMult = 0.;
+  Double_t twoPartCorrPair = 0., twoPartCorrEvSq = 0.;
+  Double_t twoPartCorrPairCent = 0., twoPartCorrEvCentSq = 0.;
+  Double_t nrPairs = 0.;
+
+  Double_t *nbins = 0x0; // Mean pT values for multiplicity bin analysis
+  Double_t *centbins = 0x0; // Mean pT values for centrality bin analysis
+
+  Double_t centralityVZERO = 0.;
+  Int_t centralityVZEROBin = 0.;
+
+  Int_t centBin = 0.;
+  Double_t evMptCent = 0.;
+
+
+//   Printf("MC: %d",fMC);
+
+
+// --- Mean pT values ---
+
+// !!!!! Have to be calculated in a first run - for each sample, which should be analysed, separately! !!!!!
+
+
+// -- Mean pT values for multiplicity bins (first bin is for nrTracks = 0 and has to be 0) --
+
+  // MC PbPb 2.76 ATeV 11a10a (Hijing) (Pt-Range: 0.15 - 2)
+  Double_t nbinsMC276[32] = {0.000, 0.525, 0.535, 0.539, 0.541, 0.543, 0.544, 0.545, 0.545, 0.546, 0.547, 0.548, 0.548, 0.549, 0.549, 0.550, 0.551, 0.551, 0.552, 0.552, 0.553, 0.554, 0.554, 0.555, 0.556, 0.556, 0.557, 0.558, 0.549, 0.000, 0.000, 0.000};
+
+  // Data PbPb 2.76 ATeV 10h.pass1 (Pt-Range: 0.15 - 2)
+  Double_t nbinsData276[32] = {0.000, 0.580, 0.605, 0.618, 0.626, 0.632, 0.636, 0.640, 0.642, 0.645, 0.647, 0.648, 0.649, 0.650, 0.651, 0.652, 0.653, 0.653, 0.654, 0.654, 0.654, 0.654, 0.654, 0.654, 0.655, 0.656, 0.658, 0.660, 0.661, 0.659, 0.660, 0.656};
+
+// -- End mean pT values for multiplicity bins --
+
+
+// -- Mean pT values for centrality bins --
+
+  // MC PbPb 2.76 ATeV 11a10a (Hijing) (Pt-Range: 0.15 - 2)
+  Double_t centbinsMC276[21] = {0.555, 0.553, 0.551, 0.549, 0.548, 0.547, 0.545, 0.544, 0.543, 0.541, 0.539, 0.537, 0.536, 0.533, 0.530, 0.528, 0.524, 0.518, 0.505, 0.000, 0.534};
+
+  // Data PbPb 2.76 ATeV 10h.pass1 (Pt-Range: 0.15 - 2)
+  Double_t centbinsData276[21] = {0.654, 0.654, 0.652, 0.650, 0.647, 0.644, 0.640, 0.635, 0.629, 0.623, 0.616, 0.609, 0.601, 0.594, 0.586, 0.579, 0.568, 0.551, 0.535, 0.000, 0.621};
+
+// -- End mean pT values for centrality bins --
+
+
+// -- Selection of MC/Data; whole sample values --
+
+if (fMC) { // - MC -
+
+//   Printf(" -- MC, 2.76 ATeV -- ");
+
+  nbins = nbinsMC276;
+  centbins = centbinsMC276;
+
+} // - End MC -
+else { // - Data -
+
+//   Printf(" -- Data, 2.76 ATeV -- ");
+
+  nbins = nbinsData276;
+  centbins = centbinsData276;
+
+} // - End data -
+
+// -- End selection of MC/Data; whole sample values --
+
+// --- End mean pT values ---
+
+
+  nrESDTracks = fESD->GetNumberOfTracks();
+//   if( nrESDTracks ) { printf("Found event with %i tracks.\n",nrESDTracks); }
+
+
+  if(fUseCentrality != 0) {
+
+    AliCentrality *esdCentrality = fESD->GetCentrality();
+
+    //V0
+    if (fUseCentrality == 1){
+
+      centralityVZERO = esdCentrality->GetCentralityPercentile("V0M");
+
+      if      ( centralityVZERO >   0. && centralityVZERO <   5.) centralityVZEROBin =  0;
+      else if ( centralityVZERO >=  5. && centralityVZERO <  10.) centralityVZEROBin =  5;
+      else if ( centralityVZERO >= 10. && centralityVZERO <  15.) centralityVZEROBin = 10;
+      else if ( centralityVZERO >= 15. && centralityVZERO <  20.) centralityVZEROBin = 15;
+      else if ( centralityVZERO >= 20. && centralityVZERO <  25.) centralityVZEROBin = 20;
+      else if ( centralityVZERO >= 25. && centralityVZERO <  30.) centralityVZEROBin = 25;
+      else if ( centralityVZERO >= 30. && centralityVZERO <  35.) centralityVZEROBin = 30;
+      else if ( centralityVZERO >= 35. && centralityVZERO <  40.) centralityVZEROBin = 35;
+      else if ( centralityVZERO >= 40. && centralityVZERO <  45.) centralityVZEROBin = 40;
+      else if ( centralityVZERO >= 45. && centralityVZERO <  50.) centralityVZEROBin = 45;
+      else if ( centralityVZERO >= 50. && centralityVZERO <  55.) centralityVZEROBin = 50;
+      else if ( centralityVZERO >= 55. && centralityVZERO <  60.) centralityVZEROBin = 55;
+      else if ( centralityVZERO >= 60. && centralityVZERO <  65.) centralityVZEROBin = 60;
+      else if ( centralityVZERO >= 65. && centralityVZERO <  70.) centralityVZEROBin = 65;
+      else if ( centralityVZERO >= 70. && centralityVZERO <  75.) centralityVZEROBin = 70;
+      else if ( centralityVZERO >= 75. && centralityVZERO <  80.) centralityVZEROBin = 75;
+      else if ( centralityVZERO >= 80. && centralityVZERO <  85.) centralityVZEROBin = 80;
+      else if ( centralityVZERO >= 85. && centralityVZERO <  90.) centralityVZEROBin = 85;
+      else if ( centralityVZERO >= 90. && centralityVZERO <  95.) centralityVZEROBin = 90;
+      else if ( centralityVZERO >= 95. && centralityVZERO <  99.) centralityVZEROBin = 95;
+      else if ( centralityVZERO >= 99. ) centralityVZEROBin = 100;
+      else if ( centralityVZERO <= 0.  ) centralityVZEROBin = 100;
+
+      centBin = centralityVZEROBin / 5;
+      evMptCent = centbins[centBin];
+    }
+
+//     Printf("\nRaw mult: %i   CentVZERO: %f   CentVZERObin: %i   centBin: %i   evMptCent: %.3f ",nrESDTracks,centralityVZERO,centralityVZEROBin,centBin,evMptCent);
+
+  }
+
+
+  // --- Loop over all tracks of one event ---
+
+  for (Int_t iTracks = 0; iTracks < nrESDTracks; iTracks++) {
+    AliESDtrack* track = fESD->GetTrack(iTracks);
+    if (!track) {
+      Printf("ERROR: Could not receive track %d\n", iTracks);
+      continue;
+    }
+
+      if(!fESDTrackCuts->AcceptTrack(track))continue;
+
+       trackPt = track->Pt();
+       fPtSpec->Fill(trackPt);
+       tracks[nrTracks] = trackPt;
+
+       trackEta = track->Eta();
+       fEta->Fill(trackEta);
+
+       trackPhi = track->Phi();
+
+       if (trackEta > 0.) {
+         fEtaPhiPlus->Fill(trackPhi);
+       }
+       else if (trackEta < 0.) {
+         fEtaPhiMinus->Fill(trackPhi);
+       }
+
+       sumPt = sumPt + trackPt;
+       nrTracks++;
+
+//     printf("Track Pt = %.3f   Track Eta = %.3f   Track Phi = %3f\n",trackPt,trackEta,trackPhi);
+
+  } // --- End track loop ---
+
+
+  // --- Calculation of various values and filling of histograms (for remaining events with N_acc > 0) ---
+
+  if(nrTracks != 0) {
+
+    // VtxZ after all track cuts
+    fVtxZTrackCuts->Fill(vtxZ);
+
+    // Multiplicity distributions
+    fMult->Fill(nrTracks);
+    fMultSum->Fill(nrTracks,nrTracks);
+
+    // Number of pairs in event
+    nrPairs = 0.5 * nrTracks * (nrTracks-1);
+    fMultNrPairs->Fill(nrTracks,nrPairs);
+
+    // Calculation of mean Pt and mean Pt Squared
+    eventMeanPt = sumPt / nrTracks;
+    eventMeanPtSq = eventMeanPt * eventMeanPt;
+
+    // Mean-Pt and Mean-Pt squared
+    fEventMeanPt->Fill(eventMeanPt);
+    fEventMeanPtSq->Fill(eventMeanPtSq);
+
+    // Mean-Pt and Mean-Pt squared depending on multiplicity
+    fMultEventMeanPt->Fill(nrTracks,eventMeanPt);
+    fMultEventMeanPtSq->Fill(nrTracks,eventMeanPtSq);
+    fMultSumPt->Fill(nrTracks,sumPt);
+
+//     printf("nrTracks: %i   sumPt: %.8f   meanPt: %.8f   meanPtSq: %.8f\n",nrTracks,sumPt,eventMeanPt,eventMeanPtSq);
+
+
+    // Centrality V0
+    if (fUseCentrality == 1) {
+
+      // Centrality distributions
+      fCent->Fill(centralityVZEROBin);
+      fCentSum->Fill(centralityVZEROBin,nrTracks);
+
+      // Number of pairs in event
+      fCentNrPairs->Fill(centralityVZEROBin,nrPairs);
+
+      // Mean-Pt and Mean-Pt squared depending on centrality
+      fCentEventMeanPt->Fill(centralityVZEROBin,eventMeanPt);
+      fCentEventMeanPtSq->Fill(centralityVZEROBin,eventMeanPtSq);
+      fCentSumPt->Fill(centralityVZEROBin,sumPt);
+
+//       printf("CentV0bin: %i   NrTracks: %i   NrPairs: %.0f   SumPt: %.1f\n",centralityVZEROBin,nrTracks,nrPairs,sumPt);
+    }
+
+
+    // --- Two-particle correlator ---
+
+    // -- Analysis in multiplicity bins --
+
+    for (int k=1; k<32; k++) {
+      if (nrTracks > 100 * (k-1) && nrTracks <= 100 * k) {
+       evMptMult = nbins[k];
+       break;
+      }
+    }
+//     printf("nrTracks = %3i   evMptMult = %.3f \n",nrTracks,evMptMult);
+
+    for (int i=0; i<nrTracks; i++) {
+      for (int j=i+1; j<nrTracks; j++) {
+       twoPartCorrPair = (tracks[i] - evMptMult) * (tracks[j] - evMptMult);
+       twoPartCorrEv = twoPartCorrEv + twoPartCorrPair;
+      }
+    }
+
+    twoPartCorrEvSq = twoPartCorrEv * twoPartCorrEv;
+    fTwoPartCorrEv->Fill(nrTracks,twoPartCorrEv);
+    fTwoPartCorrEvSq->Fill(nrTracks,twoPartCorrEvSq);
+
+//     printf("twoPartCorrEv = %.3f   twoPartCorrEvSq = %.3f \n\n",twoPartCorrEv,twoPartCorrEvSq);
+
+    // -- End analysis in multiplicity bins --
+
+
+    // -- Analysis in centrality bins --
+
+    // Centrality V0
+    if (fUseCentrality == 1) {
+
+//       printf("CentV0bin: %i  ",centralityVZEROBin);
+
+      if (centralityVZEROBin < 91) {
+
+//     printf("CentV0bin: %i   ",centralityVZEROBin);
+
+       for (int i=0; i<nrTracks; i++) {
+         for (int j=i+1; j<nrTracks; j++) {
+           twoPartCorrPairCent = (tracks[i] - evMptCent) * (tracks[j] - evMptCent);
+           twoPartCorrEvCent = twoPartCorrEvCent + twoPartCorrPairCent;
+         }
+       }
+
+       twoPartCorrEvCentSq = twoPartCorrEvCent * twoPartCorrEvCent;
+       fTwoPartCorrEvCent->Fill(centralityVZEROBin,twoPartCorrEvCent);
+       fTwoPartCorrEvCentSq->Fill(centralityVZEROBin,twoPartCorrEvCentSq);
+
+//     printf("CentV0bin: %i   evMptCent: %.3f   CorrEvCent: %.2f   CorrEvCentSq: %.2f\n",centralityVZEROBin,evMptCent,twoPartCorrEvCent,twoPartCorrEvCentSq);
+      }
+    }
+
+    // -- End analysis in centrality bins --
+
+    // --- End two-particle correlator ---
+
+
+  } // --- End calculation of various values and filling of histograms ---
+
+
+  // Post output data
+  PostData(1, fOutputList);
+}
+
+void AliAnalysisTaskPtFlucPbPb::Terminate(Option_t *)
+{
+  // Draw result to the screen
+  // Called once at the end of the query
+
+  fOutputList = dynamic_cast<TList*> (GetOutputData(1));
+  if (!fOutputList) {
+    Printf("ERROR: fOutputList not available\n");
+    return;
+  }
+
+}
diff --git a/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.h b/PWG2/EBYE/MeanPtFluctuations/AliAnalysisTaskPtFlucPbPb.h
new file mode 100644 (file)
index 0000000..193f033
--- /dev/null
@@ -0,0 +1,80 @@
+#ifndef AliAnalysisTaskPtFlucPbPb_cxx\r
+#define AliAnalysisTaskPtFlucPbPb_cxx\r
+\r
+// Analysis of Pt FLuctuations (PbPb)\r
+// Author: Stefan Heckel\r
+// Version of PbPb task: 5.0, 18.04.2011\r
+\r
+\r
+class TList;\r
+class TH1F;\r
+class TProfile;\r
+\r
+class AliESDEvent;\r
+class AliESDtrack;\r
+class AliESDtrackCuts;\r
+\r
+\r
+#include "AliAnalysisTaskSE.h"\r
+\r
+\r
+class AliAnalysisTaskPtFlucPbPb : public AliAnalysisTaskSE {\r
+  public:\r
+    AliAnalysisTaskPtFlucPbPb(const char *name = "AliAnalysisTaskPtFlucPbPb");\r
+    virtual ~AliAnalysisTaskPtFlucPbPb();\r
+\r
+    virtual void   UserCreateOutputObjects();\r
+    virtual void   UserExec(Option_t *option);\r
+    virtual void   Terminate(Option_t *);\r
+\r
+    void SetAliESDtrackCuts(AliESDtrackCuts* esdTrackCuts) {fESDTrackCuts = esdTrackCuts;}\r
+    void SetMaxVertexZ(Float_t vZ) {fMaxVertexZ = vZ;}\r
+    void SetNContributors(Int_t nCont) {fNContributors = nCont;}\r
+    void SetUseCentrality(Int_t cent) {fUseCentrality = cent;}\r
+    void SetMC(Bool_t bMC) {fMC = bMC;}\r
+\r
+\r
+  private:\r
+    AliESDEvent             *fESD;             // ESD object\r
+    TList*           fOutputList;      // List where all the output files are stored\r
+       TH1F*           fPtSpec;        // Pt spectrum\r
+       TH1F*           fMult;          // Multiplicity distribution\r
+       TH1F*           fMultSum;       // Sum of number of tracks of all events for each multiplicity bin\r
+       TH1F*           fMultSumPt;     // Sum of pTs for multiplicity bins\r
+       TH1F*           fMultNrPairs;   // Sum of number of pairs in mult. bins\r
+       TH1F*           fCent;          // Centrality distribution\r
+       TH1F*           fCentSum;       // Sum of number of tracks of all events for each centrality bin\r
+       TH1F*           fCentSumPt;     // Sum of pTs for centrality bins\r
+       TH1F*           fCentNrPairs;   // Sum of number of pairs in cent. bins\r
+       TH1F*           fEta;           // Eta distribution\r
+       TH1F*           fEtaPhiPlus;    // Phi distribution for positive eta\r
+       TH1F*           fEtaPhiMinus;   // Phi distribution for negative eta\r
+       TH1F*           fVtxZ;          // Vertex Z distribution after physics selection before any further cuts\r
+       TH1F*           fVtxZCut;       // Vertex Z dist. after vertex Z cut\r
+       TH1F*           fVtxZCont;      // Vertex Z dist. after vertex cut on nContributors\r
+       TH1F*           fVtxZTrackCuts; // Vertex Z dist. after all event and track cuts\r
+       TH1F*           fEventMeanPt;           // Event mean pT distribution\r
+       TH1F*           fEventMeanPtSq;         // Event mean pT squared dist.\r
+       TH1F*           fMultEventMeanPt;       // Event mean pT for multiplicity bins\r
+       TH1F*           fMultEventMeanPtSq;     // Event mean pT squared for mult. bins\r
+       TH1F*           fCentEventMeanPt;       // Event mean pT for centrality bins\r
+       TH1F*           fCentEventMeanPtSq;     // Event mean pT squared for cent. bins\r
+       TH1F*           fTwoPartCorrEv;         // Two-particle correlator for multiplicity bins\r
+       TH1F*           fTwoPartCorrEvSq;       // Two-part. corr. squared for mult. bins\r
+       TH1F*           fTwoPartCorrEvCent;     // Two-particle correlator for centrality bins\r
+       TH1F*           fTwoPartCorrEvCentSq;   // Two-part. corr. squared for cent. bins\r
+\r
+    AliESDtrackCuts* fESDTrackCuts;    // Esd track cuts\r
+    Float_t          fMaxVertexZ;      // Maximum value for Vertex Z position\r
+    Int_t            fNContributors;   // Minimum contributors to the vertex\r
+    Int_t            fUseCentrality;   // Use centrality (0=off, 1=VZERO, 2=SPD(not yet implemented))\r
+    Bool_t           fMC;              // Check for MC\r
+\r
+    AliAnalysisTaskPtFlucPbPb(const AliAnalysisTaskPtFlucPbPb&); // not implemented\r
+    AliAnalysisTaskPtFlucPbPb& operator=(const AliAnalysisTaskPtFlucPbPb&); // not implemented\r
+\r
+    ClassDef(AliAnalysisTaskPtFlucPbPb, 1);\r
+\r
+};\r
+\r
+#endif\r
diff --git a/PWG2/EBYE/macros/AddTaskPtFluc.C b/PWG2/EBYE/macros/AddTaskPtFluc.C
new file mode 100644 (file)
index 0000000..0a9ebed
--- /dev/null
@@ -0,0 +1,108 @@
+AliAnalysisTask *AddTaskPtFluc(){
+  //get the current analysis manager
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+    Error("AddTask_sheckel_PtFluc", "No analysis manager found.");
+    return 0;
+  }
+
+
+  // --- Check for MC ---
+  AliMCEventHandler  *mcH = static_cast<AliMCEventHandler*>(mgr->GetMCtruthEventHandler());
+
+
+  //============= Set Task Name ===================
+  TString taskName=("AliAnalysisTaskPtFluc.cxx+");
+  //===============================================
+  //            Load the task
+  gROOT->LoadMacro(taskName.Data());
+  if (gProof){
+    TString taskSO=gSystem->pwd();
+    taskSO+="/";
+    taskSO+=taskName(0,taskName.First('.'))+"_cxx.so";
+    gProof->Exec(Form("gSystem->Load(\"%s\")",taskSO.Data()),kTRUE);
+  }
+
+
+  //========= Add task to the ANALYSIS manager =====
+  AliAnalysisTaskPtFluc *task = new AliAnalysisTaskPtFluc("sheckelTaskPtFluc");
+
+
+  // --- Set ESD track Cuts ---
+
+  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("esdTrackCuts"); // My pp track cuts (ITS+TPC)
+  AliESDtrackCuts *esdCuts = new AliESDtrackCuts("esdCuts"); // TPC-only track cuts
+
+  // -- My pp track cuts (ITS+TPC)
+  // TPC
+  esdTrackCuts->SetRequireTPCStandAlone(kFALSE);
+  esdTrackCuts->SetMinNClustersTPC(70);
+  esdTrackCuts->SetMaxChi2PerClusterTPC(4.0);
+  esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  // ITS
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kAny);
+  // DCA to vertex
+  esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01"); // 7*(0.0026+0.0050/pt^1.01)
+  esdTrackCuts->SetMaxDCAToVertexZ(2.);
+  esdTrackCuts->SetDCAToVertex2D(kFALSE);
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  // Pt range and acceptance in eta
+  esdTrackCuts->SetPtRange(0.15,2.);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+
+
+  // -- TPC-only track cuts
+  // TPC
+  esdCuts->SetRequireTPCStandAlone(kTRUE);
+  esdCuts->SetMinNClustersTPC(70);
+  esdCuts->SetMaxChi2PerClusterTPC(4.0);
+  esdCuts->SetAcceptKinkDaughters(kFALSE);
+  esdCuts->SetRequireTPCRefit(kFALSE);
+  // ITS
+  esdCuts->SetRequireITSRefit(kFALSE);
+  // DCA to vertex
+  esdCuts->SetMaxDCAToVertexXY(2.4); // cm
+  esdCuts->SetMaxDCAToVertexZ(3.2);  // cm
+  esdCuts->SetDCAToVertex2D(kTRUE);
+  esdCuts->SetRequireSigmaToVertex(kFALSE);
+  // Pt range and acceptance in eta
+  esdCuts->SetPtRange(0.15,2.);
+  esdCuts->SetEtaRange(-0.8,0.8);
+
+  // --- End track Cuts ---
+
+
+  task->SetAliESDtrackCuts(esdTrackCuts);
+//   task->SetAliESDtrackCuts(esdCuts);
+
+  task->SelectCollisionCandidates(AliVEvent::kMB);
+  task->SetMaxVertexZ(10.);    // cm
+  task->SetNContributors(1);
+
+  if (mcH) task->SetMC(kTRUE);
+
+
+  mgr->AddTask(task);
+
+
+  //================================================
+  //              data containers
+  //================================================
+
+  // input container
+  AliAnalysisDataContainer *cinput  = mgr->GetCommonInputContainer();
+
+  // output container
+  AliAnalysisDataContainer *coutput =
+      mgr->CreateContainer("sheckel_PtFluc", TList::Class(),
+                           AliAnalysisManager::kOutputContainer,"sheckel_PtFluc.root");
+
+  // connect containers
+  mgr->ConnectInput  (task,  0, cinput );
+  mgr->ConnectOutput (task,  1, coutput );
+
+  return task;
+
+}
diff --git a/PWG2/EBYE/macros/AddTaskPtFlucPbPb.C b/PWG2/EBYE/macros/AddTaskPtFlucPbPb.C
new file mode 100644 (file)
index 0000000..27b6e8a
--- /dev/null
@@ -0,0 +1,109 @@
+AliAnalysisTask *AddTaskPtFlucPbPb(){
+  //get the current analysis manager
+  AliAnalysisManager *mgr = AliAnalysisManager::GetAnalysisManager();
+  if (!mgr) {
+    Error("AddTask_sheckel_PtFlucPbPb", "No analysis manager found.");
+    return 0;
+  }
+
+
+  // --- Check for MC ---
+  AliMCEventHandler  *mcH = static_cast<AliMCEventHandler*>(mgr->GetMCtruthEventHandler());
+
+
+  //============= Set Task Name ===================
+  TString taskName=("AliAnalysisTaskPtFlucPbPb.cxx+");
+  //===============================================
+  //            Load the task
+  gROOT->LoadMacro(taskName.Data());
+  if (gProof){
+    TString taskSO=gSystem->pwd();
+    taskSO+="/";
+    taskSO+=taskName(0,taskName.First('.'))+"_cxx.so";
+    gProof->Exec(Form("gSystem->Load(\"%s\")",taskSO.Data()),kTRUE);
+  }
+
+
+  //========= Add task to the ANALYSIS manager =====
+  AliAnalysisTaskPtFlucPbPb *task = new AliAnalysisTaskPtFlucPbPb("sheckelTaskPtFlucPbPb");
+
+
+  // --- Set ESD track Cuts ---
+
+  AliESDtrackCuts *esdTrackCuts = new AliESDtrackCuts("esdTrackCuts"); // My pp track cuts (ITS+TPC)
+  AliESDtrackCuts *esdCuts = new AliESDtrackCuts("esdCuts"); // TPC-only track cuts
+
+  // -- My pp track cuts (ITS+TPC)
+  // TPC
+  esdTrackCuts->SetRequireTPCStandAlone(kFALSE);
+  esdTrackCuts->SetMinNClustersTPC(70);
+  esdTrackCuts->SetMaxChi2PerClusterTPC(4.0);
+  esdTrackCuts->SetAcceptKinkDaughters(kFALSE);
+  esdTrackCuts->SetRequireTPCRefit(kTRUE);
+  // ITS
+  esdTrackCuts->SetRequireITSRefit(kTRUE);
+  esdTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kAny);
+  // DCA to vertex
+  esdTrackCuts->SetMaxDCAToVertexXYPtDep("0.0182+0.0350/pt^1.01"); // 7*(0.0026+0.0050/pt^1.01)
+  esdTrackCuts->SetMaxDCAToVertexZ(2.);
+  esdTrackCuts->SetDCAToVertex2D(kFALSE);
+  esdTrackCuts->SetRequireSigmaToVertex(kFALSE);
+  // Pt range and acceptance in eta
+  esdTrackCuts->SetPtRange(0.15,2.);
+  esdTrackCuts->SetEtaRange(-0.8,0.8);
+
+
+  // -- TPC-only track cuts
+  // TPC
+  esdCuts->SetRequireTPCStandAlone(kTRUE);
+  esdCuts->SetMinNClustersTPC(70);
+  esdCuts->SetMaxChi2PerClusterTPC(4.0);
+  esdCuts->SetAcceptKinkDaughters(kFALSE);
+  esdCuts->SetRequireTPCRefit(kFALSE);
+  // ITS
+  esdCuts->SetRequireITSRefit(kFALSE);
+  // DCA to vertex
+  esdCuts->SetMaxDCAToVertexXY(2.4); // cm
+  esdCuts->SetMaxDCAToVertexZ(3.2);  // cm
+  esdCuts->SetDCAToVertex2D(kTRUE);
+  esdCuts->SetRequireSigmaToVertex(kFALSE);
+  // Pt range and acceptance in eta
+  esdCuts->SetPtRange(0.15,2.);
+  esdCuts->SetEtaRange(-0.8,0.8);
+
+  // --- End track Cuts ---
+
+
+//   task->SetAliESDtrackCuts(esdTrackCuts);
+  task->SetAliESDtrackCuts(esdCuts);
+
+  task->SelectCollisionCandidates(AliVEvent::kMB);
+  task->SetMaxVertexZ(10.);    // cm
+  task->SetNContributors(1);
+  task->SetUseCentrality(1);   // 0=off, 1=VZERO, 2=SPD(not implemented)
+
+  if (mcH) task->SetMC(kTRUE);
+
+
+  mgr->AddTask(task);
+
+
+  //================================================
+  //              data containers
+  //================================================
+
+  // input container
+  AliAnalysisDataContainer *cinput  = mgr->GetCommonInputContainer();
+
+  // output container
+  AliAnalysisDataContainer *coutput =
+      mgr->CreateContainer("sheckel_PtFlucPbPb", TList::Class(),
+                           AliAnalysisManager::kOutputContainer,"sheckel_PtFlucPbPb.root");
+
+  // connect containers
+  mgr->ConnectInput  (task,  0, cinput );
+  mgr->ConnectOutput (task,  1, coutput );
+
+  return task;
+
+}
index 28320a2..14a7cad 100644 (file)
@@ -17,4 +17,7 @@
 
 #pragma link C++ class AliEbyEFluctuationAnalysisTask+;
 
+#pragma link C++ class AliAnalysisTaskPtFluc+;
+#pragma link C++ class AliAnalysisTaskPtFlucPbPb+;
+
 #endif