AliPhysics  cdeda5a (cdeda5a)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliJetEmbeddingFromPYTHIATask.cxx
Go to the documentation of this file.
1 // $Id$
2 //
3 // Jet embedding from PYTHIA task.
4 //
5 // Author: S.Aiola, C.Loizides
6 
8 
9 #include <TFile.h>
10 #include <TMath.h>
11 #include <TString.h>
12 #include <TRandom.h>
13 #include <TParameter.h>
14 #include <TH1I.h>
15 #include <TGrid.h>
16 #include <THashTable.h>
17 #include <TSystem.h>
18 
19 #include "AliVEvent.h"
20 #include "AliLog.h"
21 
23 
24 //________________________________________________________________________
27  fPYTHIAPath(),
28  fPtHardBinScaling(),
29  fLHC11hAnchorRun(kTRUE),
30  fAnchorRun(-1),
31  fFileTable(0),
32  fUseAsVetoTable(kTRUE),
33  fMinEntriesPerPtHardBin(1),
34  fCurrentPtHardBin(-1),
35  fPtHardBinParam(0),
36  fPtHardBinCount(0),
37  fHistPtHardBins(0)
38 {
39  // Default constructor.
40  SetSuffix("PYTHIAEmbedding");
41  fTotalFiles = 140;
42  fRandomAccess = kTRUE;
43  SetAODMC(kTRUE);
44 }
45 
46 //________________________________________________________________________
48  AliJetEmbeddingFromAODTask(name, drawqa),
49  fPYTHIAPath("alien:///alice/sim/2012/LHC12a15e_fix/%d/%d/AOD149/%04d/AliAOD.root"),
50  fPtHardBinScaling(),
51  fLHC11hAnchorRun(kTRUE),
52  fAnchorRun(-1),
53  fFileTable(0),
54  fUseAsVetoTable(kTRUE),
55  fMinEntriesPerPtHardBin(1),
56  fCurrentPtHardBin(-1),
57  fPtHardBinParam(0),
58  fPtHardBinCount(0),
59  fHistPtHardBins(0)
60 {
61  // Standard constructor.
62  SetSuffix("PYTHIAEmbedding");
63  fTotalFiles = 140;
64  fRandomAccess = kTRUE;
65  SetAODMC(kTRUE);
66 }
67 
68 //________________________________________________________________________
70 {
71  // Destructor
72 }
73 
74 //________________________________________________________________________
76 {
77  if (!fQAhistos)
78  return;
79 
81 
82  fHistPtHardBins = new TH1F("fHistPtHardBins", "fHistPtHardBins", 11, 0, 11);
83  fHistPtHardBins->GetXaxis()->SetTitle("p_{T} hard bin");
84  fHistPtHardBins->GetYaxis()->SetTitle("total events");
86 
87  const Int_t ptHardLo[11] = { 0, 5,11,21,36,57, 84,117,152,191,234};
88  const Int_t ptHardHi[11] = { 5,11,21,36,57,84,117,152,191,234,1000000};
89 
90  for (Int_t i = 1; i < 12; i++)
91  fHistPtHardBins->GetXaxis()->SetBinLabel(i, Form("%d-%d",ptHardLo[i-1],ptHardHi[i-1]));
92 
93  fHistEmbeddingQA = new TH1F("fHistEmbeddingQA", "fHistEmbeddingQA", 2, 0, 2);
94  fHistEmbeddingQA->GetXaxis()->SetTitle("Event state");
95  fHistEmbeddingQA->GetYaxis()->SetTitle("counts");
96  fHistEmbeddingQA->GetXaxis()->SetBinLabel(1, "OK");
97  fHistEmbeddingQA->GetXaxis()->SetBinLabel(2, "Not embedded");
99 
100  fHistRejectedEvents = new TH1F("fHistRejectedEvents", "fHistRejectedEvents", 500, -0.5, 499.5);
101  fHistRejectedEvents->GetXaxis()->SetTitle("# of rejected events");
102  fHistRejectedEvents->GetYaxis()->SetTitle("counts");
104 
105  PostData(1, fOutput);
106 }
107 
108 //________________________________________________________________________
110 {
111  if (fPtHardBinScaling.GetSize() > 0) {
112  Double_t sum = 0;
113  for (Int_t i = 0; i < fPtHardBinScaling.GetSize(); i++)
114  sum += fPtHardBinScaling[i];
115 
116  if (sum == 0) {
117  AliWarning("No hard pt bin scaling!");
118  sum = fPtHardBinScaling.GetSize();
119  for (Int_t i = 0; i < fPtHardBinScaling.GetSize(); i++)
120  fPtHardBinScaling[i] /= sum;
121  }
122  else {
123  for (Int_t i = 0; i < fPtHardBinScaling.GetSize(); i++)
124  fPtHardBinScaling[i] /= sum;
125  }
126  }
127 
128  fPtHardBinParam = static_cast<TParameter<int>*>(InputEvent()->FindListObject("PYTHIAPtHardBin"));
129  if (!fPtHardBinParam) {
130  fPtHardBinParam = new TParameter<int>("PYTHIAPtHardBin", 0);
131  AliDebug(3,"Adding pt hard bin param object to the event list...");
132  InputEvent()->AddObject(fPtHardBinParam);
133  }
134 
136 }
137 
138 //________________________________________________________________________
140 {
142  fPtHardBinCount = 0;
143 
144  Int_t newPtHard = GetRandomPtHardBin();
145 
146  if (newPtHard != fCurrentPtHardBin) {
147  fPtHardBinParam->SetVal(newPtHard);
148  fCurrentPtHardBin = newPtHard;
149  if (!OpenNextFile()) return kFALSE;
150  }
151  }
152 
153  fPtHardBinCount++;
154  if (fHistPtHardBins) fHistPtHardBins->SetBinContent(fCurrentPtHardBin+1, fHistPtHardBins->GetBinContent(fCurrentPtHardBin+1)+1);
155 
157 }
158 
159 //________________________________________________________________________
161 {
162  static Int_t order[20]={-1};
163 
164  if (order[0] == -1)
165  TMath::Sort(fPtHardBinScaling.GetSize(), fPtHardBinScaling.GetArray(), order);
166 
167  Double_t rnd = gRandom->Rndm();
168  Double_t sum = 0;
169  Int_t ptHard = -1;
170  for (Int_t i = 0; i < fPtHardBinScaling.GetSize(); i++) {
171  sum += fPtHardBinScaling[order[i]];
172  if (sum >= rnd) {
173  ptHard = order[i];
174  break;
175  }
176  }
177 
178  return ptHard;
179 }
180 
181 //________________________________________________________________________
183 {
184  if (!fLHC11hAnchorRun)
185  return kTRUE;
186 
187  Int_t runNumber = InputEvent()->GetRunNumber();
188 
189  Int_t semiGoodRunList[32] = {169975, 169981, 170038, 170040, 170083, 170084, 170085,
190  170088, 170089, 170091, 170152, 170155, 170159, 170163,
191  170193, 170195, 170203, 170204, 170228, 170230, 170268,
192  170269, 170270, 170306, 170308, 170309, 169238, 169160,
193  169156, 169148, 169145, 169144};
194 
195  fAnchorRun = 169838; // Assume it is a good run
196 
197  for (Int_t i = 0; i < 32; i++) {
198  if (runNumber == semiGoodRunList[i]) { // If it is semi good, change the anchor run
199  fAnchorRun = 170040;
200  break;
201  }
202  }
203 
204  return kTRUE;
205 }
206 
207 //________________________________________________________________________
209 {
210  fCurrentAODFileID = TMath::Nint(gRandom->Rndm()*(fTotalFiles-1))+1;
211 
212  if (fMinEntriesPerPtHardBin < 0) {
215  }
216 
218 
219  if (fAnchorRun>0)
220  fileName = Form(fPYTHIAPath.Data(), fAnchorRun, fCurrentPtHardBin, fCurrentAODFileID);
221  else
222  fileName = Form(fPYTHIAPath.Data(), fCurrentPtHardBin, fCurrentAODFileID);
223 
224  if (fFileTable && fFileTable->GetEntries() > 0) {
225  TObject* obj = fFileTable->FindObject(fileName);
226  if (obj != 0 && fUseAsVetoTable) {
227  AliWarning(Form("File %s found in the vetoed file table. Skipping...", fileName.Data()));
228  return 0;
229  }
230  if (obj == 0 && !fUseAsVetoTable) {
231  AliWarning(Form("File %s not found in the allowed file table. Skipping...", fileName.Data()));
232  return 0;
233  }
234  }
235 
236  if (fileName.BeginsWith("alien://") && !gGrid) {
237  AliInfo("Trying to connect to AliEn ...");
238  TGrid::Connect("alien://");
239  }
240 
241  TString baseFileName(fileName);
242  if (baseFileName.Contains(".zip#")) {
243  Ssiz_t pos = baseFileName.Last('#');
244  baseFileName.Remove(pos);
245  }
246 
247  if (gSystem->AccessPathName(baseFileName)) {
248  AliError(Form("File %s does not exist!", baseFileName.Data()));
249  return 0;
250  }
251 
252  AliDebug(3,Form("Trying to open file %s...", fileName.Data()));
253  TFile *file = TFile::Open(fileName);
254 
255  if (!file || file->IsZombie()) {
256  AliError(Form("Unable to open file: %s!", fileName.Data()));
257  return 0;
258  }
259 
260  return file;
261 }
double Double_t
Definition: External.C:58
TString fileName
TList * fOutput
! output list for QA histograms
TSystem * gSystem
TH1 * fHistEmbeddingQA
File ID not embedded.
Class for Pythia event embedding into data.
Bool_t ExecOnce()
generate a particle with random eta,phi, and correlated pt,mass values
TRandom * gRandom
TH1 * fHistPtHardBins
Number of event embedded from the current pt hard bin.
Bool_t ExecOnce()
generate a particle with random eta,phi, and correlated pt,mass values
ClassImp(AliJetEmbeddingFromPYTHIATask) AliJetEmbeddingFromPYTHIATask
int Int_t
Definition: External.C:63
TParameter< int > * fPtHardBinParam
Pt hard bin of the current open file.
Bool_t fQAhistos
draw QA histograms
TFile * file
Int_t fCurrentAODFileID
Current file being processed (via the event handler)
void SetSuffix(const char *s)
bool Bool_t
Definition: External.C:53
Class for embedding a AOD event into a data event.