Belle II Software  release-06-00-14
eclMuMuEAlgorithm Class Reference

Calibrate ecl crystals using muon pair events. More...

#include <eclMuMuEAlgorithm.h>

Inheritance diagram for eclMuMuEAlgorithm:
Collaboration diagram for eclMuMuEAlgorithm:

Public Types

enum  EResult {
  c_OK ,
  c_Iterate ,
  c_NotEnoughData ,
  c_Failure ,
  c_Undefined
}
 The result of calibration. More...
 

Public Member Functions

 eclMuMuEAlgorithm ()
 ..Constructor
 
virtual ~eclMuMuEAlgorithm ()
 ..Destructor
 
std::string getPrefix () const
 Get the prefix used for getting calibration data.
 
bool checkPyExpRun (PyObject *pyObj)
 Checks that a PyObject can be successfully converted to an ExpRun type. More...
 
Calibration::ExpRun convertPyExpRun (PyObject *pyObj)
 Performs the conversion of PyObject to ExpRun. More...
 
std::string getCollectorName () const
 Alias for prefix. More...
 
void setPrefix (const std::string &prefix)
 Set the prefix used to identify datastore objects.
 
void setInputFileNames (PyObject *inputFileNames)
 Set the input file names used for this algorithm from a Python list. More...
 
PyObject * getInputFileNames ()
 Get the input file names used for this algorithm and pass them out as a Python list of unicode strings.
 
std::vector< Calibration::ExpRun > getRunListFromAllData () const
 Get the complete list of runs from inspection of collected data.
 
RunRange getRunRangeFromAllData () const
 Get the complete RunRange from inspection of collected data.
 
IntervalOfValidity getIovFromAllData () const
 Get the complete IoV from inspection of collected data.
 
void fillRunToInputFilesMap ()
 Fill the mapping of ExpRun -> Files.
 
std::string getGranularity () const
 Get the granularity of collected data.
 
EResult execute (std::vector< Calibration::ExpRun > runs={}, int iteration=0, IntervalOfValidity iov=IntervalOfValidity())
 Runs calibration over vector of runs for a given iteration. More...
 
EResult execute (PyObject *runs, int iteration=0, IntervalOfValidity iov=IntervalOfValidity())
 Runs calibration over Python list of runs. Converts to C++ and then calls the other execute() function.
 
std::list< Database::DBImportQuery > & getPayloads ()
 Get constants (in TObjects) for database update from last execution.
 
std::list< Database::DBImportQuerygetPayloadValues ()
 Get constants (in TObjects) for database update from last execution but passed by VALUE.
 
bool commit ()
 Submit constants from last calibration into database.
 
bool commit (std::list< Database::DBImportQuery > payloads)
 Submit constants from a (potentially previous) set of payloads.
 
const std::string & getDescription () const
 Get the description of the algoithm (set by developers in constructor)
 
bool loadInputJson (const std::string &jsonString)
 Load the m_inputJson variable from a string (useful from Python interface). The rturn bool indicates success or failure.
 
const std::string dumpOutputJson () const
 Dump the JSON string of the output JSON object.
 
const std::vector< Calibration::ExpRun > findPayloadBoundaries (std::vector< Calibration::ExpRun > runs, int iteration=0)
 Used to discover the ExpRun boundaries that you want the Python CAF to execute on. This is optional and only used in some.
 
template<>
std::shared_ptr< TTree > getObjectPtr (const std::string &name, const std::vector< Calibration::ExpRun > &requestedRuns)
 Specialization of getObjectPtr<TTree>.
 

Public Attributes

int cellIDLo
 ..Parameters to control Novosibirsk fit to energy deposited in each crystal by mu+mu- events More...
 
int cellIDHi
 Last cellID to be fit.
 
int minEntries
 All crystals to be fit must have at least minEntries events in the fit range.
 
int nToRebin
 If fewer entries than this, rebin and fix eta parameter.
 
double tRatioMin
 entries/peak at low edge of fit must be greater than this

 
double tRatioMax
 entries/peak at high edge of fit must be greater than this
 
double lowerEdgeThresh
 Lower edge is where the fit = lowerEdgeThresh * peak value.
 
bool performFits
 if false, input histograms are copied to output, but no fits are done.
 
bool findExpValues
 if true, fits are used to find expected energy deposit for each crystal instead of the calibration constant
 
int storeConst
 controls which values are written to the database. More...
 

Protected Member Functions

virtual EResult calibrate () override
 ..Run algorithm on events More...
 
void setInputFileNames (std::vector< std::string > inputFileNames)
 Set the input file names used for this algorithm. More...
 
virtual bool isBoundaryRequired (const Calibration::ExpRun &)
 Given the current collector data, make a decision about whether or not this run should be the start of a payload boundary.
 
virtual void boundaryFindingSetup (std::vector< Calibration::ExpRun >, int)
 If you need to make some changes to your algorithm class before 'findPayloadBoundaries' is run, make them in this function.
 
virtual void boundaryFindingTearDown ()
 Put your algorithm back into a state ready for normal execution if you need to.
 
const std::vector< Calibration::ExpRun > & getRunList () const
 Get the list of runs for which calibration is called.
 
int getIteration () const
 Get current iteration.
 
std::vector< std::string > getVecInputFileNames () const
 Get the input file names used for this algorithm as a STL vector.
 
template<class T >
std::shared_ptr< T > getObjectPtr (const std::string &name, const std::vector< Calibration::ExpRun > &requestedRuns)
 Get calibration data object by name and list of runs, the Merge function will be called to generate the overall object.
 
template<class T >
std::shared_ptr< T > getObjectPtr (std::string name)
 Get calibration data object (for all runs the calibration is requested for) This function will only work during or after execute() has been called once.
 
template<>
shared_ptr< TTree > getObjectPtr (const string &name, const vector< ExpRun > &requestedRuns)
 We cheekily cast the TChain to TTree for the returned pointer so that the user never knows Hopefully this doesn't cause issues if people do low level stuff to the tree...
 
std::string getGranularityFromData () const
 Get the granularity of collected data.
 
void saveCalibration (TClonesArray *data, const std::string &name)
 Store DBArray payload with given name with default IOV.
 
void saveCalibration (TClonesArray *data, const std::string &name, const IntervalOfValidity &iov)
 Store DBArray with given name and custom IOV.
 
void saveCalibration (TObject *data)
 Store DB payload with default name and default IOV.
 
void saveCalibration (TObject *data, const IntervalOfValidity &iov)
 Store DB payload with default name and custom IOV.
 
void saveCalibration (TObject *data, const std::string &name)
 Store DB payload with given name with default IOV.
 
void saveCalibration (TObject *data, const std::string &name, const IntervalOfValidity &iov)
 Store DB payload with given name and custom IOV.
 
void updateDBObjPtrs (const unsigned int event, const int run, const int experiment)
 Updates any DBObjPtrs by calling update(event) for DBStore.
 
void setDescription (const std::string &description)
 Set algorithm description (in constructor)
 
void clearCalibrationData ()
 Clear calibration data.
 
Calibration::ExpRun getAllGranularityExpRun () const
 Returns the Exp,Run pair that means 'Everything'. Currently unused.
 
void resetInputJson ()
 Clears the m_inputJson member variable.
 
void resetOutputJson ()
 Clears the m_outputJson member variable.
 
template<class T >
void setOutputJsonValue (const std::string &key, const T &value)
 Set a key:value pair for the outputJson object, expected to used interally during calibrate()
 
template<class T >
const T getOutputJsonValue (const std::string &key) const
 Get a value using a key from the JSON output object, not sure why you would want to do this.
 
template<class T >
const T getInputJsonValue (const std::string &key) const
 Get an input JSON value using a key. The normal exceptions are raised when the key doesn't exist.
 
const nlohmann::json & getInputJsonObject () const
 Get the entire top level JSON object. We explicitly say this must be of object type so that we might pick.
 
bool inputJsonKeyExists (const std::string &key) const
 Test for a key in the input JSON object.
 

Protected Attributes

std::vector< Calibration::ExpRun > m_boundaries
 

Private Member Functions

std::string getExpRunString (Calibration::ExpRun &expRun) const
 Gets the "exp.run" string repr. of (exp,run)
 
std::string getFullObjectPath (const std::string &name, Calibration::ExpRun expRun) const
 constructs the full TDirectory + Key name of an object in a TFile based on its name and exprun
 

Private Attributes

int fitOK = 16
 fit is OK
 
int iterations = 8
 fit reached max number of iterations, but is useable
 
int noLowerEdge = 5
 could not determine lower edge of fit
 
int atLimit = 4
 a parameter is at the limit; fit not useable
 
int poorFit = 3
 low chi square; fit not useable
 
int noPeak = 2
 Novosibirsk component of fit is negligible; fit not useable.
 
int notFit = -1
 no fit performed
 
std::vector< std::string > m_inputFileNames
 List of input files to the Algorithm, will initially be user defined but then gets the wildcards expanded during execute()
 
std::map< Calibration::ExpRun, std::vector< std::string > > m_runsToInputFiles
 Map of Runs to input files. Gets filled when you call getRunRangeFromAllData, gets cleared when setting input files again.
 
std::string m_granularityOfData
 Granularity of input data. This only changes when the input files change so it isn't specific to an execution.
 
ExecutionData m_data
 Data specific to a SINGLE execution of the algorithm. Gets reset at the beginning of execution.
 
std::string m_description {""}
 Description of the algorithm.
 
std::string m_prefix {""}
 The name of the TDirectory the collector objects are contained within.
 
nlohmann::json m_jsonExecutionInput = nlohmann::json::object()
 Optional input JSON object used to make decisions about how to execute the algorithm code.
 
nlohmann::json m_jsonExecutionOutput = nlohmann::json::object()
 Optional output JSON object that can be set during the execution by the underlying algorithm code.
 

Static Private Attributes

static const Calibration::ExpRun m_allExpRun = make_pair(-1, -1)
 

Detailed Description

Calibrate ecl crystals using muon pair events.

Definition at line 21 of file eclMuMuEAlgorithm.h.

Member Enumeration Documentation

◆ EResult

enum EResult
inherited

The result of calibration.

Enumerator
c_OK 

Finished successfuly =0 in Python.

c_Iterate 

Needs iteration =1 in Python.

c_NotEnoughData 

Needs more data =2 in Python.

c_Failure 

Failed =3 in Python.

c_Undefined 

Not yet known (before execution) =4 in Python.

Definition at line 40 of file CalibrationAlgorithm.h.

Member Function Documentation

◆ calibrate()

CalibrationAlgorithm::EResult calibrate ( )
overrideprotectedvirtual

..Run algorithm on events


ranges of various fit parameters, and tolerance to determine that fit is at the limit

Put root into batch mode so that we don't try to open a graphics window


Write out the job parameters


Clean up existing histograms if necessary


Histograms containing the data collected by eclMuMuECollectorModule


Record the number of entries per crystal in the normalized energy histogram and calculate the average expected energy per crystal and calibration constants from Collector


Write out the basic histograms in all cases


If we have not been asked to do fits, we can quit now


Check that every crystal has enough entries.

Insufficient data. Quit if we are required to have a successful fit for every crystal


Some prep for the many fits about to follow

histograms to store results for database

Diagnostic histograms

1D summary histograms


Fits are requested and there is sufficient data. Loop over specified crystals and performs fits to the amplitude distributions

Project 1D histogram of energy in this crystal

Fit function (xmin, xmax, nparameters) for this histogram

Estimate initial parameters from the histogram. For peak, use maximum bin in the allowed range

eta is nominal values


Fit

Set the initial parameters

Fit

bins on either side of this value

look for the target value between these two points


Fit status

did not find lower edge

poor fit

parameter at limit

fill diagnostic histograms

1D summary histograms

Store histogram with fit


Interpret results of fit as expected energy or calibration constant

if the fit is not successful, set peakE and edge to -1, so that output calib = -1 * abs(input calib)

Find expected energies from MC, if requested

Otherwise, calibration constant


Write output to DB if requested and successful

Store expected energies

Store calibration constants


Write out diagnostic histograms

Histograms containing values written to DB


Clean up histograms in case Algorithm is called again


Set the return code appropriately

Implements CalibrationAlgorithm.

Definition at line 58 of file eclMuMuEAlgorithm.cc.

59 {
62  double limitTol = 0.0005; /*< tolerance for checking if a parameter is at the limit */
63  double minFitLimit = 1e-25; /*< cut off for labeling a fit as poor */
64  double peakMin(0.4), peakMax(2.2); /*< range for peak of normalized energy distribution */
65  double peakTol = limitTol * (peakMax - peakMin); /*< fit is at limit if it is within peakTol of min or max */
66  double effSigMin(0.02), effSigMax(0.2); /*< range for effective sigma of normalized energy distribution */
67  double effSigTol = limitTol * (effSigMax - effSigMin); /*< fit is at limit if it is within effSigTol of min or max */
68  double etaNom(-0.41); /*< Nominal tail parameter; fixed to this value for low statistics fits */
69  double etaMin(-1.), etaMax(0.); /*< Novosibirsk tail parameter range */
70  double etaTol = limitTol * (etaMax - etaMin); /*< fit is at limit if it is within etaTol of min or max */
71 
73  gROOT->SetBatch();
74 
77  B2INFO("eclMuMuEAlgorithm parameters:");
78  B2INFO("cellIDLo = " << cellIDLo);
79  B2INFO("cellIDHi = " << cellIDHi);
80  B2INFO("minEntries = " << minEntries);
81  B2INFO("tRatioMin = " << tRatioMin);
82  B2INFO("tRatioMax = " << tRatioMax);
83  B2INFO("lowerEdgeThresh = " << lowerEdgeThresh);
84  B2INFO("performFits = " << performFits);
85  B2INFO("findExpValues = " << findExpValues);
86  B2INFO("storeConst = " << storeConst);
87 
90  TH1F* dummy;
91  dummy = (TH1F*)gROOT->FindObject("IntegralVsCrysID");
92  if (dummy) {delete dummy;}
93  dummy = (TH1F*)gROOT->FindObject("AverageExpECrys");
94  if (dummy) {delete dummy;}
95  dummy = (TH1F*)gROOT->FindObject("AverageElecCalib");
96  if (dummy) {delete dummy;}
97  dummy = (TH1F*)gROOT->FindObject("AverageInitCalib");
98  if (dummy) {delete dummy;}
99 
102  auto TrkPerCrysID = getObjectPtr<TH1F>("TrkPerCrysID");
103  auto EnVsCrysID = getObjectPtr<TH2F>("EnVsCrysID");
104  auto ExpEvsCrys = getObjectPtr<TH1F>("ExpEvsCrys");
105  auto ElecCalibvsCrys = getObjectPtr<TH1F>("ElecCalibvsCrys");
106  auto InitialCalibvsCrys = getObjectPtr<TH1F>("InitialCalibvsCrys");
107  auto CalibEntriesvsCrys = getObjectPtr<TH1F>("CalibEntriesvsCrys");
108  auto RawDigitAmpvsCrys = getObjectPtr<TH2F>("RawDigitAmpvsCrys");
109  auto RawDigitTimevsCrys = getObjectPtr<TH2F>("RawDigitTimevsCrys");
110  auto hitCrysVsExtrapolatedCrys = getObjectPtr<TH2F>("hitCrysVsExtrapolatedCrys");
111 
115  TH1F* IntegralVsCrysID = new TH1F("IntegralVsCrysID", "Integral of EnVsCrysID for each crystal;crystal ID;Entries", 8736, 0, 8736);
116  TH1F* AverageExpECrys = new TH1F("AverageExpECrys", "Average expected E per crys from collector;Crystal ID;Energy (GeV)", 8736, 0,
117  8736);
118  TH1F* AverageElecCalib = new TH1F("AverageElecCalib", "Average electronics calib const vs crystal;Crystal ID;Calibration constant",
119  8736, 0, 8736);
120  TH1F* AverageInitCalib = new TH1F("AverageInitCalib", "Average initial calib const vs crystal;Crystal ID;Calibration constant",
121  8736, 0, 8736);
122 
123  for (int crysID = 0; crysID < 8736; crysID++) {
124  TH1D* hEnergy = EnVsCrysID->ProjectionY("hEnergy", crysID + 1, crysID + 1);
125  int Integral = hEnergy->Integral();
126  IntegralVsCrysID->SetBinContent(crysID + 1, Integral);
127 
128  double TotEntries = CalibEntriesvsCrys->GetBinContent(crysID + 1);
129 
130  double expectedE = 0.;
131  if (TotEntries > 0.) {expectedE = ExpEvsCrys->GetBinContent(crysID + 1) / TotEntries;}
132  AverageExpECrys->SetBinContent(crysID + 1, expectedE);
133 
134  double calibconst = 0.;
135  if (TotEntries > 0.) {calibconst = ElecCalibvsCrys->GetBinContent(crysID + 1) / TotEntries;}
136  AverageElecCalib->SetBinContent(crysID + 1, calibconst);
137 
138  calibconst = 0.;
139  if (TotEntries > 0.) {calibconst = InitialCalibvsCrys->GetBinContent(crysID + 1) / TotEntries;}
140  AverageInitCalib->SetBinContent(crysID + 1, calibconst);
141  }
142 
145  TFile* histfile = new TFile("eclMuMuEAlgorithm.root", "recreate");
146  TrkPerCrysID->Write();
147  EnVsCrysID->Write();
148  IntegralVsCrysID->Write();
149  AverageExpECrys->Write();
150  AverageElecCalib->Write();
151  AverageInitCalib->Write();
152  RawDigitAmpvsCrys->Write();
153  RawDigitTimevsCrys->Write();
154  hitCrysVsExtrapolatedCrys->Write();
155 
158  if (!performFits) {
159  B2RESULT("eclMuMuEAlgorithm has not been asked to perform fits; copying input histograms and quitting");
160  histfile->Close();
161  return c_NotEnoughData;
162  }
163 
166  bool sufficientData = true;
167  for (int crysID = cellIDLo - 1; crysID < cellIDHi; crysID++) {
168  if (IntegralVsCrysID->GetBinContent(crysID + 1) < minEntries) {
169  if (storeConst == 1) {B2RESULT("eclMuMuEAlgorithm: crystal " << crysID << " has insufficient statistics: " << IntegralVsCrysID->GetBinContent(crysID + 1) << ". Requirement is " << minEntries);}
170  sufficientData = false;
171  break;
172  }
173  }
174 
176  if (!sufficientData && storeConst == 1) {
177  histfile->Close();
178  return c_NotEnoughData;
179  }
180 
185  TH1F* CalibVsCrysID = new TH1F("CalibVsCrysID", "Calibration constant vs crystal ID;crystal ID;counts per GeV", 8736, 0, 8736);
186  TH1F* ExpEnergyperCrys = new TH1F("ExpEnergyperCrys", "Expected energy per crystal;Crystal ID;Peak energy (GeV)", 8736, 0, 8736);
187 
189  TH1F* PeakVsCrysID = new TH1F("PeakVsCrysID", "Peak of Novo fit vs crystal ID;crystal ID;Peak normalized energy", 8736, 0,
190  8736);
191  TH1F* EdgeVsCrysID = new TH1F("EdgeVsCrysID", "Lower edge of Novo fit vs crystal ID;crystal ID", 8736, 0,
192  8736);
193  TH1F* effSigVsCrysID = new TH1F("effSigVsCrysID", "effSigma vs crystal ID;crystal ID;sigma)", 8736, 0, 8736);
194  TH1F* etaVsCrysID = new TH1F("etaVsCrysID", "eta vs crystal ID;crystal ID;Novo eta parameter", 8736, 0, 8736);
195  TH1F* normVsCrysID = new TH1F("normVsCrysID", "Novosibirsk normalization vs crystal ID;crystal ID;normalization", 8736, 0, 8736);
196  TH1F* lowerLimitVsCrysID = new TH1F("lowerLimitVsCrysID", "fit range lower limit vs crystal ID;crystal ID;lower fit limit", 8736, 0,
197  8736);
198  TH1F* fitLimitVsCrysID = new TH1F("fitLimitVsCrysID", "fit range upper limit vs crystal ID;crystal ID;upper fit limit", 8736, 0,
199  8736);
200  TH1F* StatusVsCrysID = new TH1F("StatusVsCrysID", "Fit status vs crystal ID;crystal ID;Fit status", 8736, 0, 8736);
201 
203  TH1F* hStatus = new TH1F("hStatus", "Fit status", 25, -5, 20);
204  TH1F* hPeak = new TH1F("hPeak", "Peaks of normalized energy distributions, successful fits;Peak of Novosibirsk fit", 200, 0.8, 1.2);
205  TH1F* fracPeakUnc = new TH1F("fracPeakUnc", "Fractional uncertainty on peak uncertainty, successful fits;Uncertainty on peak", 100,
206  0, 0.1);
207 
208 
211  bool allFitsOK = true;
212  for (int crysID = cellIDLo - 1; crysID < cellIDHi; crysID++) {
213 
215  TString name = "Enormalized";
216  name += crysID;
217  TH1D* hEnergy = EnVsCrysID->ProjectionY(name, crysID + 1, crysID + 1);
218 
220  double histMin = hEnergy->GetXaxis()->GetXmin();
221  double histMax = hEnergy->GetXaxis()->GetXmax();
222  TF1* func = new TF1("eclNovoConst", eclNovoConst, histMin, histMax, 5);
223  func->SetParNames("normalization", "peak", "effSigma", "eta", "const");
224  func->SetParLimits(1, peakMin, peakMax);
225  func->SetParLimits(2, effSigMin, effSigMax);
226  func->SetParLimits(3, etaMin, etaMax);
227 
228  //..Currently not using the constant term
229  double constant = 0.;
230  func->FixParameter(4, constant);
231 
232  //..If low statistics, rebin, and fix eta
233  if (hEnergy->GetEntries() < nToRebin) {
234  hEnergy->Rebin(2);
235  func->FixParameter(3, etaNom);
236  }
237 
239  hEnergy->GetXaxis()->SetRangeUser(peakMin, peakMax);
240  double peak = hEnergy->GetMaximum();
241  int maxBin = hEnergy->GetMaximumBin();
242  double peakE = hEnergy->GetBinLowEdge(maxBin);
243  double peakEUnc = 0.;
244  double normalization = hEnergy->GetMaximum();
245  double normUnc = 0.;
246  double effSigma = hEnergy->GetRMS();
247  double sigmaUnc = 0.;
248  hEnergy->GetXaxis()->SetRangeUser(histMin, histMax);
249  double fitProb = 0.;
250 
252  double eta = etaNom;
253  double etaUnc = 0.;
254 
255  //..Will find the lower edge of normalized energy at the end of the fit
256  double lowerEnEdge = 0.;
257 
258  //..Fit range from set of bins with sufficient entries.
259  double targetY = tRatioMin * peak;
260  int iLow = maxBin;
261  while (hEnergy->GetBinContent(iLow) > targetY) {iLow--;}
262  double fitlow = hEnergy->GetBinLowEdge(iLow);
263 
264  targetY = tRatioMax * peak;
265  int iHigh = maxBin;
266  while (hEnergy->GetBinContent(iHigh) > targetY) {iHigh++;}
267  double fithigh = hEnergy->GetBinLowEdge(iHigh + 1);
268 
271  bool fitHist = IntegralVsCrysID->GetBinContent(crysID + 1) >= minEntries; /* fit only if enough events */
272  if (fitHist) {
273 
275  func->SetParameters(normalization, peakE, effSigma, eta, constant);
276 
278  hEnergy->Fit(func, "LIQ", "", fitlow, fithigh);
279  normalization = func->GetParameter(0);
280  normUnc = func->GetParError(0);
281  peakE = func->GetParameter(1);
282  peakEUnc = func->GetParError(1);
283  effSigma = func->GetParameter(2);
284  sigmaUnc = func->GetParError(2);
285  eta = func->GetParameter(3);
286  etaUnc = func->GetParError(3);
287  fitProb = func->GetProb();
288 
289  //..Lower edge of the fit function is used to find the calibration constant.
290  // Can now use the peak of the fit, instead of the bin content.
291  peak = func->Eval(peakE);
292  targetY = lowerEdgeThresh * peak;
293 
295  iHigh = hEnergy->GetXaxis()->FindBin(peakE) + 1;
296  iLow = hEnergy->GetXaxis()->FindBin(fitlow);
297  int iLast = iHigh;
298  for (int ibin = iHigh; ibin > iLow; ibin--) {
299  double xc = hEnergy->GetBinCenter(ibin);
300  if (func->Eval(xc) > targetY) {iLast = ibin;}
301  }
302  double xHigh = hEnergy->GetBinCenter(iLast);
303  double xLow = hEnergy->GetBinCenter(iLast - 1);
304 
306  if (func->Eval(xLow) < targetY and func->Eval(xHigh) > targetY) {
307  func->SetNpx(1000);
308  lowerEnEdge = func->GetX(targetY, xLow, xHigh);
309  }
310  }
311 
314  int iStatus = fitOK; // success
315 
317  if (lowerEnEdge < 0.01) {iStatus = noLowerEdge;}
318 
320  if (fitProb <= minFitLimit) {iStatus = poorFit;}
321 
323  if ((peakE < peakMin + peakTol) || (peakE > peakMax - peakTol)) {iStatus = atLimit;}
324  if ((effSigma < effSigMin + effSigTol) || (effSigma > effSigMax - effSigTol)) {iStatus = atLimit;}
325  if ((eta < etaMin + etaTol) || (eta > etaMax - etaTol)) {iStatus = atLimit;}
326 
327  //** No fit
328  if (!fitHist) {iStatus = notFit;} // not fit
329 
331  int histbin = crysID + 1;
332  PeakVsCrysID->SetBinContent(histbin, peakE);
333  PeakVsCrysID->SetBinError(histbin, peakEUnc);
334  EdgeVsCrysID->SetBinContent(histbin, lowerEnEdge);
335  EdgeVsCrysID->SetBinError(histbin, peakEUnc); // approximate
336  effSigVsCrysID->SetBinContent(histbin, effSigma);
337  effSigVsCrysID->SetBinError(histbin, sigmaUnc);
338  etaVsCrysID->SetBinContent(histbin, eta);
339  etaVsCrysID->SetBinError(histbin, etaUnc);
340  normVsCrysID->SetBinContent(histbin, normalization);
341  normVsCrysID->SetBinError(histbin, normUnc);
342  lowerLimitVsCrysID->SetBinContent(histbin, fitlow);
343  lowerLimitVsCrysID->SetBinError(histbin, 0);
344  fitLimitVsCrysID->SetBinContent(histbin, fithigh);
345  fitLimitVsCrysID->SetBinError(histbin, 0);
346  StatusVsCrysID->SetBinContent(histbin, iStatus);
347  StatusVsCrysID->SetBinError(histbin, 0);
348 
350  hStatus->Fill(iStatus);
351  if (iStatus >= iterations) {
352  hPeak->Fill(peakE);
353  fracPeakUnc->Fill(peakEUnc / peakE);
354  }
355 
357  B2INFO("cellID " << crysID + 1 << " status = " << iStatus << " fit probability = " << fitProb);
358  histfile->cd();
359  hEnergy->Write();
360 
361  } /* end of loop over crystals */
362 
365  for (int crysID = 0; crysID < 8736; crysID++) {
366  int histbin = crysID + 1;
367  double fitstatus = StatusVsCrysID->GetBinContent(histbin);
368  double peakE = PeakVsCrysID->GetBinContent(histbin);
369  double edge = EdgeVsCrysID->GetBinContent(histbin);
370  double fracpeakEUnc = PeakVsCrysID->GetBinError(histbin) / peakE;
371 
373  if (fitstatus < iterations) {
374  peakE = -1.;
375  edge = -1.;
376  fracpeakEUnc = 0.;
377  if (histbin >= cellIDLo && histbin <= cellIDHi) {
378  B2RESULT("eclMuMuEAlgorithm: cellID " << histbin << " is not a successful fit. Status = " << fitstatus);
379  allFitsOK = false;
380  }
381  }
382 
384  if (findExpValues) {
385  double inputExpE = abs(AverageExpECrys->GetBinContent(histbin));
386  ExpEnergyperCrys->SetBinContent(histbin, inputExpE * edge);
387  ExpEnergyperCrys->SetBinError(histbin, fracpeakEUnc * inputExpE * peakE);
388  } else {
389 
391  double inputCalib = abs(AverageInitCalib->GetBinContent(histbin));
392  CalibVsCrysID->SetBinContent(histbin, inputCalib / edge);
393  CalibVsCrysID->SetBinError(histbin, fracpeakEUnc * inputCalib / peakE);
394  }
395  }
396 
399  bool DBsuccess = false;
400  if (storeConst == 0 || (storeConst == 1 && allFitsOK)) {
401  DBsuccess = true;
402  if (findExpValues) {
403 
405  std::vector<float> tempE;
406  std::vector<float> tempUnc;
407  for (int crysID = 0; crysID < 8736; crysID++) {
408  tempE.push_back(ExpEnergyperCrys->GetBinContent(crysID + 1));
409  tempUnc.push_back(ExpEnergyperCrys->GetBinError(crysID + 1));
410  }
411  ECLCrystalCalib* ExpectedE = new ECLCrystalCalib();
412  ExpectedE->setCalibVector(tempE, tempUnc);
413  saveCalibration(ExpectedE, "ECLExpMuMuE");
414  B2RESULT("eclCosmicEAlgorithm: successfully stored expected energies ECLExpMuMuE");
415 
416  } else {
417 
419  std::vector<float> tempCalib;
420  std::vector<float> tempCalibUnc;
421  for (int crysID = 0; crysID < 8736; crysID++) {
422  tempCalib.push_back(CalibVsCrysID->GetBinContent(crysID + 1));
423  tempCalibUnc.push_back(CalibVsCrysID->GetBinError(crysID + 1));
424  }
425  ECLCrystalCalib* MuMuECalib = new ECLCrystalCalib();
426  MuMuECalib->setCalibVector(tempCalib, tempCalibUnc);
427  saveCalibration(MuMuECalib, "ECLCrystalEnergyMuMu");
428  B2RESULT("eclMuMuEAlgorithm: successfully stored ECLCrystalEnergyMuMu calibration constants");
429  }
430  }
431 
435  PeakVsCrysID->Write();
436  EdgeVsCrysID->Write();
437  effSigVsCrysID->Write();
438  etaVsCrysID->Write();
439  normVsCrysID->Write();
440  lowerLimitVsCrysID->Write();
441  fitLimitVsCrysID->Write();
442  StatusVsCrysID->Write();
443  hPeak->Write();
444  fracPeakUnc->Write();
445  hStatus->Write();
446 
448  if (findExpValues) {
449  ExpEnergyperCrys->Write();
450  } else {
451  CalibVsCrysID->Write();
452  }
453  histfile->Close();
454 
457  dummy = (TH1F*)gROOT->FindObject("PeakVsCrysID"); delete dummy;
458  dummy = (TH1F*)gROOT->FindObject("EdgeVsCrysID"); delete dummy;
459  dummy = (TH1F*)gROOT->FindObject("effSigVsCrysID"); delete dummy;
460  dummy = (TH1F*)gROOT->FindObject("etaVsCrysID"); delete dummy;
461  dummy = (TH1F*)gROOT->FindObject("normVsCrysID"); delete dummy;
462  dummy = (TH1F*)gROOT->FindObject("lowerLimitVsCrysID"); delete dummy;
463  dummy = (TH1F*)gROOT->FindObject("fitLimitVsCrysID"); delete dummy;
464  dummy = (TH1F*)gROOT->FindObject("StatusVsCrysID"); delete dummy;
465  dummy = (TH1F*)gROOT->FindObject("fitProbSame"); delete dummy;
466  dummy = (TH1F*)gROOT->FindObject("fracPeakUnc"); delete dummy;
467  dummy = (TH1F*)gROOT->FindObject("hStatus"); delete dummy;
468  dummy = (TH1F*)gROOT->FindObject("ExpEnergyperCrys"); delete dummy;
469  dummy = (TH1F*)gROOT->FindObject("CalibVsCrysID"); delete dummy;
470 
471 
474  if (storeConst == -1) {
475  B2RESULT("eclMuMuEAlgorithm performed fits but was not asked to store contants");
476  return c_Failure;
477  } else if (!DBsuccess) {
478  if (findExpValues) { B2RESULT("eclMuMuEAlgorithm: failed to store expected values"); }
479  else { B2RESULT("eclMuMuEAlgorithm: failed to store calibration constants"); }
480  return c_Failure;
481  }
482  return c_OK;
483 }
void saveCalibration(TClonesArray *data, const std::string &name)
Store DBArray payload with given name with default IOV.
@ c_OK
Finished successfuly =0 in Python.
@ c_NotEnoughData
Needs more data =2 in Python.
@ c_Failure
Failed =3 in Python.
General DB object to store one calibration number per ECL crystal.
void setCalibVector(const std::vector< float > &CalibConst, const std::vector< float > &CalibConstUnc)
Set vector of constants with uncertainties.
double tRatioMin
entries/peak at low edge of fit must be greater than this
int poorFit
low chi square; fit not useable
int iterations
fit reached max number of iterations, but is useable
int noLowerEdge
could not determine lower edge of fit
int storeConst
controls which values are written to the database.
bool performFits
if false, input histograms are copied to output, but no fits are done.
int cellIDHi
Last cellID to be fit.
double lowerEdgeThresh
Lower edge is where the fit = lowerEdgeThresh * peak value.
int cellIDLo
..Parameters to control Novosibirsk fit to energy deposited in each crystal by mu+mu- events
bool findExpValues
if true, fits are used to find expected energy deposit for each crystal instead of the calibration co...
int minEntries
All crystals to be fit must have at least minEntries events in the fit range.
int nToRebin
If fewer entries than this, rebin and fix eta parameter.
int atLimit
a parameter is at the limit; fit not useable
double tRatioMax
entries/peak at high edge of fit must be greater than this

◆ checkPyExpRun()

bool checkPyExpRun ( PyObject *  pyObj)
inherited

Checks that a PyObject can be successfully converted to an ExpRun type.

Checks if the PyObject can be converted to ExpRun.

Definition at line 28 of file CalibrationAlgorithm.cc.

◆ convertPyExpRun()

ExpRun convertPyExpRun ( PyObject *  pyObj)
inherited

Performs the conversion of PyObject to ExpRun.

Converts the PyObject to an ExpRun. We've preoviously checked the object so this assumes a lot about the PyObject.

Definition at line 70 of file CalibrationAlgorithm.cc.

◆ execute()

EResult execute ( std::vector< Calibration::ExpRun >  runs = {},
int  iteration = 0,
IntervalOfValidity  iov = IntervalOfValidity() 
)
inherited

Runs calibration over vector of runs for a given iteration.

You can also specify the IoV to save the database payload as. By default the Algorithm will create an IoV from your requested ExpRuns, or from the overall ExpRuns of the input data if you haven't specified ExpRuns in this function.

No checks are performed to make sure that a IoV you specify matches the data you ran over, it simply labels the IoV to commit to the database later.

◆ getCollectorName()

std::string getCollectorName ( ) const
inlineinherited

Alias for prefix.

For convenience and less writing, we say developers to set this to default collector module name in constructor of base class. One can however use the dublets of collector+algorithm multiple times with different settings. To bind these together correctly, the prefix has to be set the same for algo and collector. So we call the setter setPrefix rather than setModuleName or whatever. This getter will work out of the box for default cases -> return the name of module you have to add to your path to collect data for this algorihtm.

Definition at line 164 of file CalibrationAlgorithm.h.

◆ setInputFileNames() [1/2]

void setInputFileNames ( PyObject *  inputFileNames)
inherited

Set the input file names used for this algorithm from a Python list.

Set the input file names used for this algorithm and resolve the wildcards.

Definition at line 166 of file CalibrationAlgorithm.cc.

◆ setInputFileNames() [2/2]

void setInputFileNames ( std::vector< std::string >  inputFileNames)
protectedinherited

Set the input file names used for this algorithm.

Set the input file names used for this algorithm and resolve the wildcards.

Definition at line 194 of file CalibrationAlgorithm.cc.

Member Data Documentation

◆ cellIDLo

int cellIDLo

..Parameters to control Novosibirsk fit to energy deposited in each crystal by mu+mu- events

First cellID to be fit

Definition at line 31 of file eclMuMuEAlgorithm.h.

◆ storeConst

int storeConst

controls which values are written to the database.

0 (default): store value found by successful fits, or -|input value| otherwise; -1 : do not store values 1 : store values if every fit for [cellIDLo,cellIDHi] was successful

Definition at line 40 of file eclMuMuEAlgorithm.h.


The documentation for this class was generated from the following files: