1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
19 // AliAnalysisTriggerScalers : utility class to play with OCDB scalers.
21 // Can use it to e.g. :
23 // - get the integrated luminosity for a given period
24 // (see IntegratedLuminosity method)
26 // - plot the trigger rates (see PlotTriggerEvolution)
28 // Please note that this class is doing an OCDB "scan" (for as many run
29 // as you put in the SetRunList method), so it can be really long
30 // if you're using the raw:// OCDB. If you're planning on using this
31 // class for anything more than a brief test, you'd better think
32 // or making a local copy of the OCDB objects required by this class :
37 // and also (to get the fill and period ranges, mainly for drawing)
43 // Please note also that this class evolved from a set a quick macros to
44 // follow the luminosity and trigger rates during data taking to this stage.
46 // Now (Feb 2013) would probably be a good time to regroup a bit, think about it and
47 // make it more general, more robust and just less "had-hoc"...
49 // Ideas for improvement :
51 // - make sure the "what" parameter mean the same thing in all methods
52 // and so can take the same values in all methods
54 // - get the lumi trigger name and cross-section from e.g. OADB instead of
57 // - organize the class so that the CTP Scalers and Config can be fetched
58 // from elsewhere (e.g. from a run based container available in some
59 // distant future in the ESDs/AODs ?)
61 // - robustify the PAC fraction computation
64 // L. Aphecetche (Subatech)
67 #include "AliAnalysisTriggerScalers.h"
68 #include "AliCDBEntry.h"
69 #include "AliCDBManager.h"
70 #include "AliGRPObject.h"
72 #include "AliTimeStamp.h"
73 #include "AliTriggerBCMask.h"
74 #include "AliTriggerClass.h"
75 #include "AliTriggerCluster.h"
76 #include "AliTriggerConfiguration.h"
77 #include "AliTriggerDescriptor.h"
78 #include "AliTriggerInput.h"
79 #include "AliTriggerRunScalers.h"
80 #include "AliTriggerScalers.h"
81 #include "AliTriggerScalersESD.h"
82 #include "AliTriggerScalersRecord.h"
83 #include "Riostream.h"
87 #include "TGraphErrors.h"
89 #include "TObjArray.h"
90 #include "TObjString.h"
94 #include "AliLHCData.h"
101 ClassImp(AliAnalysisTriggerScalers)
105 //______________________________________________________________________________
106 UChar_t GetIndex(ULong64_t mask)
108 for ( Int_t i = 0; i < 64; ++i )
110 if ( mask & ( 1ull << i ) ) return i+1;
116 //______________________________________________________________________________
117 Double_t Mu(Double_t L0B, Double_t Nb)
119 /// L0B = trigger rate before any veto
120 /// Nb = number of crossing bunches
122 Double_t p0 = 1-L0B/(Nb*11245.0); // proba to get *no* collision per bunch crossing
124 return -TMath::Log(p0);
127 //______________________________________________________________________________
128 TCanvas* NewCanvas(const char* name)
130 TCanvas* c = new TCanvas(name,name);
131 c->SetLeftMargin(0.15);
135 //______________________________________________________________________________
136 void TimeAxis(TGraph* g)
138 // g->GetXaxis()->SetTimeDisplay(1);
139 // g->GetXaxis()->SetTimeFormat("%d/%m %H:%M%F2010-12-31 24:00:00");
140 // g->GetXaxis()->SetNdivisions(505);
141 g->GetYaxis()->SetTitleOffset(1.5);
143 g->GetXaxis()->SetTimeDisplay(1);
144 // g->GetXaxis()->SetTimeFormat("%d/%m %H:%M%F2010-12-31 24:00:00");
145 g->GetXaxis()->SetTimeFormat("%d/%m %H:%M");
146 g->GetXaxis()->SetTimeOffset(0,"gmt");
147 g->GetXaxis()->SetNdivisions(505);
153 //_____________________________________________________________________________
154 AliAnalysisTriggerScalers::AliAnalysisTriggerScalers(Int_t runNumber, const char* ocdbPath) : fRunList(), fOCDBPath(ocdbPath), fShouldCorrectForPileUp(kTRUE), fCrossSectionUnit("UB")
156 // ctor using a single run number
157 SetRunList(runNumber);
160 //_____________________________________________________________________________
161 AliAnalysisTriggerScalers::AliAnalysisTriggerScalers(const std::vector<int>& runs, const char* ocdbPath) : fRunList(), fOCDBPath(ocdbPath), fShouldCorrectForPileUp(kTRUE), fCrossSectionUnit("UB")
163 // ctor using a vector of run numbers
167 //_____________________________________________________________________________
168 AliAnalysisTriggerScalers::AliAnalysisTriggerScalers(const std::set<int>& runs, const char* ocdbPath) : fRunList(), fOCDBPath(ocdbPath), fShouldCorrectForPileUp(kTRUE), fCrossSectionUnit("UB")
170 // ctor using a set of run numbers
174 //_____________________________________________________________________________
175 AliAnalysisTriggerScalers::AliAnalysisTriggerScalers(const char* runlist, const char* ocdbPath) : fRunList(), fOCDBPath(ocdbPath), fShouldCorrectForPileUp(kTRUE), fCrossSectionUnit("UB")
177 // ctor from a run list (txt file)
181 //_____________________________________________________________________________
182 AliAnalysisTriggerScalers::~AliAnalysisTriggerScalers()
187 //______________________________________________________________________________
188 Bool_t AliAnalysisTriggerScalers::CheckRecord(const AliTriggerScalersRecord& record,
192 UInt_t timelapse) const
194 /// Check if a trigger scaler record should be kept for our
195 /// luminosity computations
197 const AliTriggerScalers* scaler = record.GetTriggerScalersForClass(index);
199 UInt_t a = scaler->GetLOCA() - refa;
200 UInt_t b = scaler->GetLOCB() - refb;
204 if ( b <= 2 || ( a <= 2 ) || timelapse <= 9 ) // empirical cuts
213 //______________________________________________________________________________
214 void AliAnalysisTriggerScalers::DrawFill(Int_t run1, Int_t run2, double ymin, double ymax, const char* label)
216 // Draw a yellow box to indicate a run range
218 AliDebugClass(1,Form("RUN1 %09d RUN2 %09d YMIN %e YMAX %e %s",
219 run1,run2,ymin,ymax,label));
220 TBox* b = new TBox(run1*1.0,ymin,run2*1.0,ymax);
223 TText* text = new TText((run1+run2)/2.0,ymin + (ymax-ymin)*0.85,label);
224 text->SetTextSize(0.025);
225 text->SetTextFont(42);
226 text->SetTextAlign(23);
227 text->SetTextAngle(45);
231 //_____________________________________________________________________________
232 void AliAnalysisTriggerScalers::DrawFills(Double_t ymin, Double_t ymax)
234 /// Draw the fill ranges corresponding to the list of runs
235 /// Note that this method will scan the OCDB to get the run -> fill number relationship,
236 /// so it's better in this case to use a local copy of the OCDB if you have one. Otherwise
239 std::map<int, std::pair<int,int> > fills;
241 GetFillBoundaries(fills);
243 std::map<int, std::pair<int,int> >::const_iterator it;
245 for ( it = fills.begin(); it != fills.end(); ++it )
247 const std::pair<int,int>& p = it->second;
249 fillnumber.Form("%d",it->first);
250 DrawFill(p.first,p.second,ymin,ymax,fillnumber.Data());
254 //_____________________________________________________________________________
255 void AliAnalysisTriggerScalers::GetCTPObjects(Int_t runNumber,
256 AliTriggerConfiguration*& tc,
257 AliTriggerRunScalers*& trs,
258 AliLHCData*& lhc) const
260 /// For a given run, get the CTP objects we need to do our work
262 /// FIXME: this is the method that should probably be virtualized so
263 /// we can get those objects from either OCDB or run based container in AODs/ESDs
265 tc = static_cast<AliTriggerConfiguration*>(GetOCDBObject("GRP/CTP/Config",runNumber));
266 trs = static_cast<AliTriggerRunScalers*>(GetOCDBObject("GRP/CTP/Scalers",runNumber));
267 lhc = static_cast<AliLHCData*>(GetOCDBObject("GRP/GRP/LHCData",runNumber));
270 //_____________________________________________________________________________
271 void AliAnalysisTriggerScalers::GetFillBoundaries(std::map<int, std::pair<int,int> >& fills)
273 /// Given a list of runs get the fills they are in
277 for ( std::vector<int>::const_iterator it = fRunList.begin(); it != fRunList.end(); ++it )
281 int fill = GetFillNumberFromRunNumber(runnumber);
283 if (fills.count(fill))
285 std::pair<int,int>& p = fills[fill];
286 p.first = TMath::Min(runnumber,p.first);
287 p.second = TMath::Max(runnumber,p.second);
291 fills[fill] = make_pair<int,int>(runnumber,runnumber);
297 //______________________________________________________________________________
298 Int_t AliAnalysisTriggerScalers::GetFillNumberFromRunNumber(Int_t runNumber)
300 /// Get the fill number of a run
302 AliLHCData* lhcData = static_cast<AliLHCData*>(GetOCDBObject("GRP/GRP/LHCData",runNumber));
305 Int_t fillNumber = lhcData->GetFillNumber();
306 if ( fillNumber == 0 && runNumber == 189694)
308 // manual hack because GRP info incorrect for this run ?
317 //______________________________________________________________________________
318 TObject* AliAnalysisTriggerScalers::GetOCDBObject(const char* path, Int_t runNumber) const
320 // Helper method to get an object from the OCDB
322 if ( !AliCDBManager::Instance()->IsDefaultStorageSet() )
324 AliCDBManager::Instance()->SetDefaultStorage(fOCDBPath.Data());
327 AliCDBManager::Instance()->SetRun(runNumber);
329 gErrorIgnoreLevel=kWarning; // to avoid all the TAlienFile::Open messages...
331 AliCDBEntry* e = AliCDBManager::Instance()->Get(path);
332 return e ? e->GetObject() : 0x0;
336 //______________________________________________________________________________
337 TString AliAnalysisTriggerScalers::GetLHCPeriodFromRunNumber(Int_t runNumber) const
339 /// Get a list of LHC periods from a list of run numbers
341 AliGRPObject* grp = static_cast<AliGRPObject*>(GetOCDBObject("GRP/GRP/Data",runNumber));
345 return grp->GetLHCPeriod();
348 //_____________________________________________________________________________
350 AliAnalysisTriggerScalers::GetLuminosityTriggerAndCrossSection(Int_t runNumber,
351 TString& lumiTriggerClassName,
352 Double_t& lumiTriggerCrossSection,
353 Double_t& lumiTriggerCrossSectionError) const
355 /// For one given run, get the trigger class name to be used as the luminometer,
356 /// and its corresponding cross-section
358 /// FIXME: all is hard-coded for now... (use OADB for this ?)
360 lumiTriggerClassName="";
361 lumiTriggerCrossSection=0.0;
362 lumiTriggerCrossSectionError=0.0; // FIXME
364 TString period = GetLHCPeriodFromRunNumber(runNumber);
366 if ( period.BeginsWith("LHC11") )
368 AliError("Not implemented yet");
370 else if ( period.BeginsWith("LHC12") )
372 if ( period == "LHC12h" || period == "LHC12i" )
374 // pp 8 TeV main-satellites
375 lumiTriggerClassName = "C0TVX-S-NOPF-ALLNOTRD";
376 lumiTriggerCrossSection = 28.0; // FIXME: not from a vdM yet
380 AliError("Not implemented yet");
383 else if ( period.BeginsWith("LHC13") )
385 if ( period == "LHC13g" )
388 lumiTriggerClassName = "C0TVX-B-NOPF-ALLNOTRD";
389 lumiTriggerCrossSection = 18.0; // FIXME: not from a vdM yet
393 // p-Pb or Pb-p 5.02 TeV
394 lumiTriggerClassName = "C0TVX-B-NOPF-ALLNOTRD";
395 lumiTriggerCrossSection = 0.755*2000; // FIXME: not from a vdM yet
400 AliError("Not implemented yet");
405 if ( fCrossSectionUnit=="PB" )
409 else if (fCrossSectionUnit=="NB")
413 else if ( fCrossSectionUnit=="UB" )
417 else if ( fCrossSectionUnit=="MB" )
422 lumiTriggerCrossSection *= csMult;
423 lumiTriggerCrossSectionError *= csMult;
426 //_____________________________________________________________________________
427 void AliAnalysisTriggerScalers::GetLHCPeriodBoundaries(std::map<std::string, std::pair<int,int> >& periods)
429 /// Given a list of runs get the fills they are in
433 for ( std::vector<int>::const_iterator it = fRunList.begin(); it != fRunList.end(); ++it )
437 std::string period = GetLHCPeriodFromRunNumber(runnumber).Data();
439 if (periods.count(period))
441 std::pair<int,int>& p = periods[period];
442 p.first = TMath::Min(runnumber,p.first);
443 p.second = TMath::Max(runnumber,p.second);
447 periods[period] = make_pair<int,int>(runnumber,runnumber);
453 //______________________________________________________________________________
454 Int_t AliAnalysisTriggerScalers::GetTriggerInput(Int_t runNumber, const char* inputname)
456 /// Get one trigger input for a given run
457 AliTriggerConfiguration* tc = static_cast<AliTriggerConfiguration*>(GetOCDBObject("GRP/CTP/Config",runNumber));
460 const TObjArray& inputs = tc->GetInputs();
462 AliTriggerInput* ti = static_cast<AliTriggerInput*>(inputs.FindObject(inputname));
466 return ti->GetSignature();
469 //______________________________________________________________________________
471 AliAnalysisTriggerScalers::GetPauseAndConfigCorrection(Int_t runNumber, const char* triggerClassName)
473 /// Tries to estimate the duration of the Pause and Configure phase(s) in a given run
474 /// Probably not really accurate for the moment though.
479 AliTriggerConfiguration* tc = static_cast<AliTriggerConfiguration*>(GetOCDBObject("GRP/CTP/Config",runNumber));
481 AliTriggerRunScalers* trs = static_cast<AliTriggerRunScalers*>(GetOCDBObject("GRP/CTP/Scalers",runNumber));
482 const TObjArray* scalers = trs->GetScalersRecords();
485 AliTriggerScalersRecord* record;
491 Int_t classIndex = tc->GetClassIndexFromName(triggerClassName);
493 while ( ( record = static_cast<AliTriggerScalersRecord*>(next()) ) )
495 const AliTriggerScalers* scaler = record->GetTriggerScalersForClass(classIndex);
497 const AliTimeStamp* ats = record->GetTimeStamp();
499 UInt_t seconds = ats->GetSeconds();// - TTimeStamp::GetZoneOffset();
501 TTimeStamp ts(seconds,ats->GetMicroSecs());
503 UInt_t l0b = scaler->GetLOCB() - refl0b;
504 UInt_t l0a = scaler->GetLOCA() - refl0a;
505 UInt_t timelapse = seconds - reft;
507 if ( l0b <=2 || timelapse < 9 ) continue;
510 refl0b = scaler->GetLOCB();
511 refl0a = scaler->GetLOCA();
527 return total > 0 ? nopac*1.0/total : 0.0;
530 //______________________________________________________________________________
531 AliAnalysisTriggerScalerItem*
532 AliAnalysisTriggerScalers::GetTriggerScaler(Int_t runNumber, const char* level, const char* triggerClassName)
534 // Get the scaler for a given trigger class and a given trigger level.
535 // Returned object must be deleted by the user.
537 AliTriggerConfiguration* tc = static_cast<AliTriggerConfiguration*>(GetOCDBObject("GRP/CTP/Config",runNumber));
538 AliTriggerRunScalers* trs = static_cast<AliTriggerRunScalers*>(GetOCDBObject("GRP/CTP/Scalers",runNumber));
539 AliGRPObject* grp = static_cast<AliGRPObject*>(GetOCDBObject("GRP/GRP/Data",runNumber));
541 if (!tc || !trs || !grp) return 0x0;
543 TString diCurrent(Form("L3:%5.0f;DIP:%5.0f [L3:%d;DIP:%d]",
544 grp->GetL3Current((AliGRPObject::Stats)0),
545 grp->GetDipoleCurrent((AliGRPObject::Stats)0),
546 grp->GetL3Polarity(),
547 grp->GetDipolePolarity()));
549 const TObjArray& trClasses = tc->GetClasses();
551 const TObjArray* scalers = trs->GetScalersRecords();
552 const AliTriggerScalersRecord* begin = (AliTriggerScalersRecord*)(scalers->First());
553 const AliTriggerScalersRecord* end = (AliTriggerScalersRecord*)(scalers->Last());
555 time_t duration = (end->GetTimeStamp()->GetBunchCross() - begin->GetTimeStamp()->GetBunchCross())*AliTimeStamp::fNanosecPerBC*1E-9;
557 AliTriggerClass* triggerClass = static_cast<AliTriggerClass*>(trClasses.FindObject(triggerClassName));
562 UChar_t index = GetIndex(triggerClass->GetMask());
566 const AliTriggerScalers* sbegin = begin->GetTriggerScalersForClass(index);
567 const AliTriggerScalers* send = end->GetTriggerScalersForClass(index);
569 TString swhat(level);
572 if ( swhat.BeginsWith("L0A") )
574 value = send->GetLOCA() - sbegin->GetLOCA();
576 else if ( swhat.BeginsWith("L0B") )
578 value = send->GetLOCB() - sbegin->GetLOCB();
580 else if ( swhat.BeginsWith("L1A") )
582 value = send->GetL1CA() - sbegin->GetL1CA();
584 else if ( swhat.BeginsWith("L1B") )
586 value = send->GetL1CB() - sbegin->GetL1CB();
588 else if ( swhat.BeginsWith("L2A") )
590 value = send->GetL2CA() - sbegin->GetL2CA();
592 else if ( swhat.BeginsWith("L2B") )
594 value = send->GetL2CB() - sbegin->GetL2CB();
598 // FIXME: get the downscaling factor here for all cases (only BC1 supported here so far)
599 if ( TString(triggerClassName).Contains("_B1") )
601 ds = GetTriggerInput(runNumber,"BC1");
605 swhat.ReplaceAll("RATE","");
607 return new AliAnalysisTriggerScalerItem(runNumber,swhat.Data(),diCurrent.Data(),triggerClass->GetName(),value,triggerClass->GetBCMask(),ds,duration);
610 //______________________________________________________________________________
612 AliAnalysisTriggerScalers::IntegratedLuminosityGraph(const char* triggerClassName, const char* triggerClassNameForPACEstimate)
614 // Get the integrated luminosity graph for one trigger
616 if ( fRunList.empty() )
618 AliError("Cannot work without a run list");
622 Double_t offset(0.0);
624 std::vector<double> vx;
625 std::vector<double> vy;
629 for ( std::vector<int>::size_type i = 0; i < fRunList.size(); ++i )
631 TGraph* g = IntegratedLuminosityGraph(fRunList[i],triggerClassName,triggerClassNameForPACEstimate);
635 for ( Int_t j = 0 ; j < g->GetN(); ++j )
639 vy.push_back(y+offset);
647 TGraph* g = new TGraph(vx.size(),&vx[0],&vy[0]);
654 //______________________________________________________________________________
656 AliAnalysisTriggerScalers::IntegratedLuminosityGraph(Int_t runNumber, const char* triggerClassName, const char* triggerClassNameForPACEstimate)
658 // Get the integrated luminosity graph for one trigger for the given run
660 // the triggerClassNameForPACEstimate is only used if triggerClassName is,
661 // for the given run, the same as the luminometer trigger class name,
662 // and should be the trigger class with the higher and non-zero L0A rate (except during PACs)
665 AliTriggerConfiguration* tc(0x0);
666 AliTriggerRunScalers* trs(0x0);
667 AliLHCData* lhc(0x0);
669 GetCTPObjects(runNumber,tc,trs,lhc);
671 if (!tc || !trs || !lhc)
676 const TObjArray* scalerRecords= trs->GetScalersRecords();
677 const TObjArray& triggerClasses = tc->GetClasses();
679 TString lumiTriggerClassName;
680 Double_t lumiTriggerCrossSection(0.0);
681 Double_t lumiTriggerCrossSectionError(0.0);
683 GetLuminosityTriggerAndCrossSection(runNumber,lumiTriggerClassName,lumiTriggerCrossSection,lumiTriggerCrossSectionError);
685 AliTriggerClass* lumiTriggerClass = static_cast<AliTriggerClass*>(triggerClasses.FindObject(lumiTriggerClassName));
687 if (!lumiTriggerClass)
689 AliError(Form("Could not find lumiTriggerClass=%s",lumiTriggerClassName.Data()));
693 if (lumiTriggerCrossSection==0.0)
695 AliError(Form("Could not get cross-section for lumiTriggerClass=%s",lumiTriggerClassName.Data()));
699 AliTriggerClass* triggerClass = static_cast<AliTriggerClass*>(triggerClasses.FindObject(triggerClassName));
703 AliError(Form("Could not find triggerClass=%s",triggerClassName));
707 Int_t nbcx = NumberOfInteractingBunches(*lhc);
709 if (nbcx <= 0 && ShouldCorrectForPileUp())
711 AliError("Could not get the number of bunches, so won't be able to correct for pile-up");
714 TIter nextScaler(scalerRecords);
715 AliTriggerScalersRecord* record;
718 UInt_t refl0b[] = { 0, 0 };
719 UInt_t refl0a[] = { 0, 0 };
721 UInt_t l0b[] = { 0, 0 };
722 UInt_t l0a[] = { 0, 0 };
726 UChar_t classIndices[2];
728 // class 0 = class for luminomiter
729 // class 1 = class for livetime estimate
730 // or for PAC estimate if triggerClassNameForPACEstimate
731 // is given and triggerClass is the same as the lumi class
733 classIndices[0] = GetIndex(lumiTriggerClass->GetMask());
734 classIndices[1] = GetIndex(triggerClass->GetMask());
736 Bool_t sameClass = ( classIndices[0] == classIndices[1] );
737 Bool_t pacRemoval(kFALSE);
739 if ( sameClass && strlen(triggerClassNameForPACEstimate) > 0 )
741 AliTriggerClass* triggerClassForPACEstimate = static_cast<AliTriggerClass*>(triggerClasses.FindObject(triggerClassNameForPACEstimate));
742 if (!triggerClassForPACEstimate)
744 AliError(Form("Could not find triggerClassForPACEstimate=%s. Will not correct for PAC durations",triggerClassNameForPACEstimate));
747 classIndices[1] = GetIndex(triggerClassForPACEstimate->GetMask());
748 sameClass = ( classIndices[0] == classIndices[1] );
755 ULong64_t counter(0);
757 std::vector<Double_t> vseconds;
758 std::vector<Double_t> vlivetime;
759 std::vector<Double_t> vlumi;
761 while ( ( record = static_cast<AliTriggerScalersRecord*>(nextScaler()) ) )
763 const AliTimeStamp* ats = record->GetTimeStamp();
765 UInt_t seconds = ats->GetSeconds();// - TTimeStamp::GetZoneOffset();
767 TTimeStamp ts(seconds,ats->GetMicroSecs());
769 UInt_t timelapse = seconds - reft;
771 Bool_t ok = sameClass || CheckRecord(*record,classIndices[1],refl0b[1],refl0a[1],timelapse);
773 if (ok) reft = seconds;
775 for ( Int_t i = 0; i < 2; ++i )
777 const AliTriggerScalers* scaler = record->GetTriggerScalersForClass(classIndices[i]);
781 l0b[i] = scaler->GetLOCB() - refl0b[i];
787 l0a[i] = scaler->GetLOCA() - refl0a[i];
789 refl0b[i] = scaler->GetLOCB();
790 refl0a[i] = scaler->GetLOCA();
799 if ( ShouldCorrectForPileUp() && nbcx > 0 && l0b[0] > 0 )
801 Double_t mu = Mu(l0b[0]/timelapse,nbcx);
802 Double_t factor = mu/(1-TMath::Exp(-mu));
807 vseconds.push_back(seconds);
811 if ( ok && !sameClass && !pacRemoval )
813 lt = (l0a[1]*1.0)/l0b[1];
816 counter += l0b[0]*lt;
818 vlumi.push_back( counter / lumiTriggerCrossSection );
822 if ( vseconds.empty() ) return 0x0;
824 TGraph* g = new TGraph(vseconds.size(),&vseconds[0],&vlumi[0]);
831 //______________________________________________________________________________
832 void AliAnalysisTriggerScalers::IntegratedLuminosity(const char* triggerList,
833 const char* lumiTrigger,
834 Double_t lumiCrossSection,
835 const char* csvOutputFile,
839 // Compute the luminosity for a set of triggers
841 // for T0 based lumi (end of pp 2012), use lumiTrigger = C0TVX-S-NOPF-ALLNOTRD and lumiCrossSection = 28 mb (and csUnit="nb")
843 // for T0 based lumi (pPb 2013), use lumiTrigger = C0TVX-B-NOPF-ALLNOTRD
844 // and lumiCrossSection = 0.755*2000 mb (and csUnit="mub")
846 // for T0 based lumi (pp 2.76 TeV 2013), use lumiTrigger = C0TVX-B-NOPF-ALLNOTRD
847 // and lumiCrossSection = 18 mb (and csUnit="nb")
851 std::map<std::string,float> lumiPerTrigger;
852 std::map<int, std::map<std::string,float> > lumiPerFillPerTrigger;
854 std::map<std::string,time_t> durationPerTrigger;
856 TString sTriggerList(triggerList);
858 if ( sTriggerList.Length()==0 )
860 if ( lumiCrossSection < 30 && lumiCrossSection > 20 )
862 // assuming it's pp 2012
863 sTriggerList = "CMUL8-S-NOPF-ALLNOTRD,CMUL7-S-NOPF-ALLNOTRD,CMUL8-S-NOPF-MUON,CMUL7-S-NOPF-MUON";
865 // for C0MUL-SC-NOPF-MUON must use C0TVX-SC as lumiTrigger (with same cross-section as C0TVX=28mb)
867 else if ( lumiCrossSection < 30 )
869 // assuming it's pp 2013
870 sTriggerList = "CMSL7-B-NOPF-MUON,CMSH7-B-NOPF-MUON,CMUL7-B-NOPF-MUON,CMSL7-B-NOPF-ALLNOTRD,CMSH7-B-NOPF-ALLNOTRD,CMUL7-B-NOPF-ALLNOTRD";
874 sTriggerList = "CMSL7-B-NOPF-MUON,CMSH7-B-NOPF-MUON,CMUL7-B-NOPF-MUON,CMSL7-B-NOPF-ALLNOTRD,CMSH7-B-NOPF-ALLNOTRD,CMUL7-B-NOPF-ALLNOTRD";
878 TObjArray* triggerArray = sTriggerList.Tokenize(",");
879 triggerArray->SetOwner(kTRUE);
881 std::ofstream* out(0x0);
883 if ( TString(csvOutputFile).Length() > 0 )
885 out = new std::ofstream(gSystem->ExpandPathName(csvOutputFile));
886 if (!out || out->bad())
893 TIter nextTrigger(triggerArray);
898 (*out) << "Fill number" << sep;
902 std::map<std::string,float>::const_iterator it;
904 while ( ( trigger = static_cast<TObjString*>(nextTrigger())))
906 lumiPerTrigger[trigger->String().Data()] = 0;
909 for ( it = lumiPerTrigger.begin(); it != lumiPerTrigger.end(); ++it )
911 (*out) << "lumi from " << it->first.c_str() << sep;
914 (*out) << "comments" << sep;
916 for ( it = lumiPerTrigger.begin(); it != lumiPerTrigger.end(); ++it )
918 (*out) << "recorded " << it->first.c_str() << " (integrated)" << sep;
921 (*out) << Form("LHC delivered (%s^-1) per fill ",csUnit)
922 << sep << Form("LHC delivered (%s^-1 integrated)",csUnit) << sep;
923 (*out) << "lumi tot muon" << sep << "eff (%)" << sep;
929 Int_t currentFillNumber(-1);
933 TString scsUnit(csUnit);
940 else if (scsUnit=="NB")
944 else if ( scsUnit=="UB" )
948 else if ( scsUnit=="MB" )
954 AliError(Form("Don't know how to deal with cross-section unit=%s",csUnit));
959 for ( std::vector<int>::size_type i = 0; i < fRunList.size(); ++i )
961 Int_t runNumber = fRunList[i];
962 Bool_t atLeastOneTriggerFound(kFALSE);
964 // find out which trigger classes to use
966 AliTriggerConfiguration* tc = static_cast<AliTriggerConfiguration*>(GetOCDBObject("GRP/CTP/Config",runNumber));
967 const TObjArray& trClasses = tc->GetClasses();
971 fillNumber = GetFillNumberFromRunNumber(runNumber);
972 if ( fillNumber == 0 )
974 AliError(Form("Got fillNumber = 0 for run %09d !",runNumber));
977 if ( fillNumber != currentFillNumber )
979 std::map<std::string,float>::const_iterator it;
981 for ( it = lumiPerTrigger.begin(); it != lumiPerTrigger.end(); ++it )
983 lumiPerFillPerTrigger[fillNumber][it->first.c_str()] = 0;
985 currentFillNumber = fillNumber;
989 AliTriggerRunScalers* trs = static_cast<AliTriggerRunScalers*>(GetOCDBObject("GRP/CTP/Scalers",runNumber));
990 const TObjArray* scalers = trs->GetScalersRecords();
991 const AliTriggerScalersRecord* begin = (AliTriggerScalersRecord*)(scalers->First());
992 const AliTriggerScalersRecord* end = (AliTriggerScalersRecord*)(scalers->Last());
994 time_t duration = (end->GetTimeStamp()->GetBunchCross() - begin->GetTimeStamp()->GetBunchCross())*AliTimeStamp::fNanosecPerBC*1E-9;
998 TString lumiTriggerClassName(lumiTrigger);
999 float lumiSigma = lumiCrossSection*csMult; //in csUnit
1000 AliAnalysisTriggerScalerItem* lumiB = GetTriggerScaler(runNumber,"L0B",lumiTriggerClassName.Data());
1004 AliError(Form("Did not find lumiTrigger %s for run %09d",lumiTriggerClassName.Data(),runNumber));
1008 Float_t pacCorrection(1.0);
1010 while ( ( trigger = static_cast<TObjString*>(nextTrigger()) ) )
1012 TString muTriggerClassName(trigger->String());
1015 if ( !trClasses.FindObject(muTriggerClassName.Data() ) )
1020 if ( muTriggerClassName.Contains("CMUL8") ) ++n;
1021 if ( muTriggerClassName.Contains("CMUL7") ) ++n;
1025 AliError(Form("More than 1 relevant trigger class found for run %09d ! Check that !",runNumber));
1030 AliAnalysisTriggerScalerItem* muonA = GetTriggerScaler(runNumber,"L0A",muTriggerClassName.Data());
1031 AliAnalysisTriggerScalerItem* muonB = GetTriggerScaler(runNumber,"L0B",muTriggerClassName.Data());
1033 if (!muonA || !muonB) continue;
1035 atLeastOneTriggerFound = kTRUE;
1039 if ( muonB->ValueCorrectedForDownscale() > 0 )
1041 ratio = muonA->ValueCorrectedForDownscale()/muonB->ValueCorrectedForDownscale();
1044 ratio *= lumiB->ValueCorrectedForDownscale()/lumiSigma;
1046 if ( muTriggerClassName.BeginsWith("CMUL"))
1051 if ( muTriggerClassName.Contains("CMSH") )
1053 pacCorrection = GetPauseAndConfigCorrection(runNumber,muTriggerClassName.Data());
1056 lumiPerTrigger[muTriggerClassName.Data()] += ratio;
1057 durationPerTrigger[muTriggerClassName.Data()] += duration;
1058 lumiPerFillPerTrigger[currentFillNumber][muTriggerClassName.Data()] += ratio;
1062 lumiPerTrigger[lumiTriggerClassName.Data()] += lumiB->ValueCorrectedForDownscale()/lumiSigma;
1063 durationPerTrigger[lumiTriggerClassName.Data()] += duration;
1064 lumiPerFillPerTrigger[currentFillNumber][lumiTriggerClassName.Data()] += lumiB->ValueCorrectedForDownscale()/lumiSigma;
1066 TString lumiPACCorrected(Form("%s(-PAC)",lumiTriggerClassName.Data()));
1068 lumiPerTrigger[lumiPACCorrected.Data()] += pacCorrection*lumiB->ValueCorrectedForDownscale()/lumiSigma;
1069 durationPerTrigger[lumiPACCorrected.Data()] += duration;
1070 lumiPerFillPerTrigger[currentFillNumber][lumiPACCorrected.Data()] += pacCorrection*lumiB->ValueCorrectedForDownscale()/lumiSigma;
1072 if (!atLeastOneTriggerFound && sTriggerList.Contains("CMUL") )
1074 AliError(Form("Found no relevant trigger for run %09d",runNumber));
1078 AliInfo(Form("Integrated lumi %7.4f %s^-1",intLumi,csUnit));
1080 std::map<std::string,float>::const_iterator it;
1082 for ( it = lumiPerTrigger.begin(); it != lumiPerTrigger.end(); ++it )
1084 AliInfo(Form("Trigger %30s Lumi %10.4f %s^-1 duration %10ld s",it->first.c_str(),it->second,csUnit,durationPerTrigger[it->first]));
1089 std::map<int, std::map<std::string, float> >::const_iterator fit;
1091 lumiPerTrigger.clear();
1093 for ( fit = lumiPerFillPerTrigger.begin(); fit != lumiPerFillPerTrigger.end(); ++fit )
1095 int fill = fit->first;
1096 std::map<std::string,float>::const_iterator tit;
1097 (*out) << fill << sep;
1099 for ( tit = fit->second.begin(); tit != fit->second.end(); ++tit )
1101 (*out) << Form("%e",tit->second) << sep;
1104 (*out) << sep; // comment (empty)
1106 for ( tit = fit->second.begin(); tit != fit->second.end(); ++tit )
1108 lumiPerTrigger[tit->first] += tit->second;
1110 (*out) << Form("%e",lumiPerTrigger[tit->first]) << sep;
1112 (*out) << sep << "0" << sep << "0" << sep << "0" << sep << "0" << sep << std::endl; // LHC per fill, LHC integrated, lumi tot muon , efficiency
1116 delete triggerArray;
1119 //______________________________________________________________________________
1120 TGraph* AliAnalysisTriggerScalers::MakeGraph(const std::vector<int>& vx,
1121 const std::vector<int>& vex,
1122 const std::vector<double>& vy,
1123 const std::vector<double>& vey)
1125 /// Build a graph from a set of stl vectors
1127 if ( ! ( vx.size() == vex.size() && vx.size() == vy.size() && vx.size() == vey.size() ) )
1129 std::cerr << "incompatible sizes" << std::endl;
1133 double* x = new double[vx.size()];
1134 double* ex = new double[vx.size()];
1135 double* y = new double[vx.size()];
1136 double* ey = new double[vx.size()];
1138 for ( size_t i = 0; i < vx.size(); ++i )
1146 TGraph* g = new TGraphErrors(vx.size(),x,y,ex,ey);
1158 //______________________________________________________________________________
1159 Int_t AliAnalysisTriggerScalers::NumberOfInteractingBunches(const AliLHCData& lhc) const
1161 /// Extract the number of colliding bunches from the LHC data
1163 Int_t numberOfInteractingBunches(0);
1164 Int_t numberOfInteractingBunchesMeasured(0);
1170 AliLHCDipValI* val = lhc.GetBunchConfigDeclared(beam1,0);
1172 for ( Int_t i = 0; i < val->GetSizeTotal(); ++i )
1174 if ( val->GetValue(i) < 0 ) ++numberOfInteractingBunches;
1177 AliLHCDipValI* valm = lhc.GetBunchConfigMeasured(beam1,0);
1179 for ( Int_t i = 0; i < valm->GetSizeTotal(); ++i )
1181 if ( valm->GetValue(i) < 0 ) ++numberOfInteractingBunchesMeasured;
1184 valm = lhc.GetBunchConfigMeasured(beam2,0);
1186 for ( Int_t i = 0; i < valm->GetSizeTotal(); ++i )
1188 if ( valm->GetValue(i) <= 0 ) ++nIBM2;
1191 if ( numberOfInteractingBunches != numberOfInteractingBunchesMeasured ||
1192 numberOfInteractingBunches != nIBM2 )
1194 AliWarning(Form("Got some different number of interacting bunches here ! NumberOfInteractingBunches=%3d NumberOfInteractingBunchesMeasured=%3d NIBM2=%3d",
1195 numberOfInteractingBunches,numberOfInteractingBunchesMeasured,nIBM2));
1198 return numberOfInteractingBunches;
1202 //______________________________________________________________________________
1203 TGraph* AliAnalysisTriggerScalers::PlotTrigger(const char* triggerClassName,
1206 // Plot one of the scalers (L0A,L0B,L0AOVERB,etc...) of a given triggerClass
1207 // Get one value per run.
1209 std::vector<Float_t> x;
1210 std::vector<Float_t> y;
1215 for ( std::vector<int>::size_type i = 0; i < fRunList.size(); ++i )
1217 Int_t runNumber = fRunList[i];
1219 AliAnalysisTriggerScalerItem* s = GetTriggerScaler(runNumber,what,triggerClassName);
1223 x.push_back(runNumber);
1225 Double_t value = s->ValueCorrectedForDownscale();
1227 if ( TString(what).Contains("RATE") )
1238 if ( fRunList.size() ) {
1239 mean /= fRunList.size();
1242 AliInfo(Form("Integral %e mean %e",integral,mean));
1244 return new TGraph(x.size(),&x[0],&y[0]);
1247 //______________________________________________________________________________
1248 TGraph* AliAnalysisTriggerScalers::PlotTriggerEvolution(const char* triggerClassName,
1254 /// Make a graph of "what" for a given trigger class.
1258 /// - L0B (level 0 before any veto applied)
1259 /// - L0A (level 0 after all vetoes)
1260 /// - L0AOVERB (L0A/L0B)
1261 /// - mu ( = -TMath::Log( 1 - P(0) ) where P(0) is the proba to have zero collisions in a bunch crossing)
1262 /// - pileupfactor = mu/(1-exp(-mu)) : the factor to apply to correct the cint1b count rate
1263 /// - vsnb = L0B/(NumberOfInteractingBunches*11245)
1264 /// - NBCX = NumberOfInteractingBunches
1266 TString swhat(what);
1269 if ( swhat.Contains(";"))
1271 swhat.ReplaceAll(";",",");
1272 AliWarningClass("; is not a valid separator, replaced it with ,");
1278 TObjArray* a = swhat.Tokenize(",");
1279 if (a->GetEntries()>1)
1283 TObjString* str(0x0);
1284 Double_t ymin(TMath::Limits<Double_t>::Max());
1286 Double_t xmin(TMath::Limits<Double_t>::Max());
1290 while ( ( str = static_cast<TObjString*>(next())) )
1292 g = PlotTriggerEvolution(triggerClassName,str->String().Data(),false,0x0,removeZeros);
1294 for ( Int_t i = 0; i < g->GetN(); ++i )
1296 ymin = TMath::Min(ymin,g->GetY()[i]);
1297 ymax = TMath::Max(ymax,g->GetY()[i]);
1299 xmin = TMath::Min(xmin,g->GetX()[0]);
1300 xmax = TMath::Max(xmax,g->GetX()[g->GetN()-1]);
1303 gStyle->SetOptTitle(0);
1305 AliInfoClass(Form("x %e ; %e y %e ; %e",xmin,xmax,ymin,ymax));
1306 TH2* h = new TH2F("h",triggerClassName,100,xmin,xmax,100,ymin,ymax);
1307 h->SetStats(kFALSE);
1308 h->GetXaxis()->SetTimeDisplay(1);
1309 h->GetXaxis()->SetTimeFormat("%d/%m %H:%M");
1310 h->GetXaxis()->SetTimeOffset(0,"gmt");
1311 h->GetXaxis()->SetNdivisions(505);
1314 TIter nextGraph(&graphs);
1316 while ( ( g = static_cast<TGraph*>(nextGraph())) )
1318 AliInfoClass(g->GetTitle());
1320 g->SetLineColor(color);
1321 g->SetMarkerColor(color);
1322 g->SetMarkerStyle(marker);
1331 std::vector<int> vx;
1332 std::vector<int> vex;
1333 std::vector<double> vy;
1334 std::vector<double> vey;
1340 for ( std::vector<int>::size_type iRun = 0; iRun < fRunList.size(); ++iRun )
1342 Int_t runNumber = fRunList[iRun];
1344 AliTriggerConfiguration* tc(0x0);
1345 AliTriggerRunScalers* trs(0x0);
1346 AliLHCData* lhc(0x0);
1348 GetCTPObjects(runNumber,tc,trs,lhc);
1350 Int_t numberOfInteractingBunches = NumberOfInteractingBunches(*lhc);
1352 const TObjArray* scalers = trs->GetScalersRecords();
1354 const TObjArray& trClasses = tc->GetClasses();
1355 TIter next(&trClasses);
1356 AliTriggerClass* triggerClass;
1358 while ( ( triggerClass = static_cast<AliTriggerClass*>(next()) ) )
1360 UChar_t index = GetIndex(triggerClass->GetMask());
1362 if ( !TString(triggerClass->GetName()).Contains(triggerClassName) ) continue;
1364 TIter nextScaler(scalers);
1365 AliTriggerScalersRecord* record;
1373 Bool_t first(kTRUE);
1375 while ( ( record = static_cast<AliTriggerScalersRecord*>(nextScaler()) ) )
1377 const AliTriggerScalers* scaler = record->GetTriggerScalersForClass(index);
1380 const AliTimeStamp* ats = record->GetTimeStamp();
1382 UInt_t seconds = ats->GetSeconds();// - TTimeStamp::GetZoneOffset();
1384 TTimeStamp ts(seconds,ats->GetMicroSecs());
1386 UInt_t l0b = scaler->GetLOCB() - refl0b;
1387 UInt_t l0a = scaler->GetLOCA() - refl0a;
1388 UInt_t l1b = scaler->GetL1CB() - refl1b;
1389 UInt_t l1a = scaler->GetL1CA() - refl1a;
1390 UInt_t l2b = scaler->GetL2CB() - refl2b;
1391 UInt_t l2a = scaler->GetL2CA() - refl2a;
1392 UInt_t timelapse = seconds - reft;
1394 if ( removeZeros && ( l0b <= 2 || ( l0a <= 2 && l0a != 0 ) || timelapse <= 9 ) )
1396 AliInfo(Form("Skipping point for %s l0b %d l0a %d timelapse %d ts %s",
1397 triggerClassName,l0b,l0a,timelapse,ts.AsString()));
1402 refl0b = scaler->GetLOCB();
1403 refl0a = scaler->GetLOCA();
1404 refl1b = scaler->GetL1CB();
1405 refl1a = scaler->GetL1CA();
1406 refl2b = scaler->GetL2CB();
1407 refl2a = scaler->GetL2CA();
1418 if (swhat.Contains("L0AOVERB") )
1420 value = l0a*1.0/l0b;
1423 error = value*TMath::Sqrt(1.0/l0b+1.0/l0a);
1430 else if ( swhat.Contains("L0B") )
1434 else if (swhat.Contains("L0A") )
1438 else if ( swhat.Contains("L1B") )
1442 else if (swhat.Contains("L1A") )
1446 else if ( swhat.Contains("L2B") )
1450 else if (swhat.Contains("L2A") )
1454 else if ( swhat.Contains("MU") )
1456 value = Mu(l0b/timelapse,numberOfInteractingBunches);
1457 error = 0.0; // FIXME
1459 else if ( swhat.Contains("PILEUPFACTOR") )
1461 Double_t mu = Mu(l0b/timelapse,numberOfInteractingBunches);
1462 value = mu/(1-TMath::Exp(-mu));
1463 error = -1.0; // FIXME
1465 else if ( swhat.Contains("VSNB") )
1467 value = l0b/(11245.0*numberOfInteractingBunches);
1468 error = -1.0; // FIXME
1470 else if ( swhat.Contains("NBCX"))
1472 value = numberOfInteractingBunches;
1473 error = 1.0; // FIXME
1478 AliInfo(Form("Unknown what %s",what));
1481 if ( ! swhat.Contains("OVER") && ! swhat.Contains("RATIO") &&
1482 ! swhat.Contains("MU") && ! swhat.Contains("PILEUPFACTOR") &&
1483 ! swhat.Contains("RAW") & ! swhat.Contains("NBCX") )
1488 if ( !TMath::Finite(value) ) continue;
1490 if ( !swhat.Contains("L0AOVERB") && error > 0 )
1492 error = value >0 ? 1.0/TMath::Sqrt(value) : 0.0;
1495 if (removeZeros && TMath::Abs(value) < 1E-6 )
1505 vx.push_back(seconds);
1508 vy.push_back(value);
1510 vey.push_back(error);
1518 if (mean && nvalues)
1523 if ( vx.empty() ) return 0;
1525 TGraph* g = MakeGraph(vx,vex,vy,vey);
1526 TString title(Form("TriggerEvolution%s-%s",triggerClassName,swhat.Data()));
1528 g->SetName(title.Data());
1529 g->SetTitle(title.Data());
1531 g->GetYaxis()->SetTitle(title.Data());
1535 TCanvas* c = NewCanvas(g->GetName());
1536 g->SetLineColor(color);
1537 g->SetMarkerColor(color);
1538 g->SetMarkerStyle(marker);
1540 c->SaveAs(Form("%s.png",title.Data()));
1547 //______________________________________________________________________________
1548 TGraph* AliAnalysisTriggerScalers::PlotTriggerRatio(const char* triggerClassName1,
1550 const char* triggerClassName2,
1553 // Plot the ratio of two scalers.
1554 // Get one value per run.
1556 std::vector<Float_t> x;
1557 std::vector<Float_t> y;
1559 for ( std::vector<int>::size_type i = 0; i < fRunList.size(); ++i )
1561 Int_t runNumber = fRunList[i];
1563 AliAnalysisTriggerScalerItem* s1 = GetTriggerScaler(runNumber,what1,triggerClassName1);
1564 AliAnalysisTriggerScalerItem* s2 = GetTriggerScaler(runNumber,what2,triggerClassName2);
1566 if (!s1 || !s2) continue;
1568 x.push_back(runNumber);
1572 if ( s2->ValueCorrectedForDownscale() > 0 )
1574 ratio = s1->ValueCorrectedForDownscale()/s2->ValueCorrectedForDownscale();
1579 AliDebug(1,Form("RUN %09d %20s (%s) %12llu (%5d) %20s (%s) %12llu (%5d) R %7.2f",
1581 triggerClassName1,what1,
1582 s1->Value(),s1->DownscalingFactor(),
1583 triggerClassName2,what2,
1584 s2->Value(),s2->DownscalingFactor(),
1588 return new TGraph(x.size(),&x[0],&y[0]);
1591 //______________________________________________________________________________
1592 TGraph* AliAnalysisTriggerScalers::PlotTriggerRatioEvolution(const char* triggerClassName1,
1594 const char* triggerClassName2,
1597 /// Plots the evolution of one trigger ratio
1599 Bool_t draw(kFALSE);
1600 Bool_t removeZeros(kFALSE);
1602 TGraph* g1 = PlotTriggerEvolution(triggerClassName1,what1,draw,0x0,removeZeros);
1603 TGraph* g2 = PlotTriggerEvolution(triggerClassName2,what2,draw,0x0,removeZeros);
1605 if (!g1 || !g2) return 0x0;
1607 if ( g1->GetN() != g2->GetN() )
1609 AliError("Oups. Did not get the same number of points for the 2 graphs ?!");
1615 Double_t x1err,x2err;
1616 Double_t y1err,y2err;
1618 TGraphErrors* g = new TGraphErrors(g1->GetN());
1621 for ( Int_t i = 0; i < g1->GetN(); ++i )
1623 g1->GetPoint(i,x1,y1);
1624 g2->GetPoint(i,x2,y2);
1626 x1err = g1->GetErrorX(i);
1627 x2err = g2->GetErrorX(i);
1629 y1err = g1->GetErrorY(i);
1630 y2err = g2->GetErrorY(i);
1634 AliError(Form("Points at %d don't have the same x ! ",i));
1640 if ( TMath::Abs(y2) > 1E-12 )
1647 if ( TMath::Abs(x1) > 1E-12)
1649 yerr += TMath::Sqrt(x1err*x1err/(x1*x1));
1652 if ( TMath::Abs(x2) > 1E-12)
1654 yerr += TMath::Sqrt(x2err*x2err/(x2*x2));
1659 g->SetPoint(j,x1,y);
1660 g->SetPointError(j,x1err,yerr);
1673 //______________________________________________________________________________
1674 void AliAnalysisTriggerScalers::Print(Option_t* /* opt */) const
1676 /// print our runlist
1677 AliAnalysisTriggerScalers::PrintIntegers(fRunList,',');
1680 //______________________________________________________________________________
1681 void AliAnalysisTriggerScalers::PrintIntegers(const std::vector<int>& integers,
1685 /// print a list of integers
1686 for ( std::vector<int>::size_type i = 0; i < integers.size(); ++i )
1688 out << integers[i] << sep;
1693 //______________________________________________________________________________
1694 void AliAnalysisTriggerScalers::ReadIntegers(const char* filename,
1695 std::vector<int>& integers,
1698 /// Read integers from filename, where integers are either
1699 /// separated by "," or by return carriage
1700 std::ifstream in(gSystem->ExpandPathName(filename));
1703 std::set<int> runset;
1707 for ( std::vector<int>::size_type j = 0; j < integers.size(); ++ j )
1709 runset.insert(integers[j]);
1715 in.getline(line,10000,'\n');
1717 TString sline(line);
1719 if (sline.Contains(","))
1721 TObjArray* a = sline.Tokenize(",");
1724 while ( ( s = static_cast<TObjString*>(next()) ) )
1726 runset.insert(s->String().Atoi());
1732 runset.insert(sline.Atoi());
1742 for ( std::set<int>::const_iterator it = runset.begin(); it != runset.end(); ++it )
1744 integers.push_back((*it));
1747 std::sort(integers.begin(),integers.end());
1751 //______________________________________________________________________________
1752 void AliAnalysisTriggerScalers::SetRunList(Int_t runNumber)
1754 // Make the runlist be a single run
1756 fRunList.push_back(runNumber);
1759 //______________________________________________________________________________
1760 void AliAnalysisTriggerScalers::SetRunList(const std::vector<int>& runs)
1762 // Make the runlist be a single run
1766 //______________________________________________________________________________
1767 void AliAnalysisTriggerScalers::SetRunList(const std::set<int>& runs)
1769 // Make the runlist be a single run
1771 for ( std::set<int>::const_iterator it = runs.begin(); it != runs.end(); ++it )
1773 fRunList.push_back(*it);
1777 //______________________________________________________________________________
1778 void AliAnalysisTriggerScalers::SetRunList(const char* runlist)
1780 // Read the runlist from an ASCII file or a comma separated list
1781 // or a space separated list
1785 if ( TString(runlist).Contains(",") || TString(runlist).Contains(" ") )
1787 TObjArray* runs = 0x0;
1788 if ( TString(runlist).Contains(",") )
1790 runs = TString(runlist).Tokenize(",");
1794 runs = TString(runlist).Tokenize(" ");
1798 std::set<int> runset;
1800 while ( ( s = static_cast<TObjString*>(next()) ) )
1802 runset.insert(s->String().Atoi());
1805 for ( std::set<int>::const_iterator it = runset.begin(); it != runset.end(); ++it )
1807 fRunList.push_back((*it));
1810 std::sort(fRunList.begin(),fRunList.end());
1816 ReadIntegers(runlist,fRunList);
1820 //______________________________________________________________________________
1821 //______________________________________________________________________________
1822 //______________________________________________________________________________
1823 //______________________________________________________________________________
1824 //______________________________________________________________________________
1825 //______________________________________________________________________________
1826 //______________________________________________________________________________
1828 ClassImp(AliAnalysisTriggerScalerItem)
1830 //______________________________________________________________________________
1832 AliAnalysisTriggerScalerItem::BCMaskName() const
1834 // return bc mask name
1835 if ( BCMask() ) return BCMask()->GetName(); else return "";
1838 //______________________________________________________________________________
1839 Int_t AliAnalysisTriggerScalerItem::Compare(const TObject* obj) const
1841 /// compare two scaler items (by means of their run number only)
1842 const AliAnalysisTriggerScalerItem* s = static_cast<const AliAnalysisTriggerScalerItem*>(obj);
1844 if ( s->RunNumber() < RunNumber() )
1848 else if ( s->RunNumber() > RunNumber() )
1855 //______________________________________________________________________________
1856 void AliAnalysisTriggerScalerItem::Print(Option_t* opt) const
1861 if ( RunNumber() > 0 )
1865 if ( sopt.Contains("FULL") )
1870 std::cout << Form("RUN %6d CLASS %24s (%5s %4d BCs) SCALER %s %12llu DS %6d DURATION %ld",
1871 RunNumber(),TriggerClassName(),
1873 BCMask() ? BCMask()->GetNUnmaskedBCs() : 0,
1875 Value(),DownscalingFactor(),Duration()) << std::endl;
1880 std::cout << Form("CLASS %24s SCALER %20llu NRUNS %d",
1881 TriggerClassName(),Value(),NofRuns()) << std::endl;