AliPhysics  1adf5bd (1adf5bd)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
FastSim.C
Go to the documentation of this file.
1 #ifndef FASTSIM_H
2 #define FASTSIM_H
3 #include <TSelector.h>
4 #include <TQObject.h>
5 #ifndef __CINT__
6 # include "AliGenerator.h"
7 # include "AliRunLoader.h"
8 # include "AliStack.h"
9 # include "AliHeader.h"
10 # include "AliGenEventHeader.h"
11 # include "AliRun.h"
12 # include "AliCollisionGeometry.h"
13 # include "AliGenPythiaEventHeader.h"
14 # include "AliGenDPMjetEventHeader.h"
15 # include "AliGenGeVSimEventHeader.h"
16 # include "AliGenHerwigEventHeader.h"
17 # include "AliGenHijingEventHeader.h"
18 # include "AliGenCocktailEventHeader.h"
19 # include "AliPDG.h"
20 # include <TROOT.h>
21 # include <TString.h>
22 # include <TMath.h>
23 # include <TParticle.h>
24 # include <TH1.h>
25 # include <TTree.h>
26 # include <TClonesArray.h>
27 # include <TList.h>
28 # include <TProof.h>
29 # include <TChain.h>
30 # include <TParticlePDG.h>
31 # include <TStopwatch.h>
32 # include <TFile.h>
33 # include <TProofOutputFile.h>
34 # include <TCanvas.h>
35 # include <TTimer.h>
36 # include <TRandom.h>
37 # include <TUrl.h>
38 # include <TMacro.h>
39 # include <TSystemDirectory.h>
40 # include <TFileCollection.h>
41 # include <TPRegexp.h>
42 # include <fstream>
43 # include <TLeaf.h>
44 # include <algorithm>
45 # include "FastShortHeader.C"
46 # include "FastCentEstimators.C"
47 # include "FastMonitor.C"
48 #else
49 class AliGenerator;
50 class AliRunLoader;
51 class AliStack;
52 class AliHeader;
53 class AliGenEventHeader;
54 class TH1;
55 class TTree;
56 class TClonesArray;
57 class TBrowser;
58 class TList;
59 class TFile;
60 class TProofOutputFile;
61 class TCanvas;
62 class TVirtualPad;
63 class TTimer;
64 class TUrl;
65 class TAxis;
66 class TParticle;
67 class TMacro;
68 class TLeaf;
69 class BCentEstimator;
70 class FastMonitor;
71 #endif
72 
73 
74 
76 typedef struct {
77  Double_t rproj;
78  Double_t rtarg;
79  Double_t bimpac;
80  Int_t nwtsam;
81  Int_t nwasam;
82  Int_t nwbsam;
83  Int_t nwtacc;
84  Int_t nwtaac;
85  Int_t nwtbac;
86  Int_t ncp;
87  Int_t nct;
88 } DtglcpCommon;
89 DtglcpCommon* _dtglcp = 0;
90 
91 //====================================================================
92 struct FastSimMonitor : public FastMonitor
93 {
94  FastSimMonitor(TSelector* s=0)
95  : FastMonitor(s, "FastMonitor")
96  {
97  Register("list/histograms/b", "hist", 0, false);
98  Register("list/histograms/cent", "hist", 0, false);
99  Register("list/histograms/timing", "hist", 0x18, false);
100  Register("list/histograms/dNdeta", "", 0x18, false);
101  Register("list/histograms/dNdy", "", 0x18, false);
102  Register("list/estimators/rawV0MP", "e", 0x02, false);
103  if (gProof) gProof->AddFeedback("list");
104  }
105 };
106 
107 //====================================================================
111 struct FastSim : public TSelector
112 {
123  FastSim(const char* eg="",
124  ULong_t runNo=0,
125  Double_t bMin=0,
126  Double_t bMax=20,
127  Long64_t nEvents=0,
128  Int_t monitor=0)
129  : TSelector(),
130  fEGName(eg),
131  fRunNo(runNo),
132  fBMin(bMin),
133  fBMax(bMax),
134  fGRP(0),
135  fOverrides(0),
136  fNEvents(nEvents),
137  fIsTgtA(false),
138  fIsProjA(false),
139  fGenerator(0),
140  fRunLoader(0),
141  fStack(0),
142  fHeader(0),
143  fTree(0),
144  fParticles(0),
145  fList(0),
146  fHEta(0),
147  fHY(0),
148  fHIpz(0),
149  fHType(0),
150  fHCent(0),
151  fHB(0),
152  fHPhiR(0),
153  fHTime(0),
154  fBEstimator(0),
155  fCentEstimators(0),
156  fProofFile(0),
157  fAliceFile(0),
158  fKineFile(0),
159  fFile(0),
160  fFileName(""),
161  fVerbose(true),
162  fMonitor(monitor)
163  {}
164  const char* FileName() const
165  {
166  static TString fn;
167  if (fn.IsNull()) {
168  if (!fFileName.IsNull()) fn = fFileName;
169  else {
170  TString egName = (fGenerator ? fGenerator->GetName() : "");
171  if (egName.IsNull()) egName = fEGName;
172 
173  fn = Form("%s_%09d", egName.Data(), fRunNo);
174  if (fGenerator) {
175  TString tgt, proj;
176  Int_t tgtA, tgtZ, projA, projZ;
177  fGenerator->GetTarget(tgt, tgtA, tgtZ);
178  fGenerator->GetProjectile(proj, projA, projZ);
179  fn.Append(Form("_%s%s", tgt.Data(), proj.Data()));
180  fn.Append(Form("_%05d", Int_t(fGenerator->GetEnergyCMS())));
181  }
182  else {
183  Long_t en = gROOT->ProcessLine("grp->energy");
184  // Bool_t a1 = gROOT->ProcessLine("(GRPData*)%p->beam1.IsA()",fGRP);
185  // Bool_t a2 = gROOT->ProcessLine("(GRPData*)%p->beam2.IsA()",fGRP);
186  fn.Append(Form("_%s%s", (fIsTgtA ? "A" : "p"),
187  (fIsProjA ? "A" : "p")));
188  fn.Append(Form("_%05ld", (fGRP ? en : 0)));
189  }
190 
191  if (fNEvents > 0) {
192  if (fNEvents >= 1000000)
193  fn.Append(Form("_%lldM", fNEvents/1000000));
194  else if (fNEvents >= 1000)
195  fn.Append(Form("_%lldk", fNEvents/1000));
196  else
197  fn.Append(Form("_%lld", fNEvents));
198  }
199  fn.Append(".root");
200  fFileName = fn;
201  }
202  }
203  return fn.Data();
204  }
210  const char* GetName() const { return "FastSim"; }
216  const char* GetTitle() const { return "ALICE Event Generator simulation"; }
222  virtual UShort_t GetSNN() const
223  {
224  Float_t fsNN = 0;
225  if (fGenerator) fsNN = fGenerator->GetEnergyCMS();
226  else fsNN = gROOT->ProcessLine("grp->energy");
227 
228  UShort_t sNN = (TMath::Abs(fsNN - 2760) < 10 ? 2760 :
229  TMath::Abs(fsNN - 5023) < 10 ? 5023 :
230  TMath::Abs(fsNN - 2360) < 10 ? 2360 :
231  TMath::Abs(fsNN - 900) < 10 ? 900 :
232  TMath::Abs(fsNN - 7000) < 10 ? 7000 :
233  TMath::Abs(fsNN - 8000) < 10 ? 8000 :
234  TMath::Abs(fsNN - 13000) < 10 ? 13000 :
235  0);
236  return sNN;
237  }
243  virtual const char* GetEGTitle() const
244  {
245 
246  if (fGenerator) return fGenerator->GetTitle();
247  return "";
248  }
254  Bool_t SetupOutput()
255  {
256  Printf("=== Output =============================\n"
257  " File name: %s\n"
258  "========================================", FileName());
259 
260  if (fVerbose) Info("SetupOutput", "First the file");
261  Bool_t isProof = false;
262  if (fInput && fInput->FindObject("PROOF_Ordinal"))
263  isProof = true;
264  if (isProof) {
265  Info("SetupOutput", "Making Proof File");
266  fProofFile = new TProofOutputFile(FileName(), "M");
267  // TProofOutputFile::kMerge,
268  // TProofOutputFile::kRemote);
269  fFile = fProofFile->OpenFile("RECREATE");
270  }
271  else
272  fFile = TFile::Open(FileName(), "RECREATE");
273 
274  UShort_t sNN = GetSNN();
275  TString tit = GetEGTitle();
276  if (fVerbose) Info("SetupOutput", "Making our tree: %s", tit.Data());
277  fTree = new TTree("T", tit.Data());
278  fParticles = new TClonesArray("TParticle");
279  fTree->Branch("header", &fShortHead,
280  "run/i:event:ntgt:nproj:nbin:type:"
281  "ipx/D:ipy:ipz:b:c:phir:"
282  "nsnp/i:nsnt:nspp:nspt");
283  fTree->Branch("particles", &fParticles);
284  fTree->SetAutoSave(500); // Save every on every 100 events
285  fTree->SetDirectory(fFile);
286  fTree->SetAlias("primary", "(particles.fBits&(1<<14))");
287  fTree->SetAlias("weak", "(particles.fBits&(1<<15))");
288  fTree->SetAlias("charged", "(particles.fBits&(1<<16))");
289  fTree->SetAlias("pt", "(sqrt(pow(particles.fPx,2)+"
290  /* */"pow(particles.fPy,2)))");
291  fTree->SetAlias("eta", "(particles.Pt()<1e-100?"
292  "(particles.fPz>0?100:-100):"
293  "-log(tan(atan2(particles.Pt(),particles.fPz)/2)))");
294  fTree->SetAlias("good", "(primary&&charged&&abs(eta)<1000)");
295  fTree->SetAlias("sd", "(header.fType & 0x1)");
296  fTree->SetAlias("dd", "(header.fType & 0x2)");
297  fTree->SetAlias("pion", "(abs(particles.fPdgCode)==211)");
298  fTree->SetAlias("kaon", "(abs(particles.fPdgCode)==321)");
299  fTree->SetAlias("proton", "(abs(particles.fPdgCode)==2212)");
300  fTree->SetAlias("electron","(abs(particles.fPdgCode)==11)");
301  fTree->SetAlias("other", "(!pion&&!kaon&&!proton&&!electron)");
302  fTree->SetAlias("beta", "(particles.P()/particle.Energy())");
303  fTree->SetAlias("gamma", "(1./sqrt(1-beta*beta))");
304  fTree->SetAlias("npart", "(header.ntgt+header.nproj)");
305 
306  // Let's add the title of the generator to the output. We make a
307  // histogram because that can be merged.
308  Info("SetupOutput", "Making generator histogram: %s",tit.Data());
309  TH1* egTitle = new TH1I("eg", tit.Data(), 1, 0, 1);
310  egTitle->SetDirectory(0);
311  egTitle->Fill(0.5);
312  egTitle->SetXTitle("N_{\\hbox{workers}}");
313  egTitle->SetYTitle("Count");
314  egTitle->SetFillColor(kYellow+4);
315  egTitle->SetLineColor(kYellow+4);
316  egTitle->SetMarkerColor(kYellow+4);
317  egTitle->SetFillStyle(1001);
318  egTitle->SetMarkerStyle(1);
319  egTitle->SetLineStyle(1);
320  egTitle->SetLineWidth(2);
321  egTitle->SetMarkerSize(1);
322 
323 
324  if (fVerbose) Info("SetupOutput", "Making histograms");
325  Double_t maxEta = 10;
326  Double_t dEta = 10./200;
327  fHEta = new TH1D("dNdeta", "Charged particle pseudo-rapidity density",
328  Int_t(2*maxEta/dEta+.5), -maxEta, +maxEta);
329  fHEta->Sumw2();
330  fHEta->SetXTitle("\\eta");
331  fHEta->SetYTitle("\\hbox{d}N_{\\hbox{ch}}/\\hbox{d}\\eta");
332  fHEta->SetMarkerColor(kRed+2);
333  fHEta->SetMarkerStyle(20);
334  fHEta->SetDirectory(0);
335 
336  fHY = new TH1D("dNdy", "Charged particle rapidity density",
337  Int_t(2*maxEta/dEta+.5), -maxEta, +maxEta);
338  fHY->Sumw2();
339  fHY->SetXTitle("\\mathit{y}");
340  fHY->SetYTitle("\\hbox{d}N_{\\hbox{ch}}/\\hbox{d}y");
341  fHY->SetMarkerColor(kRed+2);
342  fHY->SetMarkerStyle(20);
343  fHY->SetDirectory(0);
344 
345  fHIpz = new TH1D("ipZ", "Z-coordinate of interaction point",
346  10, -10, 10);
347  fHIpz->SetMarkerColor(kGreen+2);
348  fHIpz->SetFillColor(kGreen+2);
349  fHIpz->SetFillStyle(3001);
350  fHIpz->SetXTitle("IP_{#it{z}} [cm]");
351  fHIpz->SetDirectory(0);
352 
353  fHType = new TH1D("type", "Diffractive", 3, .5, 3.5);
354  fHType->SetMarkerColor(kOrange+2);
355  fHType->SetFillColor(kOrange+2);
356  fHType->SetFillStyle(3001);
357  fHType->SetDirectory(0);
358  fHType->GetXaxis()->SetBinLabel(1, "Non");
359  fHType->GetXaxis()->SetBinLabel(2, "Single");
360  fHType->GetXaxis()->SetBinLabel(3, "Double");
361 
362  fHCent = new TH1D("cent", "Centrality", 20, 0, 100);
363  fHCent->SetMarkerColor(kPink+2);
364  fHCent->SetFillColor(kPink+2);
365  fHCent->SetFillStyle(3001);
366  fHCent->SetDirectory(0);
367  fHCent->SetYTitle("Events");
368  fHCent->SetXTitle("Centrality [%]");
369 
370  fHB = new TH1D("b", "Impact parameter", 20, 0, 20);
371  fHB->SetMarkerColor(kYellow+2);
372  fHB->SetFillColor(kYellow+2);
373  fHB->SetFillStyle(3001);
374  fHB->SetYTitle("Events");
375  fHB->SetXTitle("#it{b} [fm]");
376  fHB->SetDirectory(0);
377 
378  fHPhiR = new TH1D("phiR", "Event plane angle", 360, 0, 360);
379  fHPhiR->SetMarkerColor(kMagenta+2);
380  fHPhiR->SetFillColor(kMagenta+2);
381  fHPhiR->SetFillStyle(3001);
382  fHPhiR->SetXTitle("#it{#Phi}_{R} [degrees]");
383  fHPhiR->SetDirectory(0);
384 
385  fHTime = new TH1D("timing", "Timing of processing", 5,0.5,5.5);
386  fHTime->SetMarkerColor(kBlue+2);
387  fHTime->SetFillColor(kBlue+2);
388  fHTime->SetFillStyle(3001);
389  fHTime->SetYTitle("seconds / event");
390  fHTime->GetXaxis()->SetBinLabel(1,"Reset");
391  fHTime->GetXaxis()->SetBinLabel(2,"Generate");
392  fHTime->GetXaxis()->SetBinLabel(3,"Header");
393  fHTime->GetXaxis()->SetBinLabel(4,"Particles");
394  fHTime->GetXaxis()->SetBinLabel(5,"Filling");
395  fHTime->SetDirectory(0);
396 
397  fList = new TList;
398  fList->SetName("list");
399  fList->Add(egTitle);
400 
401  TList* histos = new TList;
402  histos->SetName("histograms");
403  histos->SetOwner(true);
404  histos->Add(fHEta);
405  histos->Add(fHY);
406  histos->Add(fHIpz);
407  histos->Add(fHType);
408  histos->Add(fHCent);
409  histos->Add(fHB);
410  histos->Add(fHPhiR);
411  histos->Add(fHTime);
412  fList->Add(histos);
413 
414  TList* estimators = new TList;
415  estimators->SetName("estimators");
416  estimators->SetOwner(true);
417  fList->Add(estimators);
418 
419 
420  TIter next(fCentEstimators);
421  FastCentEstimator* estimator = 0;
422  while ((estimator = static_cast<FastCentEstimator*>(next()))) {
423  estimator->Setup(estimators, fTree,sNN,fIsTgtA,fIsProjA);
424  estimator->SetVerbose(fVerbose);
425  // estimator->Print("nah");
426  }
427 
428  if (fVerbose) Info("SetupOutput", "Adding list ot outputs");
429  fOutput->Add(fList);
430  // fOutput->ls();
431 
432  return true;
433  }
439  virtual void SetupSeed()
440  {
441  std::ifstream in("/dev/urandom");
442  UInt_t seed = 0;
443  in.read(reinterpret_cast<char*>(&seed), sizeof(seed));
444  in.close();
445  Printf("=== Random =============================\n"
446  " Random number seed: %d\n"
447  "========================================", seed);
448  gRandom->SetSeed(seed);
449  }
450  virtual Bool_t SetupGRP()
451  {
452  Printf(" === Setup ==============================");
453  Printf(" Run #: %6d", fRunNo);
454  Printf(" EG: %30s", fEGName.Data());
455  Printf(" B range: %5.1ffm - %5.1ffm", fBMin, fBMax);
456  Printf(" ========================================");
457  Printf("Macro path: %s", gROOT->GetMacroPath());
458 
459  // --- Check if we shoud get the GRP line ------------------------
460  if (!fGRP && fInput) {
461  fGRP = fInput->FindObject("GRP");
462  std::ofstream* pout = new std::ofstream("grp.dat");
463  if (pout) {
464  if (fVerbose)
465  Info("SetupGRP", "Writing GRP line '%s' to \"grp.dat\"",
466  fGRP->GetTitle());
467  std::ostream& out = *pout;
468  out << fGRP->GetTitle() << std::endl;
469  pout->close();
470  }
471  }
472  if (fVerbose)
473  Info("SetupGRP", "Overrides: %p Input: %p", fOverrides, fInput);
474  if (!fOverrides && fInput) {
475  fOverrides = static_cast<TList*>(fInput->FindObject("overrides"));
476  if (!fOverrides && fVerbose)
477  Info("SetupGRP", "No GRP overrides found in input:");
478  }
479 
480  // --- Load our settings -----------------------------------------
481  if (fVerbose) Info("SetupGRP", "Loading scripts");
482  // Check if we have the global "grp" already
483  if (gROOT->ProcessLine("grp") == 0)
484  gROOT->Macro(Form("GRP.C(%d)", fRunNo));
485  if (fVerbose) Info("SetupGRP", "Perhaps override");
486  OverrideGRP();
487  gROOT->ProcessLine("grp->Print()");
488  return true;
489  }
495  virtual Bool_t SetupGen()
496  {
497  if (fVerbose) Info("SetupGen", "Load base config");
498  gROOT->Macro("BaseConfig.C");
499  if (fVerbose) Info("SetupGen", "Load EG config");
500  gROOT->Macro("EGConfig.C");
501 
502  gROOT->ProcessLine(Form("VirtualEGCfg::LoadGen(\"%s\")",fEGName.Data()));
503 
504  // --- Make our generator ----------------------------------------
505  // Info("SetupGen", "Creating generator");
506  TString egMk = Form("egCfg->MakeGenerator(\"%s\",%f,%f)",
507  fEGName.Data(), fBMin, fBMax);
508  Long64_t egPtr = gROOT->ProcessLine(egMk);
509  if (egPtr == 0) {
510  Error("SetupGen", "Failed to make generator");
511  return false;
512  }
513  fGenerator = reinterpret_cast<AliGenerator*>(egPtr);
514  TString tgt, proj;
515  Int_t tgtA=0, tgtZ=0, projA=0, projZ=0;
516  fGenerator->GetTarget(tgt, tgtA, tgtZ);
517  fGenerator->GetProjectile(proj, projA, projZ);
518  fIsTgtA = !(tgtA == tgtZ && tgtA == 1);
519  fIsProjA = !(projA == projZ && projZ == 1);
520  if (fVerbose)
521  Info("SetupGen", "tgt=%s (%3d,%2d) proj=%s (%3d,%2d) CMS=%fGeV",
522  tgt.Data(), tgtA, tgtZ, proj.Data(), projA, projZ,
523  fGenerator->GetEnergyCMS());
524  Info("SetupGen", "Generator: %s", fGenerator->GetTitle());
525  if (fFileName.IsNull()) FileName();
526 
527  return true;
528  }
534  virtual Bool_t SetupRun()
535  {
536  // --- gAlice (bare ROOT) ----------------------------------------
537  if (!gAlice)
538  new AliRun("gAlice", "The ALICE Off-line framework");
539 
540  Long64_t nev = (fNEvents <= 0 ? 0xFFFFFFFF : fNEvents);
541  Printf("=== Run ================================\n"
542  " Number of events: %lld\n"
543  "========================================", nev);
544  TObject* ord = (fInput ? fInput->FindObject("PROOF_Ordinal") : 0);
545  UShort_t saveMode = 0;
546  TString post = "";
547  TString dir = "";
548  if (ord) {
549  TObject* save = fInput->FindObject("PROOF_SaveGALICE");
550  if (save && fVerbose) {
551  Info("SetupRun", "Got save option:");
552  save->Print();
553  }
554  TString optSave(save ? save->GetTitle() : "split");
555  optSave.ToLower();
556  if (optSave.EqualTo("none")) saveMode = 0;
557  else if (optSave.EqualTo("merge")) saveMode = 1;
558  else if (optSave.EqualTo("split")) saveMode = 2;
559  if (fProofFile && saveMode > 0)
560  dir = fProofFile->GetDir(true);
561  if (saveMode > 1)
562  post = Form("_%s", ord->GetTitle());
563  }
564  TString galiceName(Form("%sgalice.root",dir.Data()));
565  TString kineName(Form("%sKinematics.root",dir.Data()));
566 
567  // --- Run-loader, stack, etc -----------------------------------
568  // Info("SetupRun", "Set-up run Loader");
569  fRunLoader = AliRunLoader::Open(galiceName, "FASTRUN", "RECREATE");
570  fRunLoader->SetKineFileName(kineName);
571  fRunLoader->SetCompressionLevel(2);
572  fRunLoader->SetNumberOfEventsPerFile(nev);
573  fRunLoader->LoadKinematics("RECREATE");
574  fRunLoader->MakeTree("E");
575  gAlice->SetRunLoader(fRunLoader);
576  fRunLoader->MakeStack();
577  fStack = fRunLoader->Stack();
578  fHeader = fRunLoader->GetHeader();
579 
580  // --- Initialize generator --------------------------------------
581  // Info("SetupRun", "Initializing generator");
582  fGenerator->Init();
583  fGenerator->SetStack(fStack);
584 
585  if (saveMode < 1) {
586  if (ord)
587  Info("SetupRun", "Not saving galice.root and Kinematics.root");
588  return true;
589  }
590 
591  TString aliceOut = Form("galice%s.root", post.Data());
592  fAliceFile = new TProofOutputFile(aliceOut, "M");
593 
594  TString kineOut = Form("Kinematics%s.root", post.Data());
595  fKineFile = new TProofOutputFile(kineOut, "M");
596 
597  return true;
598  }
603  Bool_t ReadGRPLine()
604  {
605  std::ifstream* pin = new std::ifstream("grp.dat");
606  if (!pin) {
607  Warning("ReadGRPLine", "Failed to open \"grp.dat\"");
608  return false;
609  }
610  std::istream& in = *pin;
611  TString line;
612  TString env;
613  do {
614  line.ReadLine(in);
615  if (line.IsNull()) continue;
616  if (line.BeginsWith("#")) continue;
617  env = line;
618  break;
619  } while (!in.eof());
620  pin->close();
621 
622  if (env.IsNull()) {
623  Warning("ReadGRPLine", "Got no line from \"grp.dat\"");
624  return false;
625  }
626 
627  fGRP = new TNamed("GRP",env.Data());
628  if (fVerbose) Info("ReadGRPLine", "Read \"%s\"", env.Data());
629  return true;
630  }
635  void OverrideGRP()
636  {
637  Long_t ret = gROOT->ProcessLine("grp");
638  if (ret == 0) {
639  Warning("OverrideGRP", "GRP not set yet, cannot override");
640  return;
641  }
642  if (!fOverrides) {
643  if (fVerbose) Info("OverrideGRP", "No overrides defined");
644  return;
645  }
646  TIter next(fOverrides);
647  TObject* o = 0;
648  while ((o = next())) {
649  if (fVerbose)
650  Info("OverrideGRP", "Overriding GRP setting %s with %s",
651  o->GetName(), o->GetTitle());
652  gROOT->ProcessLine(Form("grp->%s = %s;",
653  o->GetName(), o->GetTitle()));
654  }
655  Info("OverrideGRP", "After overriding:");
656  gROOT->ProcessLine("grp->Print()");
657  }
664  void AddOverride(const TString& field, const TString& value)
665  {
666  if (!fOverrides) {
667  fOverrides = new TList;
668  fOverrides->SetName("overrides");
669  }
670  fOverrides->Add(new TNamed(field, value));
671  }
676  virtual void Init(TTree*)
677  {
678  }
683  virtual void Begin(TTree*)
684  {
685  // Make a monitor
686  // Info("Begin", "gProof=%p Nomonitor=%p",
687  // gProof, (gProof ? gProof->GetParameter("NOMONITOR") : 0));
688  if (fVerbose) Info("Begin", "Called for FastSim");
689 
690  if (fMonitor > 0 && !gROOT->IsBatch()) {
691  FastSimMonitor* m = new FastSimMonitor;
692  m->Connect(fMonitor);
693  }
694  gROOT->Macro(Form("GRP.C(%d)", fRunNo));
695  if (ReadGRPLine()) {
696  if(gProof) {
697  gProof->AddInput(fGRP);
698  if (fOverrides) gProof->AddInput(fOverrides);
699  }
700  }
701  if (fVerbose) Info("Begin", "Perhaps override");
702  OverrideGRP();
703  if (fVerbose) Info("Begin", "Defining centrality estimators");
704 
705  fBEstimator = new BCentEstimator;
706  fCentEstimators = new TList;
707  // fCentEstimators->Add(new V0CentEstimator(-1)); //V0C
708  // fCentEstimators->Add(new V0CentEstimator( 0)); //V0M
709  // fCentEstimators->Add(new V0CentEstimator(+1)); //V0A
710  fCentEstimators->Add(fBEstimator);
711  fCentEstimators->Add(new V0CentEstimator(-1,true)); //V0CP
712  fCentEstimators->Add(new V0CentEstimator( 0,true)); //V0MP
713  fCentEstimators->Add(new V0CentEstimator(+1,true)); //V0AP
714  // fCentEstimators->Add(new ZNCentEstimator(-1,zN|zE)); //ZNCE
715  // fCentEstimators->Add(new ZNCentEstimator(+1,zN|zE)); //ZNAE
716  fCentEstimators->Add(new ZNCentEstimator(-1,true,false)); //ZNCE
717  fCentEstimators->Add(new ZNCentEstimator(+1,true,false)); //ZNAE
718  fCentEstimators->Add(new ZNCentEstimator(-1,true,false,true)); //ZNCEP
719  fCentEstimators->Add(new ZNCentEstimator(+1,true,false,true)); //ZNAE
720  fCentEstimators->Add(new ZNCentEstimator(-1,true,true)); //ZNCS
721  fCentEstimators->Add(new ZNCentEstimator(+1,true,true)); //ZNAS
722  // fCentEstimators->Add(new ZNCentEstimator(-1,zE|zP)); //ZPCE
723  // fCentEstimators->Add(new ZNCentEstimator(+1,zE|zP)); //ZPAE
724  // fCentEstimators->Add(new ZNCentEstimator(-1,zS)); //ZPCS
725  // fCentEstimators->Add(new ZNCentEstimator(+1,zS)); //ZPAS
726  // fCentEstimators->Add(new RefMultEstimator(0.8));
727  // fCentEstimators->Add(new RefMultEstimator(0.5));
728 #if 0
729  TIter next(fCentEstimators);
730  FastCentEstimator* estimator = 0;
731  while ((estimator = static_cast<FastCentEstimator*>(next()))) {
732  estimator->Print("nah");
733  }
734 #endif
735  if (fVerbose) Info("Begin", "End of begin");
736  }
741  void SlaveBegin(TTree*)
742  {
743  if (fVerbose) Info("SlavesBegin", "Called for FastSim");
744  SetupSeed();
745  SetupGRP();
746  SetupGen();
747  SetupOutput();
748  SetupRun();
749  }
750  /* Reset internal caches etc.
751  *
752  * @param iEv Event number
753  *
754  * @return true on success
755  */
756  virtual Bool_t PreEvent(Long64_t iEv)
757  {
758  // --- Reset header ----------------------------------------------
759  fShortHead.Reset(fRunNo, iEv);
760  fParticles->Clear();
761  // --- Reset header, etc. ---------------------------------------
762  fHeader->Reset(fRunNo, iEv);
763  fRunLoader->SetEventNumber(iEv);
764  fStack->Reset();
765  fRunLoader->MakeTree("K");
766 
767  TIter next(fCentEstimators);
768  FastCentEstimator* estimator = 0;
769  while ((estimator = static_cast<FastCentEstimator*>(next())))
770  estimator->PreEvent();
771 
772  return true;
773  }
774 
780  virtual Bool_t ProcessHeader()
781  {
782  // --- Copy to short header --------------------------------------
783  fShortHead.fRunNo = fHeader->GetRun();
784  fShortHead.fEventId = fHeader->GetEvent();
785  TArrayF ip;
786  fHeader->GenEventHeader()->PrimaryVertex(ip);
787  fShortHead.fIpX = ip[0];
788  fShortHead.fIpY = ip[1];
789  fShortHead.fIpZ = ip[2];
790 
791  // --- Check header type -----------------------------------------
792  AliGenEventHeader* genHeader = fHeader->GenEventHeader();
793  AliCollisionGeometry* geometry =
794  dynamic_cast<AliCollisionGeometry*>(genHeader);
795  AliGenPythiaEventHeader* pythia =
796  dynamic_cast<AliGenPythiaEventHeader*>(genHeader);
797  AliGenDPMjetEventHeader* dpm =
798  dynamic_cast<AliGenDPMjetEventHeader*>(genHeader);
799  AliGenGeVSimEventHeader* gev =
800  dynamic_cast<AliGenGeVSimEventHeader*>(genHeader);
801  AliGenHerwigEventHeader* herwig =
802  dynamic_cast<AliGenHerwigEventHeader*>(genHeader);
803  AliGenCocktailEventHeader* cocktail =
804  dynamic_cast<AliGenCocktailEventHeader*>(genHeader);
805  AliGenHijingEventHeader* hijing =
806  dynamic_cast<AliGenHijingEventHeader*>(genHeader);
807  if (cocktail) {
808  TList* headers = cocktail->GetHeaders();
809  if (!headers) Warning("", "No headers in cocktail!");
810  TIter next(headers);
811  AliGenEventHeader* header = 0;
812  AliCollisionGeometry* geom = 0;
813  while ((header = static_cast<AliGenEventHeader*>(next()))) {
814  AliCollisionGeometry* g = dynamic_cast<AliCollisionGeometry*>(header);
815  if (g) geom = g;
816  hijing = dynamic_cast<AliGenHijingEventHeader*>(header);
817  }
818  if (geom) geometry = geom;
819  }
820  if (hijing) fShortHead.fEG = FastShortHeader::kHijing;
821  if (dpm) fShortHead.fEG = FastShortHeader::kDPMJet;
822  if (pythia) fShortHead.fEG = FastShortHeader::kPythia;
823 
824  if (geometry) {
825  fShortHead.fB = geometry->ImpactParameter();
826  fShortHead.fNtgt = geometry->TargetParticipants();
827  fShortHead.fNproj = geometry->ProjectileParticipants();
828  fShortHead.fNbin = geometry->NN();
829  fShortHead.fPhiR = geometry->ReactionPlaneAngle();
830  fShortHead.fNSpecNproj = geometry->ProjSpectatorsn();
831  fShortHead.fNSpecPproj = geometry->ProjSpectatorsp();
832  fShortHead.fNSpecNtgt = geometry->TargSpectatorsn();
833  fShortHead.fNSpecPtgt = geometry->TargSpectatorsp();
834  }
835  // --- Determine diffraction flags -------------------------------
836  Bool_t sd = false;
837  Bool_t dd = false;
838  if (pythia) {
839  Int_t type = pythia->ProcessType();
840  if (type < 100) { // pythia6
841  switch (type) {
842  case 92: case 93: sd = true; break;
843  case 94: dd = true; break;
844  }
845  }
846  else {
847  switch (type) { // Pythia8
848  case 103: case 104: sd = true; break;
849  case 105: dd = true; break;
850  }
851  }
852  fShortHead.fB = pythia->GetImpactParameter();
853  fShortHead.fNtgt = 1;
854  fShortHead.fNproj = 1;
855  fShortHead.fNbin = 1;
856  }
857  if (dpm) {
858  // Try to figure out number of spectators in either side
859  UInt_t nSpecNproj = 0;
860  UInt_t nSpecNtgt = 0;
861  UInt_t nSpecPproj = 0;
862  UInt_t nSpecPtgt = 0;
863  Int_t nPart = fStack->GetNprimary();
864  for (Int_t iPart = 0; iPart < nPart; iPart++) {
865  TParticle* particle = fStack->Particle(iPart);
866  Int_t ks = particle->GetStatusCode();
867  Int_t side = 0;
868  if (ks == 13) side = -1;
869  if (ks == 14) side = +1;
870  if (side == 0) continue;
871  Int_t kf = particle->GetPdgCode();
872  if (kf == kNeutron) {
873  if (side < 0) nSpecNproj++;
874  else nSpecNtgt++;
875  }
876  else if (kf == kProton) {
877  if (side < 0) nSpecPproj++;
878  else nSpecPtgt++;
879  }
880  }
881  fShortHead.fNSpecNtgt = nSpecNtgt;
882  fShortHead.fNSpecNproj = nSpecNproj;
883  fShortHead.fNSpecPtgt = nSpecPtgt;
884  fShortHead.fNSpecPproj = nSpecPproj;
885  // fShortHead.Print();
886 
887  Int_t type = dpm->ProcessType();
888 #ifndef NO_DPMJET_TYPE
889  switch (type) {
890  case 5: case 6: sd = true;
891  case 7: dd = true;
892  }
893 #else
894  static bool first = true;
895 
896  if (first) {
897  Func_t add = gSystem->DynFindSymbol("*", "dtglcp_");
898  if (!add)
899  Warning("", "Didn't find dtglcp_");
900  else
901  _dtglcp = (DtglcpCommon*)add;
902  }
903  // The below - or rather a different implementation with some
904  // errors - was proposed by Cvetan - I don't think it's right
905  // though. See also
906  //
907  // https://cern.ch/twiki/pub/ALICE/PAPaperCentrality/normalization.pdf
908  // https://cern.ch/twiki/bin/view/ALICE/PAMCProductionStudies
909  //
910  Int_t nsd1=0, nsd2=0, ndd=0;
911  Int_t npP = dpm->ProjectileParticipants();
912  Int_t npT = dpm->TargetParticipants();
913  // Get the numbeer of single and double diffractive participants
914  dpm->GetNDiffractive(nsd1,nsd2,ndd);
915  // Check if all partipants are single/double diffractive
916  if ((ndd == 0) && ((npP == nsd1) || (npT == nsd2))) sd = true;
917  else if (ndd == (npP + npT)) dd = true;
918  Int_t ncp = dpm->NN();
919  Int_t nct = dpm->NNw();
920  Int_t nwp = dpm->NwN();
921  Int_t nwt = dpm->NwNw();
922  Int_t nwtacc = _dtglcp->nwtacc;
923  Int_t nwtsam = _dtglcp->nwtsam;
924  if (first) {
925  Printf("@ Npp sd1 Npt sd2 dd tpe Ncp Nct Nwp Nwt acc sam");
926  first = false;
927  }
928  Printf("@ %3d %3d %3d %3d %3d %3d %3d %3d %3d %3d %3d %3d",
929  npP, nsd1, npT, nsd2, ndd, type, ncp, nct, nwp, nwt,
930  nwtacc, nwtsam);
931 #endif
932  }
933  if (gev) fShortHead.fPhiR = gev->GetEventPlane();
934  if (herwig) {
935  Int_t type = herwig->ProcessType();
936  switch (type) {
937  case 5: case 6: sd = true; break;
938  }
939  fShortHead.fNtgt = 1;
940  fShortHead.fNproj = 1;
941  fShortHead.fNbin = 1;
942  }
943  fShortHead.fType = (sd ? 0x1 : 0) | (dd ? 0x2 : 0);
944 
945  // --- Check centrality -----------------------------------------
946  Double_t b = fShortHead.fB;
947  Double_t c = fBEstimator->GetCentrality(b);
948  Info("ProcessHeader", "b=%f isProjA=%d isTgtA=%d cms=%f",
949  b, fIsProjA, fIsTgtA, fGenerator ? fGenerator->GetEnergyCMS() : -1);
950  if (c >= 0) fShortHead.fC = c;
951 
952  // --- Check if within vertex cut -------------------------------
953  Bool_t selected = (fShortHead.fIpZ <= fHIpz->GetXaxis()->GetXmax() &&
954  fShortHead.fIpZ >= fHIpz->GetXaxis()->GetXmin());
955 
956  // --- Only update histograms if within IPz cut ------------------
957  if (selected) {
958  fHPhiR->Fill(fShortHead.fPhiR*TMath::RadToDeg());
959  fHB->Fill(fShortHead.fB);
960  fHIpz->Fill(fShortHead.fIpZ);
961  if (dd) fHType->Fill(3);
962  if (sd) fHType->Fill(2);
963  if (!dd && !sd) fHType->Fill(1);
964  fHCent->Fill(c);
965  fHEta ->AddBinContent(0); // Count events
966  fHY ->AddBinContent(0); // Count events
967  fHTime->AddBinContent(0); // Count events
968  // fShortHead.Print();
969  }
970  TIter next(fCentEstimators);
971  FastCentEstimator* estimator = 0;
972  while ((estimator = static_cast<FastCentEstimator*>(next())))
973  estimator->ProcessHeader(fShortHead);
974 
975  return selected;
976  }
984  virtual Bool_t ProcessParticles(Bool_t selected)
985  {
986  Int_t nPart = fStack->GetNprimary();
987  for (Int_t iPart = 0; iPart < nPart; iPart++) {
988  TParticle* particle = fStack->Particle(iPart);
989  TParticlePDG* pdg = particle->GetPDG();
990  Bool_t primary = fStack->IsPhysicalPrimary(iPart);
991  Bool_t weakDecay = fStack->IsSecondaryFromWeakDecay(iPart);
992  Bool_t charged = (pdg && TMath::Abs(pdg->Charge()) > 0);
993  if (primary) particle->SetBit(BIT(14));
994  if (weakDecay) particle->SetBit(BIT(15));
995  if (charged) particle->SetBit(BIT(16));
996 
997  new ((*fParticles)[iPart]) TParticle(*particle);
998 
999  TIter next(fCentEstimators);
1000  FastCentEstimator* estimator = 0;
1001  while ((estimator = static_cast<FastCentEstimator*>(next())))
1002  estimator->Process(particle);
1003 
1004  if (!selected || !charged || !primary) continue;
1005 
1006  Double_t y = particle->Y();
1007  if (y > fHY->GetXaxis()->GetXmin() &&
1008  y < fHY->GetXaxis()->GetXmax())
1009  // Avoid filling under/overflow bins
1010  fHY->Fill(y);
1011 
1012  Double_t pT = particle->Pt();
1013  if (pT < 1e-10) continue;
1014  Double_t pZ = particle->Pz();
1015  Double_t theta = TMath::ATan2(pT, pZ);
1016  Double_t eta = -TMath::Log(TMath::Tan(theta/2));
1017  if (eta > fHEta->GetXaxis()->GetXmin() &&
1018  eta < fHEta->GetXaxis()->GetXmax())
1019  // Avoid filling under/overflow bins
1020  fHEta->Fill(eta);
1021  }
1022  return true;
1023  }
1028  virtual void PostEvent()
1029  {
1030  fHeader->SetNprimary(fStack->GetNprimary());
1031  fHeader->SetNtrack(fStack->GetNtrack());
1032 
1033  TIter next(fCentEstimators);
1034  FastCentEstimator* estimator = 0;
1035  while ((estimator = static_cast<FastCentEstimator*>(next())))
1036  estimator->PostEvent();
1037 
1038  fTree->Fill();
1039 
1040  fStack->FinishEvent();
1041  fHeader->SetStack(fStack);
1042 
1043  fRunLoader->TreeE()->Fill();
1044  fRunLoader->WriteKinematics("OVERWRITE");
1045  }
1046  virtual void Generate() { fGenerator->Generate(); }
1047 
1055  virtual Bool_t Process(Long64_t iEv)
1056  {
1057  // --- The stopwatch ---------------------------------------------
1058  TStopwatch timer;
1059  timer.Start();
1060  PreEvent(iEv);
1061  fHTime->Fill(1, timer.RealTime());
1062 
1063  // --- Generate event --------------------------------------------
1064  timer.Start();
1065  Generate();
1066  fHTime->Fill(2, timer.RealTime());
1067 
1068  // --- Process the header ----------------------------------------
1069  timer.Start();
1070  Bool_t selected = ProcessHeader();
1071  fHTime->Fill(3, timer.RealTime());
1072 
1073  // --- Loop over particles ---------------------------------------
1074  timer.Start();
1075  ProcessParticles(selected);
1076  fHTime->Fill(4, timer.RealTime());
1077 
1078  // --- Do final stuff --------------------------------------------
1079  timer.Start();
1080  PostEvent();
1081  fHTime->Fill(5, timer.RealTime());
1082 
1083  return true;
1084  }
1085  virtual void FinishRun()
1086  {
1087  fGenerator->FinishRun();
1088  fRunLoader->WriteHeader("OVERWRITE");
1089  fGenerator->Write();
1090  fRunLoader->Write();
1091 
1092  }
1097  void SlaveTerminate()
1098  {
1099  FinishRun();
1100  if (fFile) {
1101  if (fProofFile) {
1102  if (fVerbose) fProofFile->Print();
1103  fOutput->Add(fProofFile);
1104  fOutput->Add(new TH1F("filename", fFileName.Data(),1,0,1));
1105  }
1106  // Flush out tree
1107  fFile->cd();
1108  fTree->Write(0, TObject::kOverwrite);
1109  fFile->Close();
1110  fFile->Delete();
1111  fFile = 0;
1112  }
1113  if (fAliceFile) {
1114  TFile* galice = GetGAlice();
1115  if (galice) {
1116  if (fVerbose) fAliceFile->Print();
1117  fAliceFile->AdoptFile(galice);
1118  fAliceFile->SetOutputFileName(fAliceFile->GetName());
1119  fOutput->Add(fAliceFile);
1120  galice->Write();
1121  }
1122  }
1123  if (fKineFile) {
1124  TFile* kine = GetKine();
1125  if (kine) {
1126  if (fVerbose) fKineFile->Print();
1127  fKineFile->AdoptFile(kine);
1128  fKineFile->SetOutputFileName(fKineFile->GetName());
1129  fOutput->Add(fKineFile);
1130  kine->Write();
1131  }
1132  }
1133 
1134  if (fVerbose) {
1135  Info("SlaveTerminate", "Content of output list");
1136  gROOT->IncreaseDirLevel();
1137  fOutput->ls();
1138  gROOT->DecreaseDirLevel();
1139 
1140  gSystem->Exec("echo \"Content of working directory\"");
1141  gSystem->Exec("ls -l1 | sed 's/^/ /'");
1142  }
1143  }
1151  void FlushList(TCollection* c, TDirectory* dir)
1152  {
1153  dir->cd();
1154  TIter next(c);
1155  TObject* o = 0;
1156  while ((o = next())) {
1157  if (o->IsA()->InheritsFrom(TCollection::Class())) {
1158  if (fVerbose) Info("FlushList", "Got collection: %s", c->GetName());
1159  TDirectory* cur = dir->mkdir(o->GetName());
1160  FlushList(static_cast<TCollection*>(o), cur);
1161  dir->cd();
1162  continue;
1163  }
1164  o->Write();
1165  }
1166  dir->cd();
1167  }
1168 
1173  void Terminate()
1174  {
1175  if (gProof) gProof->ClearFeedback();
1176 
1177  if (!fList)
1178  fList = static_cast<TList*>(fOutput->FindObject("histograms"));
1179  if (!fList) {
1180  Error("Terminate", "No output list");
1181  return;
1182  }
1183 
1184  if (!fProofFile) {
1185  TObject* fn = fOutput->FindObject("filename");
1186  if (fn) fFileName = fn->GetTitle();
1187  fProofFile =
1188  static_cast<TProofOutputFile*>(fOutput->FindObject(FileName()));
1189  }
1190  if (fProofFile)
1191  fFile = fProofFile->OpenFile("UPDATE");
1192  if (!fFile)
1193  fFile = TFile::Open(FileName(),"UPDATE");
1194 
1195  TList* estimators = static_cast<TList*>(fList->FindObject("estimators"));
1196  TList* histos = static_cast<TList*>(fList->FindObject("histograms"));
1197  if (!histos) {
1198  Warning("Terminate", "No histogram list found in output");
1199  fList->ls();
1200  }
1201  TIter next(fCentEstimators);
1202  FastCentEstimator* estimator = 0;
1203  while ((estimator = static_cast<FastCentEstimator*>(next())))
1204  estimator->Terminate(estimators);
1205 
1206  fHEta = static_cast<TH1*>(histos->FindObject("dNdeta"));
1207  fHY = static_cast<TH1*>(histos->FindObject("dNdy"));
1208  fHIpz = static_cast<TH1*>(histos->FindObject("ipZ"));
1209  fHType = static_cast<TH1*>(histos->FindObject("type"));
1210  fHCent = static_cast<TH1*>(histos->FindObject("cent"));
1211  fHB = static_cast<TH1*>(histos->FindObject("b"));
1212  fHPhiR = static_cast<TH1*>(histos->FindObject("phiR"));
1213  fHTime = static_cast<TH1*>(histos->FindObject("timing"));
1214 
1215  if (!(fHEta && fHY && fHIpz && fHType && fHB && fHPhiR && fHTime)) {
1216  Warning("Terminate", "Missing histograms (%p,%p,%p,%p,%p,%p,%p)",
1217  fHEta, fHY, fHIpz, fHType, fHB, fHPhiR, fHTime);
1218  return;
1219  }
1220 
1221  Int_t nTotal = fHIpz->GetEntries();
1222  fHEta ->Scale(1./nTotal, "width");
1223  fHY ->Scale(1./nTotal, "width");
1224  fHB ->Scale(1./nTotal, "width");
1225  fHPhiR->Scale(1./nTotal, "width");
1226  fHTime->Scale(1./nTotal, "width");
1227 
1228  if (!fFile){
1229  Warning("Terminate", "No file to write to");
1230  return;
1231  }
1232 
1233  FlushList(fList, fFile); // ->Write();
1234 
1235  fTree = static_cast<TTree*>(fFile->Get("T"));
1236  if (!fTree) Warning("Terminate", "No tree");
1237 
1238  if (fVerbose) fFile->ls();
1239  fFile->Close();
1240 
1241  MoveAliceFiles();
1242  }
1248  TFile* GetGAlice()
1249  {
1250  if (!fRunLoader) return 0;
1251  TString galiceName = fRunLoader->GetFileName();
1252  TFile* file = gROOT->GetFile(galiceName);
1253  if (!file) {
1254  Warning("GetGAlice", "Didn't find galice file \"%s\"", galiceName.Data());
1255  gROOT->GetListOfFiles()->ls();
1256  return 0;
1257  }
1258  return file;
1259  }
1265  TFile* GetKine()
1266  {
1267  if (!fRunLoader) return 0;
1268  TString kineName = "Kinematics.root";
1269  TString galiceName = fRunLoader->GetFileName();
1270  TString dir = gSystem->DirName(galiceName);
1271  if (dir.EqualTo(".")) dir = "";
1272  if (!dir.IsNull() && dir[dir.Length()-1] != '/') dir.Append("/");
1273  kineName.Prepend(dir);
1274 
1275  TFile* file = gROOT->GetFile(kineName);
1276  if (!file) {
1277  Warning("GetKine", "Didn't find kinematics file \"%s\"", kineName.Data());
1278  gROOT->GetListOfFiles()->ls();
1279  return 0;
1280  }
1281  return file;
1282  }
1289  void MoveAliceFiles()
1290  {
1291  if (!fInput) return;
1292 
1293  TObject* save = fInput->FindObject("PROOF_SaveGALICE");
1294  if (!save) return;
1295 
1296  TString sMode = save->GetTitle();
1297  if (!sMode.EqualTo("split", TString::kIgnoreCase)) return;
1298 
1299  TList* lst = new TList;
1300  TSystemDirectory* dir = new TSystemDirectory(".",
1301  gSystem->WorkingDirectory());
1302  TList* files = dir->GetListOfFiles();
1303  TSystemFile* file = 0;
1304  TIter next(files);
1305  while ((file = static_cast<TSystemFile*>(next()))) {
1306  if (file->IsDirectory()) continue;
1307  TString fn(file->GetName());
1308  if (!fn.BeginsWith("galice") && !fn.BeginsWith("Kinematics"))
1309  continue;
1310 
1311  TPRegexp regex("(.*)_([^_]+)\\.root");
1312  TObjArray* matches = regex.MatchS(fn);
1313  if (matches->GetEntriesFast() < 3) {
1314  delete matches;
1315  continue;
1316  }
1317  TString ord = matches->At(2)->GetName();
1318  TString bse = matches->At(1)->GetName();
1319 
1320  if (gSystem->AccessPathName(ord,kFileExists))
1321  gSystem->MakeDirectory(ord);
1322 
1323  if (fVerbose)
1324  Info("MoveAliceFiles", "Moving %s to %s/%s.root",
1325  fn.Data(), ord.Data(), bse.Data());
1326  file->Move(Form("%s/%s.root", ord.Data(), bse.Data()));
1327 
1328  if (!bse.EqualTo("galice")) continue;
1329  TObjString* url = new TObjString(Form("file://%s/%s/%s.root?#TE",
1330  file->GetTitle(),
1331  ord.Data(),
1332  bse.Data()));
1333  if (fVerbose)
1334  Info("MoveAliceFiles", "Adding \"%s\" to file list",
1335  url->GetName());
1336  lst->Add(url);
1337  }
1338  if (lst->GetEntries() <= 0) return;
1339  if (fVerbose) lst->ls();
1340 
1341  TFile* out = TFile::Open("index.root","RECREATE");
1342  lst->Write("TE",TObject::kSingleKey);
1343  out->Write();
1344  out->Close();
1345 
1346  }
1352  Int_t Version() const { return 1; }
1353 
1358  TString fEGName; // Name of event generator
1359  Int_t fRunNo; // Run to simulate
1360  Double_t fBMin; // Least impact parameter
1361  Double_t fBMax; // Largest impact parameter
1362  TObject* fGRP;
1363  TList* fOverrides;
1364  Long64_t fNEvents; // Number of requested events
1365  Bool_t fIsTgtA;
1366  Bool_t fIsProjA;
1367  /* @} */
1372  AliGenerator* fGenerator;
1373  AliRunLoader* fRunLoader;
1374  AliStack* fStack;
1375  AliHeader* fHeader;
1376  /* @} */
1381  TTree* fTree;
1382  TClonesArray* fParticles;
1383 
1387  TList* fList;
1388  TH1* fHEta;
1389  TH1* fHY;
1390  TH1* fHIpz;
1391  TH1* fHType;
1392  TH1* fHCent;
1393  TH1* fHB;
1394  TH1* fHPhiR;
1395  TH1* fHTime;
1396  /* @} */
1401  BCentEstimator* fBEstimator; // Always present
1402  TList* fCentEstimators; // Centrality estimators
1403  /* @} */
1408  TProofOutputFile* fProofFile;
1409  TProofOutputFile* fAliceFile;
1410  TProofOutputFile* fKineFile;
1411  TFile* fFile;
1412  mutable TString fFileName;
1413  /* @} */
1414  Bool_t fVerbose; // Verbosity
1415  Int_t fMonitor;
1416 
1417 #ifndef __CINT__
1418  FastShortHeader fShortHead;
1419 #endif
1420 
1434  static Bool_t LocalRun(Long64_t nev,
1435  UInt_t run,
1436  const TString& gen,
1437  Double_t bMin,
1438  Double_t bMax,
1439  Int_t monitor,
1440  Bool_t verbose,
1441  const TString& overrides="")
1442 
1443  {
1444  FastSim* sim = new FastSim(gen,run,bMin,bMax,nev,monitor);
1445  SetOverrides(sim, overrides);
1446  sim->fVerbose = verbose;
1447  sim->Begin(0);
1448  sim->SlaveBegin(0);
1449 
1450  for (Long64_t i=0; i <nev; i++) {
1451  Printf("=== Event # %6lld/%6lld ==========================",
1452  i+1, nev);
1453  sim->Process(i);
1454  }
1455  sim->SlaveTerminate();
1456  sim->Terminate();
1457 
1458  return true;
1459  }
1463  static void ProofLoadLibs()
1464  {
1465  if (!gProof) return;
1466 
1467  // Remember to copy changes to RunFast.C
1468  TList clsLib;
1469  clsLib.Add(new TNamed("TVirtualMC", "libVMC"));
1470  clsLib.Add(new TNamed("TLorentzVector", "libPhysics"));
1471  clsLib.Add(new TNamed("TLinearFitter", "libMinuit"));
1472  clsLib.Add(new TNamed("TTree", "libTree"));
1473  clsLib.Add(new TNamed("TProof", "libProof"));
1474  clsLib.Add(new TNamed("TGFrame", "libGui"));
1475  clsLib.Add(new TNamed("TSAXParser", "libXMLParser"));
1476  clsLib.Add(new TNamed("AliVEvent", "libSTEERBase"));
1477  clsLib.Add(new TNamed("AliESDEvent", "libESD"));
1478  clsLib.Add(new TNamed("AliAODEvent", "libAOD"));
1479  clsLib.Add(new TNamed("AliAnalysisManager", "libANALYSIS"));
1480  clsLib.Add(new TNamed("AliCDBManager", "libCDB"));
1481  clsLib.Add(new TNamed("AliRawVEvent", "libRAWDatabase"));
1482  clsLib.Add(new TNamed("AliHit", "libSTEER"));
1483  clsLib.Add(new TNamed("AliGenMC", "libEVGEN"));
1484  clsLib.Add(new TNamed("AliFastEvent", "libFASTSIM"));
1485 
1486  TIter next(&clsLib);
1487  TObject* obj = 0;
1488  while ((obj = next())) {
1489  gProof->Exec(Form("gROOT->LoadClass(\"%s\",\"%s\");",
1490  obj->GetName(), obj->GetTitle()));
1491  }
1492  }
1510  static Bool_t ProofRun(const TUrl& url,
1511  Long64_t nev,
1512  UInt_t run,
1513  const TString& gen,
1514  Double_t bMin,
1515  Double_t bMax,
1516  Int_t monitor=-1,
1517  Bool_t verbose=false,
1518  const TString& overrides="",
1519  const TString& save="none",
1520  const char* opt="")
1521  {
1522  TProof::Reset(url.GetUrl());
1523  TProof::Open(url.GetUrl());
1524  gProof->ClearCache();
1525 
1526  TString phy = gSystem->ExpandPathName("$(ALICE_PHYSICS)");
1527  TString ali = gSystem->ExpandPathName("$(ALICE_ROOT)");
1528  // TString fwd = gSystem->ExpandPathName("$ANA_SRC");
1529  TString fwd = phy + "/PWGLF/FORWARD/analysis2";
1530 
1531  gProof->AddIncludePath(Form("%s/include", ali.Data()));
1532  gProof->AddIncludePath(Form("%s/include", phy.Data()));
1533  ProofLoadLibs();
1534  gProof->Load(Form("%s/sim/GRP.C",fwd.Data()), true);
1535  gProof->Load(Form("%s/sim/BaseConfig.C",fwd.Data()), true);
1536  gProof->Load(Form("%s/sim/EGConfig.C",fwd.Data()), true);
1537 
1538  // gROOT->ProcessLine("gProof->SetLogLevel(5);");
1539  gProof->Load(Form("%s/sim/FastShortHeader.C", fwd.Data()));
1540  gProof->Load(Form("%s/sim/FastCentEstimators.C+%s",fwd.Data(),opt));
1541  gProof->Load(Form("%s/sim/FastMonitor.C+%s",fwd.Data(),opt));
1542  gProof->Load(Form("%s/sim/FastSim.C+%s", fwd.Data(), opt),true);
1543  gProof->SetParameter("PROOF_SaveGALICE", save);
1544 
1545  FastSim* sim = new FastSim(gen,run,bMin,bMax,nev,monitor);
1546  SetOverrides(sim, overrides);
1547  sim->fVerbose = verbose;
1548  gProof->Process(sim, nev, "");
1549 
1550  return true; // status >= 0;
1551  }
1562  static Bool_t Str2KeyVal(const TString& in,
1563  TString& key,
1564  TString& val,
1565  const char sep='=')
1566  {
1567  Int_t idx = in.Index(sep);
1568  if (idx == kNPOS) return false;
1569 
1570  key = in(0,idx);
1571  val = in(idx+1, in.Length()-idx-1);
1572  return true;
1573  }
1574  static void SetOverrides(FastSim* sim, const TString& override)
1575  {
1576  if (override.IsNull()) return;
1577 
1578  const char* valid[] = { "beamEnergy", // UInt_t [GeV]
1579  "energy", // UInt_t [GeV]
1580  "period", // String
1581  "run", // UInt_t
1582  "beam1.a", // UInt_t
1583  "beam1.z", // UInt_t
1584  "beam2.a", // UInt_t
1585  "beam2.z", // UInt_t
1586  0 };
1587  TObjArray* tokens = override.Tokenize(",");
1588  TObjString* token = 0;
1589  TIter next(tokens);
1590  while ((token = static_cast<TObjString*>(next()))) {
1591  TString& str = token->String();
1592  if (str.IsNull()) continue;
1593 
1594  TString key, val;
1595  if (!Str2KeyVal(str,key,val, ':')) {
1596  Printf("Warning: FastSim::Run: incomplete override '%s'",str.Data());
1597  continue;
1598  }
1599  const char** pvalid = valid;
1600  while (*pvalid) {
1601  if (key.EqualTo(*pvalid, TString::kIgnoreCase)) {
1602  break;
1603  }
1604  pvalid++;
1605  }
1606  if (!*pvalid) {
1607  Printf("Warning: FastSim::Run: Invalid override '%s'", key.Data());
1608  continue;
1609  }
1610  // Special case for a string
1611  if (key.EqualTo("period",TString::kIgnoreCase))
1612  val = Form("\"%s\"", val.Data());
1613  sim->AddOverride(*pvalid, val);
1614  }
1615  // delete tokens;
1616  }
1648  static Bool_t Run(const char* url, const char* opt="")
1649  {
1650  Printf("Will run fast simulation with:\n\n\t%s\n\n",url);
1651  Long64_t nev = 10000;
1652  UInt_t run = 0;
1653  TString eg = "default";
1654  TString override= "";
1655  TString save = "none";
1656  Double_t bMin = 0;
1657  Double_t bMax = 20;
1658  Int_t monitor = -1;
1659  Bool_t verbose = false;
1660  TUrl u(url);
1661  TString out;
1662  TObjArray* opts = TString(u.GetOptions()).Tokenize("&");
1663  TObjString* token = 0;
1664  TIter nextToken(opts);
1665  while ((token = static_cast<TObjString*>(nextToken()))) {
1666  TString& str = token->String();
1667  if (str.IsNull()) continue;
1668 
1669  if (str.EqualTo("verbose")) { verbose = true; str = ""; }
1670 
1671  if (str.IsNull()) continue;
1672 
1673  TString key, val;
1674  if (!Str2KeyVal(str,key,val)) {
1675  if (!out.IsNull()) out.Append("&");
1676  out.Append(str);
1677  continue;
1678  }
1679 
1680  if (key.EqualTo("events")) nev = val.Atoll();
1681  else if (key.EqualTo("run")) run = val.Atoi();
1682  else if (key.EqualTo("eg")) eg = val;
1683  else if (key.EqualTo("override")) override = val;
1684  else if (key.EqualTo("save")) save = val;
1685  else if (key.EqualTo("monitor")) monitor = val.Atoi();
1686  else if (key.EqualTo("b")) {
1687  TString min, max;
1688  if (Str2KeyVal(val, min, max, '-')) {
1689  bMin = min.Atof();
1690  bMax = max.Atof();
1691  }
1692  }
1693  else {
1694  if (!out.IsNull()) out.Append("&");
1695  out.Append(str);
1696  }
1697  }
1698  opts->Delete();
1699  u.SetOptions(out);
1700  if (!u.IsValid()) {
1701  Printf("Error: FastSim::Run: URL %s is invalid", u.GetUrl());
1702  return false;
1703  }
1704 
1705  Bool_t isLocal = TString(u.GetProtocol()).EqualTo("local");
1706 
1707  Printf("Run %s for %lld events anchored at %d\n"
1708  " Impact paramter range: %5.1f-%5.1f fm\n"
1709  " Monitor frequency: %d sec\n"
1710  " Execution url: %s",
1711  eg.Data(), nev, run, bMin, bMax, monitor, u.GetUrl());
1712 
1713 
1714  TStopwatch timer;
1715  timer.Start();
1716 
1717  Bool_t ret = false;
1718  if (isLocal)
1719  ret = LocalRun(nev, run, eg, bMin, bMax, monitor, verbose, override);
1720  else
1721  ret = ProofRun(u, nev, run, eg, bMin, bMax,
1722  monitor, verbose, override, save, opt);
1723  timer.Print();
1724 
1725  return ret;
1726  }
1727 
1728  ClassDef(FastSim,3);
1729 };
1730 
1731 
1732 struct EPosSim : public FastSim
1733 {
1734  EPosSim(UInt_t run=0, Int_t monitor=0)
1735  : FastSim("epos", run, 0, 20, 100000, monitor),
1736  fInTree(0),
1737  fInNTot(0),
1738  fInB(0),
1739  fInPDG(0),
1740  fInStatus(0),
1741  fInPx(0),
1742  fInPy(0),
1743  fInPz(0),
1744  fInE(0),
1745  fInM(0),
1746  fInNcollH(0),
1747  fInNpartP(0),
1748  fInNpartT(0),
1749  fInNcoll(0),
1750  fInNSpcPN(0),
1751  fInNSpcTN(0),
1752  fInNSpcPP(0),
1753  fInNSpcTP(0),
1754  fInPhiR(0)
1755  {
1756  }
1757  Bool_t SetupBranches()
1758  {
1759  fInNTot = fInTree->GetLeaf("nPart");
1760  fInB = fInTree->GetLeaf("ImpactParameter");
1761  fInPDG = fInTree->GetLeaf("pdgid");
1762  fInStatus = fInTree->GetLeaf("status");
1763  fInPx = fInTree->GetLeaf("px");
1764  fInPy = fInTree->GetLeaf("py");
1765  fInPz = fInTree->GetLeaf("pz");
1766  fInE = fInTree->GetLeaf("E");
1767  fInM = fInTree->GetLeaf("m");
1768  // These are probably EPOS-LHC specific
1769  fInNcollH = fInTree->GetLeaf("Ncoll_hard");
1770  fInNpartP = fInTree->GetLeaf("Npart_proj");
1771  fInNpartT = fInTree->GetLeaf("Npart_targ");
1772  fInNcoll = fInTree->GetLeaf("Ncoll");
1773  fInNSpcPN = fInTree->GetLeaf("Nspec_proj_neut");
1774  fInNSpcTN = fInTree->GetLeaf("Nspec_targ_neut");
1775  fInNSpcPP = fInTree->GetLeaf("Nspec_proj_prot");
1776  fInNSpcTP = fInTree->GetLeaf("Nspec_targ_prot");
1777  fInPhiR = fInTree->GetLeaf("phiR");
1778 
1779 
1780 
1781  return (fInNTot && fInB && fInPDG && fInPx && fInPy && fInPz && fInE);
1782  }
1783  void Begin(TTree* tree)
1784  {
1785  FastSim::Begin(tree);
1786  TIter next(fCentEstimators);
1787  FastCentEstimator* estimator = 0;
1788  while ((estimator = static_cast<FastCentEstimator*>(next()))) {
1789  if (!estimator->IsA()->InheritsFrom(V0CentEstimator::Class()))
1790  continue;
1791  V0CentEstimator* v = static_cast<V0CentEstimator*>(estimator);
1792  v->Flip(false); // Flip detector acceptance of V0A/C
1793  // v->Print("nah");
1794  }
1795  }
1796  void Init(TTree* tree)
1797  {
1798  Info("Init", "Initializing with tree %p (%s)",
1799  tree, (tree ? tree->ClassName() : ""));
1800  if (!tree) return;
1801 
1802  TFile* file = tree->GetCurrentFile();
1803  Info("Init", "Current file: (%p) %s", file,
1804  (file ? file->GetName() : ""));
1805 
1806  fInTree = tree;
1807  if (!SetupBranches())
1808  Fatal("Init", "Failed to set-up branches");
1809  // if (!SetupEstimator())
1810  // Fatal("Init", "Failed to set-up estimator");
1811  }
1817  Bool_t Notify()
1818  {
1819  if (!fInTree) {
1820  Warning("Notify", "No tree set yet!");
1821  return false;
1822  }
1823  TFile* file = fInTree->GetCurrentFile();
1824  Info("Notify", "processing file: (%p) %s", file,
1825  (file ? file->GetName() : ""));
1826  if (!file) return true;
1827  if (!SetupBranches()) {
1828  Warning("Notify", "Failed to set-up branches");
1829  return false;
1830  }
1831  return true;
1832  }
1833  const char* GetEGTitle() const { return "EPOS-LHC"; }
1834  void SetupSeed() {}
1835  Bool_t SetupGen()
1836  {
1837  fIsTgtA = gROOT->ProcessLine("grp->beam1.IsA()");
1838  fIsProjA = gROOT->ProcessLine("grp->beam2.IsA()");
1839  return true;
1840 
1841  }
1842  Bool_t SetupRun() { return true; }
1843  Bool_t PreEvent(Long64_t iEv)
1844  {
1845  Int_t read = fInTree->GetTree()->GetEntry(iEv);
1846  if (read <= 0) return false;
1847 
1848  // --- Reset header ----------------------------------------------
1849  fShortHead.Reset(fRunNo, iEv);
1850  fParticles->Clear();
1851  // Reset input
1852 
1853  TIter next(fCentEstimators);
1854  FastCentEstimator* estimator = 0;
1855  while ((estimator = static_cast<FastCentEstimator*>(next())))
1856  estimator->PreEvent();
1857 
1858  return true;
1859  }
1860  Bool_t ProcessHeader()
1861  {
1862  fShortHead.fIpX = 0;
1863  fShortHead.fIpY = 0;
1864  fShortHead.fIpZ = 0;
1865  fShortHead.fB = fInB->GetValue();
1866  fShortHead.fNtgt = (fInNpartT ? fInNpartT->GetValue() : 0);
1867  fShortHead.fNproj = (fInNpartP ? fInNpartP->GetValue() : 0);
1868  fShortHead.fNbin = (fInNcoll ? fInNcoll ->GetValue() : 0);
1869  fShortHead.fPhiR = (fInPhiR ? fInPhiR ->GetValue() : 0);
1870  fShortHead.fNSpecNproj = (fInNSpcPN ? fInNSpcPN->GetValue() : 0);
1871  fShortHead.fNSpecNtgt = (fInNSpcTN ? fInNSpcTN->GetValue() : 0);
1872  fShortHead.fNSpecPproj = (fInNSpcPP ? fInNSpcPP->GetValue() : 0);
1873  fShortHead.fNSpecPtgt = (fInNSpcTP ? fInNSpcTP->GetValue() : 0);
1874  fShortHead.fEG = FastShortHeader::kEPOS;
1875 
1876  Double_t c = fBEstimator->GetCentrality(fShortHead.fB);
1877  if (c >= 0) fShortHead.fC = c;
1878 
1879  // --- Check if within vertex cut -------------------------------
1880  Bool_t selected = (fShortHead.fIpZ <= fHIpz->GetXaxis()->GetXmax() &&
1881  fShortHead.fIpZ >= fHIpz->GetXaxis()->GetXmin());
1882 
1883  // --- Only update histograms if within IPz cut ------------------
1884  if (selected) {
1885  fHPhiR->Fill(fShortHead.fPhiR*TMath::RadToDeg());
1886  fHB->Fill(fShortHead.fB);
1887  fHIpz->Fill(fShortHead.fIpZ);
1888  fHCent->Fill(c);
1889  // fShortHead.Print();
1890  }
1891  TIter next(fCentEstimators);
1892  FastCentEstimator* estimator = 0;
1893  while ((estimator = static_cast<FastCentEstimator*>(next())))
1894  estimator->ProcessHeader(fShortHead);
1895  return selected;
1896  }
1897  virtual Bool_t ProcessParticles(Bool_t selected)
1898  {
1899  Int_t nTot = fInNTot->GetValue();
1900  for (Int_t iPart = 0; iPart < nTot; iPart++) {
1901  Int_t status = fInStatus->GetValue(iPart);
1902  Int_t pdg = fInPDG->GetValue(iPart);
1903  Double_t px = fInPx->GetValue(iPart);
1904  Double_t py = fInPy->GetValue(iPart);
1905  Double_t pz = fInPz->GetValue(iPart);
1906  // Double_t pz = -fInPz->GetValue(iPart); // flip sign on pZ
1907  Double_t e = fInE->GetValue(iPart);
1908  // Double_t m = fInM->GetValue(iPart);
1909 
1910  TParticle* particle =
1911  new ((*fParticles)[iPart]) TParticle(pdg, status,-1,-1,-1,-1,
1912  px, py, pz, e, 0, 0, 0, 0);
1913  TParticlePDG* pdgP = particle->GetPDG();
1914  Bool_t primary = status == 1;
1915  Bool_t weakDecay = false;
1916  Bool_t charged = (pdgP && TMath::Abs(pdgP->Charge()) > 0);
1917  if (primary) particle->SetBit(BIT(14));
1918  if (weakDecay) particle->SetBit(BIT(15));
1919  if (charged) particle->SetBit(BIT(16));
1920 
1921  TIter next(fCentEstimators);
1922  FastCentEstimator* estimator = 0;
1923  while ((estimator = static_cast<FastCentEstimator*>(next())))
1924  estimator->Process(particle);
1925 
1926  if (!selected || !charged || !primary) continue;
1927  fHY ->Fill(particle->Y());
1928  Double_t pT = particle->Pt();
1929  if (pT < 1e-10) continue;
1930  Double_t pZ = particle->Pz();
1931  Double_t theta = TMath::ATan2(pT, pZ);
1932  Double_t eta = -TMath::Log(TMath::Tan(theta/2));
1933  fHEta->Fill(eta);
1934  }
1935  return true;
1936  }
1937  void PostEvent()
1938  {
1939  fTree->Fill();
1940 
1941  TIter next(fCentEstimators);
1942  FastCentEstimator* estimator = 0;
1943  while ((estimator = static_cast<FastCentEstimator*>(next())))
1944  estimator->PostEvent();
1945  }
1946  void Generate() {}
1947  void FinishRun() {}
1948  TTree* fInTree;
1949  TLeaf* fInNTot;
1950  TLeaf* fInB;
1951  TLeaf* fInPDG;
1952  TLeaf* fInStatus;
1953  TLeaf* fInPx;
1954  TLeaf* fInPy;
1955  TLeaf* fInPz;
1956  TLeaf* fInE;
1957  TLeaf* fInM;
1958  TLeaf* fInNcollH;
1959  TLeaf* fInNpartP;
1960  TLeaf* fInNpartT;
1961  TLeaf* fInNcoll;
1962  TLeaf* fInNSpcPN;
1963  TLeaf* fInNSpcTN;
1964  TLeaf* fInNSpcPP;
1965  TLeaf* fInNSpcTP;
1966  TLeaf* fInPhiR;
1967 
1978  static Bool_t LocalRun(Long64_t nev,
1979  UInt_t run,
1980  Int_t monitor,
1981  Bool_t verbose)
1982 
1983  {
1984  EPosSim* sim = new EPosSim(run, monitor);
1985  sim->fVerbose = verbose;
1986  sim->Begin(0);
1987  sim->SlaveBegin(0);
1988 
1989  for (Long64_t i=0; i <nev; i++) {
1990  Printf("=== Event # %6lld/%6lld ==========================",
1991  i+1, nev);
1992  sim->Process(i);
1993  }
1994  sim->SlaveTerminate();
1995  sim->Terminate();
1996 
1997  return true;
1998  }
2007  static Bool_t SetupProof(const TUrl& url,
2008  const char* opt="")
2009  {
2010  TProof::Reset(url.GetUrl());
2011  TProof::Open(url.GetUrl());
2012  gProof->ClearCache();
2013 
2014  TString phy = gSystem->ExpandPathName("$(ALICE_PHYSICS)");
2015  TString ali = gSystem->ExpandPathName("$(ALICE_ROOT)");
2016  // TString fwd = gSystem->ExpandPathName("$ANA_SRC");
2017  TString fwd = phy + "/PWGLF/FORWARD/analysis2";
2018 
2019  gProof->AddIncludePath(Form("%s/include", ali.Data()));
2020  gProof->AddIncludePath(Form("%s/include", phy.Data()));
2021  ProofLoadLibs();
2022  gProof->Load(Form("%s/sim/GRP.C",fwd.Data()), true);
2023  gProof->Load(Form("%s/sim/BaseConfig.C",fwd.Data()), true);
2024  gProof->Load(Form("%s/sim/EGConfig.C",fwd.Data()), true);
2025 
2026  // gROOT->ProcessLine("gProof->SetLogLevel(5);");
2027  gProof->Load(Form("%s/sim/FastMonitor.C+%s",fwd.Data(),opt));
2028  gProof->Load(Form("%s/sim/FastShortHeader.C", fwd.Data()));
2029  gProof->Load(Form("%s/sim/FastCentEstimators.C+%s",fwd.Data(),opt));
2030  gProof->Load(Form("%s/sim/FastSim.C+%s", fwd.Data(), opt),true);
2031 
2032  return true; // status >= 0;
2033  }
2065  static Bool_t Run(const char* url, const char* opt="")
2066  {
2067  Printf("Will run fast simulation with:\n\n\t%s\n\n",url);
2068  UInt_t run = 0;
2069  Long64_t nev = -1;
2070  Int_t monitor = -1;
2071  Bool_t verbose = false;
2072  TUrl u(url);
2073  TString out;
2074  TObjArray* opts = TString(u.GetOptions()).Tokenize("&");
2075  TObjString* token = 0;
2076  TIter nextToken(opts);
2077  while ((token = static_cast<TObjString*>(nextToken()))) {
2078  TString& str = token->String();
2079  if (str.IsNull()) continue;
2080 
2081  if (str.EqualTo("verbose")) { verbose = true; str = ""; }
2082 
2083  if (str.IsNull()) continue;
2084 
2085  TString key, val;
2086  if (!Str2KeyVal(str,key,val)) {
2087  if (!out.IsNull()) out.Append("&");
2088  out.Append(str);
2089  continue;
2090  }
2091 
2092  if (key.EqualTo("run")) run = val.Atoi();
2093  else if (key.EqualTo("monitor")) monitor = val.Atoi();
2094  else if (key.EqualTo("events")) nev = val.Atoi();
2095  else {
2096  if (!out.IsNull()) out.Append("&");
2097  out.Append(str);
2098  }
2099  }
2100  opts->Delete();
2101  u.SetOptions(out);
2102  if (!u.IsValid()) {
2103  Printf("Error: FastSim::Run: URL %s is invalid", u.GetUrl());
2104  return false;
2105  }
2106 
2107  Printf("Run EPos anchored at %d\n"
2108  " Monitor frequency: %d sec\n"
2109  " Execution url: %s",
2110  run, monitor, u.GetUrl());
2111 
2112  TString treeName = u.GetAnchor();
2113  if (treeName.IsNull()) treeName = "Particle";
2114  TFile* file = TFile::Open(u.GetFile(), "READ");
2115  if (!file) {
2116  Printf("Error: FastAnalysis::Run: Failed to open %s",
2117  u.GetFile());
2118  return false;
2119  }
2120 
2121  TChain* chain = new TChain(treeName, treeName);
2122  TObject* o = file->Get(treeName);
2123  if (!o) {
2124  Printf("Error: FastAnalysis::Run: Couldn't get %s from %s",
2125  treeName.Data(), u.GetFile());
2126  file->Close();
2127  return false;
2128  }
2129  Int_t cret = 0;
2130  if (o->IsA()->InheritsFrom(TChain::Class()))
2131  cret = chain->Add(static_cast<TChain*>(o));
2132  else if (o->IsA()->InheritsFrom(TTree::Class()))
2133  cret = chain->AddFile(u.GetFile());
2134  else if (o->IsA()->InheritsFrom(TCollection::Class()))
2135  cret = chain->AddFileInfoList(static_cast<TCollection*>(o));
2136  else if (o->IsA()->InheritsFrom(TFileCollection::Class()))
2137  cret = chain->AddFileInfoList(static_cast<TFileCollection*>(o)
2138  ->GetList());
2139  file->Close();
2140  if (cret <= 0 || chain->GetListOfFiles()->GetEntries() <= 0) {
2141  Printf("Error: FastAnalysis::Run: Failed to create chain");
2142  return false;
2143  }
2144 
2145  TString proto = u.GetProtocol();
2146  Bool_t isProof = (proto.EqualTo("proof") || proto.EqualTo("lite"));
2147  if (isProof) {
2148  if (!SetupProof(u,opt)) return false;
2149  chain->SetProof();
2150  }
2151 
2152  EPosSim* sim = new EPosSim(run, monitor);
2153  sim->fVerbose = verbose;
2154  if (nev < 0) nev = TChain::kBigNumber;
2155 
2156  TStopwatch timer;
2157  timer.Start();
2158 
2159  Long64_t ret = chain->Process(sim, "", nev, 0);
2160  timer.Print();
2161 
2162  return ret > 0;
2163  }
2164  Int_t Version() const { return 2; }
2165  ClassDef(EPosSim,1);
2166 };
2167 
2168 #endif
2169 //
2170 // EOF
2171 //
EPosSim::PreEvent
Bool_t PreEvent(Long64_t iEv)
Definition: FastSim.C:1843
FastSim::ProofLoadLibs
static void ProofLoadLibs()
Definition: FastSim.C:1463
FastShortHeader::Reset
void Reset(UInt_t runNo, UInt_t eventNo)
Definition: FastShortHeader.C:46
FastSim::SlaveBegin
void SlaveBegin(TTree *)
Definition: FastSim.C:741
pdg
Int_t pdg
Definition: charmCutsOptimization.C:34
FastSim::fBMin
Double_t fBMin
Definition: FastSim.C:1360
str
return jsonbuilder str().c_str()
EPosSim::ClassDef
ClassDef(EPosSim, 1)
FastCentEstimators.C
FastSim::SetupGRP
virtual Bool_t SetupGRP()
Definition: FastSim.C:450
FastSim::ClassDef
ClassDef(FastSim, 3)
FastSim::fIsProjA
Bool_t fIsProjA
True if target beam is nuclei.
Definition: FastSim.C:1366
EPosSim::GetEGTitle
const char * GetEGTitle() const
Definition: FastSim.C:1833
EPosSim::Run
static Bool_t Run(const char *url, const char *opt="")
Definition: FastSim.C:2065
Double_t
double Double_t
Definition: External.C:58
FastSim::Init
virtual void Init(TTree *)
Definition: FastSim.C:676
FastSim::Process
virtual Bool_t Process(Long64_t iEv)
Definition: FastSim.C:1055
EPosSim::fInB
TLeaf * fInB
Definition: FastSim.C:1950
url
const char * url
Definition: ExtractNoiseGain.C:9
EPosSim::fInNSpcPN
TLeaf * fInNSpcPN
Definition: FastSim.C:1962
DtglcpCommon::nct
Int_t nct
Definition: FastSim.C:87
FastSim::Begin
virtual void Begin(TTree *)
Definition: FastSim.C:683
FastSim::FlushList
void FlushList(TCollection *c, TDirectory *dir)
Definition: FastSim.C:1151
Long64_t
long long Long64_t
Definition: External.C:43
FastSim::GetGAlice
TFile * GetGAlice()
Definition: FastSim.C:1248
EPosSim::SetupProof
static Bool_t SetupProof(const TUrl &url, const char *opt="")
Definition: FastSim.C:2007
FastSim::MoveAliceFiles
void MoveAliceFiles()
Definition: FastSim.C:1289
EPosSim::fInNcoll
TLeaf * fInNcoll
Definition: FastSim.C:1961
DtglcpCommon::nwbsam
Int_t nwbsam
Definition: FastSim.C:82
DtglcpCommon::ncp
Int_t ncp
Definition: FastSim.C:86
gSystem
TSystem * gSystem
FastCentEstimator::PostEvent
virtual void PostEvent()
Definition: FastCentEstimators.C:82
FastSim::fBMax
Double_t fBMax
Definition: FastSim.C:1361
FastSim::fStack
AliStack * fStack
Loader of trees.
Definition: FastSim.C:1374
EPosSim::ProcessHeader
Bool_t ProcessHeader()
Definition: FastSim.C:1860
EPosSim::LocalRun
static Bool_t LocalRun(Long64_t nev, UInt_t run, Int_t monitor, Bool_t verbose)
Definition: FastSim.C:1978
FastCentEstimator::Terminate
virtual void Terminate(TCollection *out)=0
FastSim::fGenerator
AliGenerator * fGenerator
Definition: FastSim.C:1372
DtglcpCommon::rproj
Double_t rproj
Definition: FastSim.C:77
EPosSim::fInNTot
TLeaf * fInNTot
Definition: FastSim.C:1949
c
TCanvas * c
Definition: TestFitELoss.C:172
FastSim::fAliceFile
TProofOutputFile * fAliceFile
Proof output file.
Definition: FastSim.C:1409
FastSim::fOverrides
TList * fOverrides
GRP in one line.
Definition: FastSim.C:1363
EPosSim::Init
void Init(TTree *tree)
Definition: FastSim.C:1796
EPosSim::FinishRun
void FinishRun()
Definition: FastSim.C:1947
FastCentEstimator::ProcessHeader
virtual void ProcessHeader(FastShortHeader &)
Definition: FastCentEstimators.C:70
EPosSim::fInTree
TTree * fInTree
Definition: FastSim.C:1948
FastSim::GetTitle
const char * GetTitle() const
Definition: FastSim.C:216
FastSim::fProofFile
TProofOutputFile * fProofFile
Definition: FastSim.C:1408
FastSim::fTree
TTree * fTree
Definition: FastSim.C:1381
DtglcpCommon::nwtbac
Int_t nwtbac
Definition: FastSim.C:85
DtglcpCommon::nwtacc
Int_t nwtacc
Definition: FastSim.C:83
EPosSim::Notify
Bool_t Notify()
Definition: FastSim.C:1817
gRandom
TRandom * gRandom
EPosSim::SetupBranches
Bool_t SetupBranches()
Definition: FastSim.C:1757
FastSim::OverrideGRP
void OverrideGRP()
Definition: FastSim.C:635
TUrl
Definition: External.C:92
FastSim::fVerbose
Bool_t fVerbose
Output file name.
Definition: FastSim.C:1414
FastMonitor::Connect
void Connect(Int_t freq=-1)
Definition: FastMonitor.C:136
EPosSim::fInPhiR
TLeaf * fInPhiR
Definition: FastSim.C:1966
EPosSim::fInPy
TLeaf * fInPy
Definition: FastSim.C:1954
FastCentEstimator::Print
virtual void Print(Option_t *option="") const
Definition: FastCentEstimators.C:89
FastSim::FileName
const char * FileName() const
Definition: FastSim.C:164
FastSim::Str2KeyVal
static Bool_t Str2KeyVal(const TString &in, TString &key, TString &val, const char sep='=')
Definition: FastSim.C:1562
FastSim::Version
Int_t Version() const
Definition: FastSim.C:1352
FastSim::fHB
TH1 * fHB
Event type histogram.
Definition: FastSim.C:1393
FastSim::fIsTgtA
Bool_t fIsTgtA
Definition: FastSim.C:1365
EPosSim::Begin
void Begin(TTree *tree)
Definition: FastSim.C:1783
FastSim::fMonitor
Int_t fMonitor
Definition: FastSim.C:1415
EPosSim::fInNcollH
TLeaf * fInNcollH
Definition: FastSim.C:1958
EPosSim::fInPz
TLeaf * fInPz
Definition: FastSim.C:1955
FastSim::FinishRun
virtual void FinishRun()
Definition: FastSim.C:1085
FastSim::GetKine
TFile * GetKine()
Definition: FastSim.C:1265
Int_t
int Int_t
Definition: External.C:63
FastSim::fBEstimator
BCentEstimator * fBEstimator
Definition: FastSim.C:1401
EPosSim::SetupGen
Bool_t SetupGen()
Definition: FastSim.C:1835
TH1I
Definition: External.C:204
UInt_t
unsigned int UInt_t
Definition: External.C:33
EPosSim::fInNpartP
TLeaf * fInNpartP
Definition: FastSim.C:1959
FastSim::fGRP
TObject * fGRP
Definition: FastSim.C:1362
Float_t
float Float_t
Definition: External.C:68
kNeutron
Definition: FastCentEstimators.C:22
FastSim::fCentEstimators
TList * fCentEstimators
Definition: FastSim.C:1402
V0CentEstimator::Flip
void Flip(Bool_t onlySign=true)
Definition: FastCentEstimators.C:565
FastSim::fShortHead
FastShortHeader fShortHead
Definition: FastSim.C:1418
_dtglcp
DtglcpCommon * _dtglcp
Definition: FastSim.C:89
FastShortHeader.C
EPosSim::Generate
void Generate()
Definition: FastSim.C:1946
FastSim::fHIpz
TH1 * fHIpz
dN/dy
Definition: FastSim.C:1390
FastSim::fHPhiR
TH1 * fHPhiR
B histogram.
Definition: FastSim.C:1394
BCentEstimator::GetCentrality
Double_t GetCentrality(Double_t b) const
Definition: FastCentEstimators.C:473
EPosSim::SetupRun
Bool_t SetupRun()
Definition: FastSim.C:1842
EPosSim::fInNpartT
TLeaf * fInNpartT
Definition: FastSim.C:1960
EPosSim::SetupSeed
void SetupSeed()
Definition: FastSim.C:1834
TH1D
Definition: External.C:212
DtglcpCommon::nwtaac
Int_t nwtaac
Definition: FastSim.C:84
FastSim::GetName
const char * GetName() const
Definition: FastSim.C:210
FastSim::GetSNN
virtual UShort_t GetSNN() const
Definition: FastSim.C:222
FastSim::FastSim
FastSim(const char *eg="", ULong_t runNo=0, Double_t bMin=0, Double_t bMax=20, Long64_t nEvents=0, Int_t monitor=0)
Definition: FastSim.C:123
FastSim::fHType
TH1 * fHType
IPz histogram.
Definition: FastSim.C:1391
FastSim::SetupSeed
virtual void SetupSeed()
Definition: FastSim.C:439
kProton
Definition: FastCentEstimators.C:23
ULong_t
unsigned long ULong_t
Definition: External.C:38
DtglcpCommon::rtarg
Double_t rtarg
Definition: FastSim.C:78
FastSim::fKineFile
TProofOutputFile * fKineFile
Definition: FastSim.C:1410
FastSim::fHeader
AliHeader * fHeader
Stack of particles.
Definition: FastSim.C:1375
FastSim::fHEta
TH1 * fHEta
List of outputs.
Definition: FastSim.C:1388
FastSim::ReadGRPLine
Bool_t ReadGRPLine()
Definition: FastSim.C:603
FastSim::Generate
virtual void Generate()
Definition: FastSim.C:1046
FastSim::fRunNo
Int_t fRunNo
Definition: FastSim.C:1359
EPosSim::fInPDG
TLeaf * fInPDG
Definition: FastSim.C:1951
FastSim::fList
TList * fList
Definition: FastSim.C:1387
FastSim::Run
static Bool_t Run(const char *url, const char *opt="")
Definition: FastSim.C:1648
FastSim::SetupGen
virtual Bool_t SetupGen()
Definition: FastSim.C:495
FastSim::fHTime
TH1 * fHTime
Reaction plane.
Definition: FastSim.C:1395
EPosSim::Version
Int_t Version() const
Definition: FastSim.C:2164
FastSim::fFile
TFile * fFile
Definition: FastSim.C:1411
EPosSim::fInStatus
TLeaf * fInStatus
Definition: FastSim.C:1952
EPosSim::fInM
TLeaf * fInM
Definition: FastSim.C:1957
fwd
const char * fwd
Definition: ExtractNoiseGain.C:7
FastCentEstimator::SetVerbose
void SetVerbose(Bool_t verb)
Definition: FastCentEstimators.C:52
DtglcpCommon::nwasam
Int_t nwasam
Definition: FastSim.C:81
GetListOfFiles
TList * GetListOfFiles(const char *input=".")
Definition: RunFinalQA.C:114
FastCentEstimator::Process
virtual void Process(const TParticle *p)=0
FastMonitor::Register
void Register(TObject *descr, Bool_t proof=true)
Definition: FastMonitor.C:186
FastShortHeader::fEG
enum FastShortHeader::EG_t fEG
FastSim::Terminate
void Terminate()
Definition: FastSim.C:1173
FastSim::PreEvent
virtual Bool_t PreEvent(Long64_t iEv)
Definition: FastSim.C:756
EPosSim::fInNSpcTN
TLeaf * fInNSpcTN
Definition: FastSim.C:1963
FastMonitor.C
EPosSim::fInNSpcTP
TLeaf * fInNSpcTP
Definition: FastSim.C:1965
nEvents
Float_t nEvents[nProd]
Definition: DrawProductionComparison.C:36
DtglcpCommon::bimpac
Double_t bimpac
Definition: FastSim.C:79
EPosSim::EPosSim
EPosSim(UInt_t run=0, Int_t monitor=0)
Definition: FastSim.C:1734
FastSim::ProofRun
static Bool_t ProofRun(const TUrl &url, Long64_t nev, UInt_t run, const TString &gen, Double_t bMin, Double_t bMax, Int_t monitor=-1, Bool_t verbose=false, const TString &overrides="", const TString &save="none", const char *opt="")
Definition: FastSim.C:1510
EPosSim::fInPx
TLeaf * fInPx
Definition: FastSim.C:1953
FastSimMonitor::FastSimMonitor
FastSimMonitor(TSelector *s=0)
Definition: FastSim.C:94
FastSim::fRunLoader
AliRunLoader * fRunLoader
Event generator.
Definition: FastSim.C:1373
EPosSim::PostEvent
void PostEvent()
Definition: FastSim.C:1937
file
TFile * file
Definition: DrawAnaCaloTrackQA.C:19
FastSim::PostEvent
virtual void PostEvent()
Definition: FastSim.C:1028
FastSim::fHY
TH1 * fHY
dN/deta
Definition: FastSim.C:1389
FastSim::fFileName
TString fFileName
Output file.
Definition: FastSim.C:1412
FastSim::SetupOutput
Bool_t SetupOutput()
Definition: FastSim.C:254
FastSim::LocalRun
static Bool_t LocalRun(Long64_t nev, UInt_t run, const TString &gen, Double_t bMin, Double_t bMax, Int_t monitor, Bool_t verbose, const TString &overrides="")
Definition: FastSim.C:1434
UShort_t
unsigned short UShort_t
Definition: External.C:28
FastSim::fParticles
TClonesArray * fParticles
Custom tree.
Definition: FastSim.C:1382
DtglcpCommon::nwtsam
Int_t nwtsam
Definition: FastSim.C:80
FastSim::SetupRun
virtual Bool_t SetupRun()
Definition: FastSim.C:534
Bool_t
bool Bool_t
Definition: External.C:53
FastSim::fEGName
TString fEGName
Definition: FastSim.C:1358
FastCentEstimator::PreEvent
virtual void PreEvent()
Definition: FastCentEstimators.C:69
FastSim::AddOverride
void AddOverride(const TString &field, const TString &value)
Definition: FastSim.C:664
FastSim::GetEGTitle
virtual const char * GetEGTitle() const
Definition: FastSim.C:243
FastSim::ProcessHeader
virtual Bool_t ProcessHeader()
Definition: FastSim.C:780
FastSim::ProcessParticles
virtual Bool_t ProcessParticles(Bool_t selected)
Definition: FastSim.C:984
EPosSim::ProcessParticles
virtual Bool_t ProcessParticles(Bool_t selected)
Definition: FastSim.C:1897
FastSim::SlaveTerminate
void SlaveTerminate()
Definition: FastSim.C:1097
FastCentEstimator::Setup
virtual void Setup(TCollection *out, TTree *tree, UShort_t sNN, Bool_t tgtA, Bool_t projA)=0
FastSim::fHCent
TH1 * fHCent
Event type histogram.
Definition: FastSim.C:1392
TH1
Definition: External.C:196
EPosSim::fInNSpcPP
TLeaf * fInNSpcPP
Definition: FastSim.C:1964
FastSim::fNEvents
Long64_t fNEvents
GRP setting to override.
Definition: FastSim.C:1364
FastSim::SetOverrides
static void SetOverrides(FastSim *sim, const TString &override)
Definition: FastSim.C:1574
EPosSim::fInE
TLeaf * fInE
Definition: FastSim.C:1956