AliPhysics  8d00e07 (8d00e07)
EventTimeTask.C
Go to the documentation of this file.
1 // EventTimeTask.C
2 #ifndef __CINT__
3 # include <TTree.h>
4 # include <TError.h>
5 # include <TChain.h>
6 #else
7 class TTree;
8 // Force load of libGui
9 class TGWindow;
10 #endif
11 class AliESDEvent;
12 
19 {
23  enum {
24  kBCWidth = 12,
25  kOrbitWidth = 24,
27  };
29  ULong64_t fFull;
37  UChar_t fGUID[42];
43  void CreateBranch(TTree* tree)
44  {
45  tree->Branch("event", &(this->fFull),
46  "full/l:time/i:detector:type/s:guid/C");
47  }
53  void ReadBranch(TTree* tree)
54  {
55  tree->SetBranchAddress("event", &(this->fFull));
56  }
66  static ULong64_t EncodeFull(ULong64_t period,
67  ULong64_t orbit,
68  ULong64_t bc)
69  {
70  const ULong64_t bcMask = (1 << kBCWidth) - 1;
71  const ULong64_t orbitMask = (1 << kOrbitWidth) - 1;
72  const ULong64_t periodMask = (1 << kPeriodWidth) - 1;
73  ULong64_t ret = ((bc & bcMask) |
74  (orbit & orbitMask) << kBCWidth |
75  (period & periodMask) << (kBCWidth+kOrbitWidth));
76  return ret;
77  }
85  void Fill(AliESDEvent* esd, UInt_t dets, const TString& guid);
86 };
87 
88 #ifndef NO_TASK
89 # ifndef __CINT__
90 # include <AliESDEvent.h>
91 # endif
93  const TString& guid)
94 {
95  ULong64_t period = esd->GetPeriodNumber();
96  ULong64_t orbit = esd->GetOrbitNumber();
97  ULong64_t bc = esd->GetBunchCrossNumber();
98  fType = esd->GetEventType();
99  fTime = esd->GetTimeStamp();//LDC time
100  fDetectors = dets; // esd->GetDAQDetectorPattern();
101  fFull = EncodeFull(period, orbit, bc);
102  Int_t i = 0;
103  for (i = 0; i < guid.Length() && i < 42; i++) fGUID[i] = guid[i];
104  for (; i < 42; i++) fGUID[i] = '\0';
105  fGUID[41] = '\0';
106 }
107 # else
108 inline void EventTimeData::Fill(AliESDEvent*, UInt_t, const TString&)
109 {
110  Warning("Fill", "Calling empty method - shouldn't happen");
111 }
112 #endif
113 
114 #ifndef NO_TASK
115 # ifndef __CINT__
116 # include <AliAnalysisManager.h>
117 # include <AliVEventHandler.h>
118 # include <AliESDEvent.h>
119 # include <TTree.h>
120 # include <TH2.h>
121 # include <TList.h>
122 # include <AliCDBManager.h>
123 # include <AliCDBEntry.h>
124 # include <AliTriggerConfiguration.h>
125 # include <AliTriggerClass.h>
126 # include <AliTriggerCluster.h>
127 # include <AliDAQ.h>
128 # include <TObjArray.h>
129 # include <TDirectory.h>
130 # include <TUrl.h>
131 # else
132 class AliAnalysisManager;
133 class TTree;
134 class AliTimeStamp;
135 class TList;
136 class TH2;
137 # endif
138 # include <AliAnalysisTaskSE.h>
139 # include <vector>
140 
149 {
150 public:
151  enum {
152  kListSlot = 1,
153  kTreeSlot = 2
154  };
159  : AliAnalysisTaskSE(),
160  fTree(0),
161  fHistograms(0),
162  fDetVsType(0)
163  {
164  }
170  EventTimeTask(const char* name)
171  : AliAnalysisTaskSE(name),
172  fTree(0),
173  fHistograms(0),
174  fDetVsType(0),
175  fGUID("")
176  {
177  DefineOutput(kListSlot, TList::Class());
178  DefineOutput(kTreeSlot, TTree::Class());
179  fBranchNames = "ESD:AliESDRun.,AliESDHeader.";
180  }
187  {
188  Printf("Creating tree and histogram");
189  fGUID = "";
190  fHistograms = new TList();
191  fHistograms->SetOwner();
192  fHistograms->SetName("L");
193 
194  fDetVsType = new TH2D("detVsType", "Detector vs type",
195  16, -.5, 15.5, 31, -.5, 30.5);
196  fDetVsType->SetXTitle("Type");
197  fDetVsType->SetYTitle("Detector");
198  fDetVsType->SetDirectory(0);
199  fHistograms->Add(fDetVsType);
200  Printf("Histogram (%d,%f,%f)x(%d,%f,%f)",
201  fDetVsType->GetXaxis()->GetNbins(),
202  fDetVsType->GetXaxis()->GetXmin(),
203  fDetVsType->GetXaxis()->GetXmax(),
204  fDetVsType->GetYaxis()->GetNbins(),
205  fDetVsType->GetYaxis()->GetXmin(),
206  fDetVsType->GetYaxis()->GetXmax());
207 
208  // TDirectory* savdir = gDirectory;
209  // Printf("Opening file at slot %d", kTreeSlot);
210  // OpenFile(kTreeSlot);
211  Printf("Make tree and disassociate from file");
212  fTree = new TTree("T", "T");
213  fTree->SetDirectory(0);
214  Printf("Create branch");
215  fData.CreateBranch(fTree);
216  // savdir->cd();
217 
218 
219  PostData(kListSlot, fHistograms);
220  PostData(kTreeSlot, fTree);
221  }
227  {
228  static Bool_t first = true;
229  LoadBranches();
230 
231  AliESDEvent* esd = dynamic_cast<AliESDEvent*>(InputEvent());
232  if (!esd) return;
233  if (!esd->GetHeader()) return;
234 
235  if (first) {
236  LoadTriggerConfig(esd->GetRunNumber());
237  first = false;
238  }
239  ULong64_t mask = esd->GetTriggerMask();
240  UInt_t dets = 0;
241  for (UShort_t i = 0; i < fDets.size(); i++) {
242  if ((1 << i) & mask) dets |= fDets[i];
243  }
244  // Printf("Event mask 0x%016llx -> 0x%08x", mask, dets);
245  fData.Fill(esd, dets, fGUID);
246  fTree->Fill();
247 
248  UInt_t type = esd->GetEventType();
249  UInt_t detectors = dets;
250  for (UInt_t i = 0; i < 31; i++) {
251  if ((1 << i) & detectors) fDetVsType->Fill(type, i);
252  }
253 
254  PostData(kListSlot, fHistograms);
255  PostData(kTreeSlot, fTree);
256  }
258  {
259  fGUID = "";
260  AliAnalysisManager* mgr = AliAnalysisManager::GetAnalysisManager();
261  AliVEventHandler* inp = mgr->GetInputEventHandler();
262  if (!inp) {
263  Warning("UserNotify", "No input handler");
264  return true;
265  }
266  TTree* tree = inp->GetTree();
267  if (!tree) {
268  Warning("UserNotify", "No input tree");
269  return true;
270  }
271  TFile* file = tree->GetCurrentFile();
272  if (!file) {
273  Warning("UserNotify", "No current file for tree");
274  return true;
275  }
276  const TUrl* url = file->GetEndpointUrl();
277  if (!url) {
278  Warning("UserNotify", "No end point for file");
279  return false;
280  }
281  fGUID = gSystem->BaseName(url->GetFile());
282  Printf("Got GUID=%s from %s", fGUID.Data(), url->GetUrl());
283 
284  return true;
285  }
287  {
288  Printf("Loading trigger configuration for run %d", runNo);
289  // --- Connect to CDB --------------------------------------------
290  AliCDBManager* cdb = AliCDBManager::Instance();
291  cdb->SetDefaultStorageFromRun(runNo);
292  cdb->SetRun(runNo);
293 
294  // --- Get entry -------------------------------------------------
295  AliCDBEntry* entry = cdb->Get("GRP/CTP/Config");
296  if (!entry || !entry->GetObject()) {
297  Warning("LoadTriggerConfig", "Couldn't get trigger configuration");
298  return;
299  }
300  AliTriggerConfiguration* config =
301  static_cast<AliTriggerConfiguration*>(entry->GetObject());
302 
303  // --- Get the classes, and resize cache -------------------------
304  const TObjArray& clss = config->GetClasses();
305  fDets.resize(clss.GetEntries());
306 
307  // --- Loop over configurations ----------------------------------
308  TIter next(&clss);
309  AliTriggerClass* cls = 0;
310  while ((cls = static_cast<AliTriggerClass*>(next()))) {
311  Int_t mask = cls->GetMask();
312  AliTriggerCluster* cluster = cls->GetCluster();
313  if (!cluster) {
314  Warning("LoadTriggerConfig",
315  "Failed to get trigger cluster for %s", cls->GetName());
316  continue;
317  }
318  TString names = cluster->GetDetectorsInCluster();
319  if (names.IsNull()) {
320  Warning("LoadTriggerConfig", "No detectors for cluster %s",
321  cls->GetName());
322  continue;
323  }
324  UInt_t dets = AliDAQ::DetectorPattern(names);
325  UInt_t id = UInt_t(TMath::Log2(mask));
326  fDets[id] = dets;
327  Printf("Trigger mask 0x%08x (%3d): 0x%08x (%s)",
328  mask, id, dets, names.Data());
329  }
330  // cdb->Destroy();
331  }
338  {
339  mgr->AddTask(this);
340  mgr->ConnectInput(this, 0, mgr->GetCommonInputContainer());
341  AliAnalysisDataContainer* contTree =
342  mgr->CreateContainer("T", TTree::Class(),
343  AliAnalysisManager::kOutputContainer,
344  "time.root");
345  AliAnalysisDataContainer* contList =
346  mgr->CreateContainer("L", TList::Class(),
347  AliAnalysisManager::kOutputContainer,
348  "hist.root");
349  mgr->ConnectOutput(this, kTreeSlot, contTree);
350  mgr->ConnectOutput(this, kListSlot, contList);
351  }
356  static void Create()
357  {
358  EventTimeTask* task = new EventTimeTask("eventTime");
359  task->Connect(AliAnalysisManager::GetAnalysisManager());
360  }
361  TTree* fTree; // Our tree
362  EventTimeData fData; // Our data
363  TList* fHistograms; // List
364  TH2D* fDetVsType; // Histogram
365  std::vector<UInt_t> fDets; // Per-trigger-bit detector mask
367 
368  ClassDef(EventTimeTask,3);
369 };
370 #endif // NO_TASK
371 
372 #ifndef NO_MAP
373 #include <utility>
374 #include <map>
375 
376 typedef std::pair<ULong64_t,ULong64_t> EventTimeMapPair;
377 
383 struct EventTimeMap : public TObject
384 {
386  typedef std::map<ULong64_t,ULong64_t> Map;
388  typedef Map::key_type Key;
390  typedef Map::mapped_type Value;
392  typedef Map::value_type Pair;
394  typedef Map::iterator Iterator;
396  typedef Map::const_iterator ConstIterator;
400  EventTimeMap() : TObject(), fMap() {}
404  virtual ~EventTimeMap() {}
410  const char* GetName() const { return "eventTimeMap"; }
419  Value& operator[](const Key& k)
420  {
421  return fMap[k];
422  }
430  Iterator Find(const Key& k)
431  {
432  return fMap.find(k);
433  }
441  ConstIterator Find(const Key& k) const
442  {
443  return fMap.find(k);
444  }
450  Iterator Begin()
451  {
452  return fMap.begin();
453  }
459  ConstIterator Begin() const
460  {
461  return fMap.begin();
462  }
468  Iterator End()
469  {
470  return fMap.end();
471  }
477  ConstIterator End() const
478  {
479  return fMap.end();
480  }
481  enum {
482  kInvalidTime = 0xFFFFFFFFFFFFFFFF
483  };
492  Value Get(const Key& timestamp) const
493  {
494  ConstIterator i = Find(timestamp);
495  if (i == End()) return kInvalidTime;
496  return i->second;
497  }
504  void Set(const Key& timestamp, const Value& diff)
505  {
506  this->operator[](timestamp) = diff;
507  }
508  ULong64_t Size() const { return fMap.size(); }
510  Map fMap;
511  ClassDef(EventTimeMap,1)
512 };
513 #endif // NO_MAP
514 
515 #ifndef NO_SORTER
516 # ifndef __CINT__
517 # include <TFile.h>
518 # include <TArrayL64.h>
519 # include <TMath.h>
520 # include <TParameter.h>
521 # include <TCanvas.h>
522 # include <TH1.h>
523 # else
524 class TFile;
525 class TTree;
526 class TCanvas;
527 # endif
528 
535 {
538  ULong64_t fMin;
539  ULong64_t fMax;
540 
544  EventTimeSorter() : fTree(0), fData(), fMin(0xFFFFFFFFFFFFFFFF), fMax(0) {}
551  void Progress(Long64_t cur, Long64_t total) const
552  {
553  static UShort_t old = 0;
554  if (cur == 0) old = 0;
555  UShort_t now = UShort_t(100 * (cur + 1 == total ? 1 :
556  Double_t(cur)/total));
557  if (now > old)
558  std::cout << "\rLooping over " << total << " entries: ... "
559  << now << '%' << std::flush;
560  if (now >= 100) std::cout << std::endl;
561  old = now;
562  }
571  Bool_t OpenInput(const char* inputName, const char* treeName)
572  {
573  CloseInput();
574 
575  // --- Get input -------------------------------------------------
576  TChain* chain = new TChain(treeName);
577  chain->Add(inputName);
578  if (chain->GetListOfFiles()->GetEntries() < 1) {
579  Error("Run", "Failed to add \"%s\" to chain", inputName);
580  return false;
581  }
582  fTree = chain;
583 
584  // --- Set branch address ---------------------------------------
585  fData.ReadBranch(fTree);
586 
587  return true;
588  }
592  void CloseInput()
593  {
594  if (!fTree) return;
595  TFile* file = fTree->GetCurrentFile();
596  if (file) file->Close();
597  fTree = 0;
598  }
608  Bool_t Run(const char* inputName, const char* outputName,
609  const char* treeName="T")
610  {
611  if (!OpenInput(inputName, treeName)) return false;
612 
613  // --- Loop over the data ----------------------------------------
614  Long64_t nEnt = fTree->GetEntries();
615  Long64_t n = 0;
616  ULong64_t* values = new ULong64_t[nEnt];
617  UInt_t* dets = new UInt_t[nEnt];
618  UShort_t* types = new UShort_t[nEnt];
619  UInt_t* times = new UInt_t[nEnt];
620  for (Long64_t iEnt = 0; iEnt < nEnt; iEnt++) {
621  Progress(iEnt, nEnt);
622 
623  fTree->GetEntry(iEnt);
624 
625  if (!(fData.fDetectors & 0x1000))
626  // No FMD read-out
627  continue;
628  if (fData.fType != 7) {
629  // Ignore all but physics events
630  Warning("Run", "Non-PHYSICS event (%d) with FMD in it", fData.fType);
631  // continue;
632  }
633 
634  // Store values
635  Int_t j = n;
636  values[j] = fData.fFull;
637  dets[j] = fData.fDetectors;
638  types[j] = fData.fType;
639  times[j] = fData.fTime;
640  n++;
641  }
642 
643  // --- Now sort all values ---------------------------------------
644  TArrayL64 index(n);
645  TMath::Sort(n, values, index.fArray, false);
646 
647  // Maps the unique time to the distance to previous event
648  EventTimeMap eventTimeMap;
649  ULong64_t key = values[index[0]];
650  ULong64_t prev = 0;
651  ULong64_t dt = key-prev;
652  ULong64_t min = 0xFFFFFFFFFFFFFFFF;
653  ULong64_t max = 0;
654  eventTimeMap[key] = dt;
655  for (Long64_t i = 1; i < n; i++) {
656  Long64_t j = index[i];
657  Long64_t k = index[i-1];
658  key = values[j];
659  prev = values[k];
660  dt = key - prev;
661  if (dt <= 0) {
662  Printf("0x%016llx==0x%016llx -> dt=%10llu [%10lld %10lld] "
663  "(det: 0x%08x 0x%08x type: 0x%04x 0x%04x time: %08d %08d)",
664  key, prev, dt, j, k, dets[j], dets[k],
665  types[j], types[k], times[j], times[k]);
666  // continue;
667  }
668  eventTimeMap[key] = dt;
669  min = TMath::Min(dt, min);
670  max = TMath::Max(dt, max);
671  }
672  std::cout << "Range is [" << min << "," << max << ']' << std::endl;
673 
674  TFile* output = TFile::Open(outputName, "RECREATE");
675  if (!output) {
676  Error("Run", "Failed to create output file \"%s\"", outputName);
677  return false;
678  }
679  fMin = min;
680  fMax = max;
681  eventTimeMap.Write();
682  output->Write();
683  output->Close();
684 
685  delete [] values;
686 
687  CloseInput();
688 
689  return true;
690  }
691  Bool_t Test(const char* inputName, const char* outputName,
692  const char* treeName="T")
693  {
694  if (!OpenInput(inputName, treeName)) return false;
695 
696  // --- Get our map -----------------------------------------------
697  TFile* output = TFile::Open(outputName, "UPDATE");
698  if (!output) {
699  Error("Test", "Failed to open \"%s\"", outputName);
700  return false;
701  }
702 
703  EventTimeMap* eventTimeMap =
704  static_cast<EventTimeMap*>(output->Get("eventTimeMap"));
705  if (!eventTimeMap) {
706  Error("Test", "Failed to get \"%s\" from \"%s\"",
707  "eventTimeMap", outputName);
708  return false;
709  }
710 
711  // --- Histograms --------------------------------------------------
712  ULong64_t mmin = TMath::Min(25*fMin, 900000ULL);
713  ULong64_t mmax = TMath::Min(25*fMax, 110000000ULL);
714  TH1* diff = new TH1D("diff", "Time-to-last-event (10#mus bins w/FMD)",
715  (mmax-mmin)/10000, mmin, mmax);
716  diff->SetStats(0);
717  diff->SetXTitle("#Deltat [ns]");
718  diff->SetFillColor(kRed+1);
719  diff->SetFillStyle(3001);
720  diff->SetDirectory(0);
721 
722 
723  TH1* ldiff = new TH1D("ldiff", "log(#Deltat) - Events w/FMD",
724  300, 0, 15);
725  ldiff->SetStats(0);
726  ldiff->SetXTitle("log_{10}(#Deltat) [ns]");
727  ldiff->SetFillColor(kGreen+1);
728  ldiff->SetFillStyle(3001);
729  ldiff->SetDirectory(0);
730 
731  // --- Loop over the data ----------------------------------------
732  Long64_t nEnt = fTree->GetEntries();
733  Long64_t nZero = 0;
734  Long64_t nMiss = 0;
735  Long64_t nGood = 0;
736  Long64_t nNoFMD = 0;
737  for (Long64_t iEnt = 0; iEnt < /*10*/ nEnt; iEnt++) {
738  Progress(iEnt, nEnt);
739  fTree->GetEntry(iEnt);
740 
741  if (!(fData.fDetectors & 0x1000)) {
742  // No FMD read-out
743  nNoFMD++;
744  continue;
745  }
746  if (fData.fType != 7) {
747  // Ignore all but physics events
748  Warning("Run", "Non-PHYSICS event (%d) with FMD in it", fData.fType);
749  // continue;
750  }
751 
752  // Look-up the timestamp
753  ULong64_t value = fData.fFull;
754  ULong64_t dT = eventTimeMap->Get(value);
755  if (dT == EventTimeMap::kInvalidTime) {
756  Warning("Test", "Value %llu not found", value);
757  ldiff->Fill(1);
758  nMiss++;
759  continue;
760  }
761  if (dT <= 0) {
762 #if 0
763  Warning("Test", "Impossible dt=%llu found for %llu (0x%0x %2d)",
764  dT, value, fData.fDetectors, fData.fType);
765 #endif
766  ldiff->Fill(1);
767  nZero++;
768  continue;
769  }
770  nGood++;
771  Double_t dt = 25.*dT;
772  diff->Fill(dt);
773  Double_t logDt = TMath::Log10(dt);
774  ldiff->Fill(logDt);
775  }
776  CloseInput();
777  Printf("missing: %llu, bad: %llu, good: %llu, no FMD: %llu, all: %llu",
778  nMiss,nZero,nGood,nNoFMD,nEnt);
779 
780  // --- Draw results ----------------------------------------------
781  TCanvas* c = new TCanvas("c", "C");
782  c->SetTopMargin(0.01);
783  c->SetRightMargin(0.01);
784  c->Divide(2,1); // ,0,0);
785  TVirtualPad* p = c->cd(2);
786  // p->SetRightMargin(0.10);
787  p->SetLogy();
788  ldiff->DrawCopy();
789 
790  p = c->cd(1);
791  p->SetLogy();
792  p->SetLogx();
793  diff->DrawCopy();
794 
795  c->Print("dt.png");
796  c->Print("dt.C");
797 
798  // --- Write histogram to file -----------------------------------
799  output->cd();
800  diff->Write();
801  output->Write();
802  output->Close();
803 
804  return true;
805  }
806 };
807 #endif // NO_SORTER
808 #ifdef __MAKECINT__
809 # ifndef NO_MAP
810 # pragma link C++ class std::pair<ULong64_t,ULong64_t>+;
811 # pragma link C++ class std::map<ULong64_t,ULong64_t>+;
812 # endif
813 # ifndef NO_TASK
814 # pragma link C++ class std::vector<UInt_t>+;
815 # endif
816 #endif
817 
818 // EOF
819 
820 
void ReadBranch(TTree *tree)
Definition: EventTimeTask.C:53
double Double_t
Definition: External.C:58
const char * url
long long Long64_t
Definition: External.C:43
std::vector< UInt_t > fDets
void CreateBranch(TTree *tree)
Definition: EventTimeTask.C:43
TSystem * gSystem
UChar_t fGUID[42]
Definition: EventTimeTask.C:37
Bool_t UserNotify()
ULong64_t fFull
Definition: EventTimeTask.C:29
static ULong64_t EncodeFull(ULong64_t period, ULong64_t orbit, ULong64_t bc)
Definition: EventTimeTask.C:66
TCanvas * c
Definition: TestFitELoss.C:172
Bool_t OpenInput(const char *inputName, const char *treeName)
Int_t types
ConstIterator Begin() const
void Connect(AliAnalysisManager *mgr)
Bool_t Run(const char *inputName, const char *outputName, const char *treeName="T")
Definition: External.C:92
Bool_t Test(const char *inputName, const char *outputName, const char *treeName="T")
void LoadTriggerConfig(Int_t runNo)
Map::const_iterator ConstIterator
int Int_t
Definition: External.C:63
static void Create()
unsigned int UInt_t
Definition: External.C:33
TList * fHistograms
void Set(const Key &timestamp, const Value &diff)
std::map< ULong64_t, ULong64_t > Map
void UserCreateOutputObjects()
void Fill(AliESDEvent *esd, UInt_t dets, const TString &guid)
Definition: EventTimeTask.C:92
Definition: External.C:228
Definition: External.C:212
Iterator Find(const Key &k)
EventTimeData fData
ULong64_t Size() const
std::pair< ULong64_t, ULong64_t > EventTimeMapPair
Map::mapped_type Value
Map::value_type Pair
ConstIterator Find(const Key &k) const
Map::iterator Iterator
ConstIterator End() const
Definition: External.C:220
TFile * file
TList with histograms for a given trigger.
UShort_t fType
Definition: EventTimeTask.C:35
const char * dets[]
Definition: PlotSysInfo.C:138
void Progress(Long64_t cur, Long64_t total) const
unsigned short UShort_t
Definition: External.C:28
const char Option_t
Definition: External.C:48
UInt_t fDetectors
Definition: EventTimeTask.C:33
EventTimeTask(const char *name)
Value & operator[](const Key &k)
EventTimeData fData
Iterator Begin()
bool Bool_t
Definition: External.C:53
const char * GetName() const
virtual ~EventTimeMap()
Map::key_type Key
Value Get(const Key &timestamp) const
void UserExec(Option_t *)
Definition: External.C:196
Iterator End()