ATLAS Offline Software
Public Member Functions | Public Attributes | Private Types | Private Member Functions | Static Private Member Functions | Private Attributes | List of all members
ZDCNLCalibration Class Reference

#include <ZDCNLCalibration.h>

Collaboration diagram for ZDCNLCalibration:

Public Member Functions

 ZDCNLCalibration (const std::string &file, int maxNLPower=3, bool useGRL=true, int debugLevel=0)
 
virtual ~ZDCNLCalibration ()
 
std::array< float, 4 > FindSNPeaks (size_t LBLow, size_t LBHigh, size_t side)
 
std::pair< float, float > FindSNRange (size_t LBLow, size_t LBHigh, size_t side)
 
std::pair< std::pair< float, float >, std::pair< float, float > > FindSNTwoNRanges (size_t LBLow, size_t LBHigh, size_t side)
 
void SetDefaultCalibration (size_t side, const CalibData &calib)
 
void Calibrate (size_t side, const std::string &calibInput, const std::string &calibOutput, size_t LBLow, size_t LBHigh, std::array< int, 4 > maxPowerModule, const std::vector< std::pair< double, double > > &nNeutERange, bool excludeHE, float heSumThresh, float HEDeweight)
 
void TestCalibration (int side, const std::string &calibName)
 
TH1 * GetTestSNHist ()
 
TH1 * GetTestFracHist (size_t module)
 
TTree * GetTestTree ()
 

Public Attributes

std::array< std::map< std::string, CalibData >, 2 > m_calibrations
 
bool m_haveTest {}
 
TH1D * m_testCalibSNHist {}
 
std::array< TH1D *, 4 > m_testCalibHEFracHist
 
std::array< TH1D *, 4 > m_testCalibEnergyHist
 
TTree * m_testTree {}
 

Private Types

typedef std::multimap< unsigned int, std::pair< unsigned int, unsigned int > > LBEvtMap
 

Private Member Functions

void FillLumiBlockEvtMap ()
 
void FillMinimizationData (TMatrixD &minimMatrix, TVectorD &minimVector, std::array< int, 4 > maxPowerModule, float HEDeweight, const std::vector< std::vector< double > > &sums1DVec, const std::vector< double > &sumsHE, const std::vector< std::vector< double > > &sums2DVec, const std::vector< double > &sumsHE2D)
 
void AddCalibration (size_t side, const std::string &tag, const CalibData &calib)
 
CalibData GetCalibration (size_t side, const std::string &tag)
 
void AddToSums (std::vector< double > &sums1D, std::vector< double > &sums2D, float *amps)
 

Static Private Member Functions

static double CalculateEnergy (const float *moduleAmps, const CalibData &calib)
 

Private Attributes

TFile * m_TFile
 
TTree * m_tree
 
size_t m_maxNLPower
 
bool m_useGRL
 
int m_debugLevel
 
const float m_SNEnergy
 
const std::vector< float > m_HEFraction
 
UInt_t runNumber
 
UInt_t eventNumber
 
UInt_t lumiBlock
 
UInt_t bcid
 
UInt_t passBits
 
Float_t zdc_ZdcAmp [2]
 
UInt_t zdc_ZdcModuleMask
 
Float_t zdc_ZdcModuleAmp [2][4]
 
Bool_t L1_ZDC_A
 
Bool_t L1_ZDC_C
 
Bool_t L1_ZDC_AND
 
Bool_t L1_ZDC_A_C
 
TBranch * b_runNumber
 
TBranch * b_eventNumber
 
TBranch * b_lumiBlock
 
TBranch * b_bcid
 
TBranch * b_passBits
 
TBranch * b_zdc_ZdcAmp
 
TBranch * b_zdc_ZdcModuleMask
 
TBranch * b_zdc_ZdcModuleAmp
 
TBranch * b_L1_ZDC_A
 
TBranch * b_L1_ZDC_C
 
TBranch * b_L1_ZDC_AND
 
TBranch * b_L1_ZDC_A_C
 
LBEvtMap m_LumiBlockEvtMap
 

Detailed Description

Definition at line 55 of file ZDCNLCalibration.h.

Member Typedef Documentation

◆ LBEvtMap

typedef std::multimap<unsigned int, std::pair<unsigned int, unsigned int> > ZDCNLCalibration::LBEvtMap
private

Definition at line 107 of file ZDCNLCalibration.h.

Constructor & Destructor Documentation

◆ ZDCNLCalibration()

ZDCNLCalibration::ZDCNLCalibration ( const std::string &  file,
int  maxNLPower = 3,
bool  useGRL = true,
int  debugLevel = 0 
)

◆ ~ZDCNLCalibration()

virtual ZDCNLCalibration::~ZDCNLCalibration ( )
inlinevirtual

Definition at line 122 of file ZDCNLCalibration.h.

122 {}

Member Function Documentation

◆ AddCalibration()

void ZDCNLCalibration::AddCalibration ( size_t  side,
const std::string &  tag,
const CalibData calib 
)
private

Definition at line 51 of file ZDCNLCalibration.cxx.

52 {
53  //cppcheck-suppress mismatchingContainers
55  if (iter != m_calibrations[side].end()) (*iter).second = calib;//*iter = std::pair<std::string, CalibData>(tag, calib);
56  else {
57  m_calibrations[side].insert(std::pair<std::string, CalibData>(tag, calib));
58  }
59 }

◆ AddToSums()

void ZDCNLCalibration::AddToSums ( std::vector< double > &  sums1D,
std::vector< double > &  sums2D,
float *  amps 
)
inlineprivate

Definition at line 161 of file ZDCNLCalibration.h.

162  {
163  size_t index1D = 0;
164  size_t index2D = 0;
165 
166  for (size_t module1 : {0, 1, 2, 3}) {
167  double amp1 = amps[module1];
168  double amp1Pow = amp1;
169 
170  for (size_t power1 = 0; power1 < m_maxNLPower; power1++) {
171  for (size_t module2 : {0, 1, 2, 3}) {
172  double amp2 = amps[module2];
173  double amp2Pow = amp2;
174 
175  for (size_t power2 = 0; power2 < m_maxNLPower; power2++) {
176  sums2D.at(index2D++) += amp1Pow*amp2Pow;
177  amp2Pow *= amp2;
178  }
179  }
180 
181  sums1D.at(index1D++) += amp1Pow;
182  amp1Pow *= amp1;
183  }
184  }
185  }

◆ CalculateEnergy()

static double ZDCNLCalibration::CalculateEnergy ( const float *  moduleAmps,
const CalibData calib 
)
inlinestaticprivate

Definition at line 187 of file ZDCNLCalibration.h.

188  {
189  double energy = 0;
190 
191  for (size_t module : {0, 1, 2, 3}) {
192  float amp = moduleAmps[module];
193  float ampPow = amp;
194 
195  for (size_t power = 0; power < calib.weights[module].size(); power++) {
196  energy += calib.weights[module][power]*ampPow;
197  ampPow *= amp;
198  }
199  }
200 
201  return energy;
202  }

◆ Calibrate()

void ZDCNLCalibration::Calibrate ( size_t  side,
const std::string &  calibInput,
const std::string &  calibOutput,
size_t  LBLow,
size_t  LBHigh,
std::array< int, 4 >  maxPowerModule,
const std::vector< std::pair< double, double > > &  nNeutERange,
bool  excludeHE,
float  heSumThresh,
float  HEDeweight 
)

Definition at line 197 of file ZDCNLCalibration.cxx.

201 {
202  std::for_each(maxPowerModule.begin(), maxPowerModule.end(),
203  [](int power) {
204  if (power == 0) {
205  std::cout << "Found max power = 0, quitting" << std::endl;
206  return;
207  }
208  } );
209 
210  // Set up container for sums that we will use to calcualting the weights
211  //
212  size_t numWeights = std::accumulate(maxPowerModule.begin(), maxPowerModule.end(), 0);
213  size_t arrayDim = 4*m_maxNLPower;
214 
215  size_t arrayDim2D = arrayDim*arrayDim;
216 
217  size_t numNeutMax = nNeutERange.size();
218 
219  std::vector<double> sumsHE(arrayDim, 0);
220  std::vector<double> sums2DHE(arrayDim2D, 0);
221 
222  std::vector<std::vector<double> > sumsNeutVec(numNeutMax, std::vector<double>(arrayDim, 0));
223  std::vector<std::vector<double> > sums2DNeutVec(numNeutMax, std::vector<double>(arrayDim2D, 0));
224 
225  LBEvtMap::const_iterator LBStart = m_LumiBlockEvtMap.lower_bound(LBLow);
226  LBEvtMap::const_iterator LBEnd = m_LumiBlockEvtMap.upper_bound(LBHigh);
227 
228  // Select which trigger we're going to look at based on the side
229  //
230  const bool& triggerOpp = (side == 1 ? L1_ZDC_C : L1_ZDC_A);
231  const bool& triggerSame = (side == 1 ? L1_ZDC_A : L1_ZDC_C);
232 
233  std::vector<int> numNeutEvent(numNeutMax, 0);
234  int numHEEvent = 0;
235 
236  //
237  //
238  std::string calibName = (calibInput != "" ? calibInput : "default");
239  CalibData calib = GetCalibration(side, calibName);
240 
241  for ( LBEvtMap::const_iterator iter = LBStart; iter != LBEnd; ++iter) {
242  unsigned int entryStart = iter->second.first;
243  unsigned int entryEnd = iter->second.second;
244 
245  for (unsigned int entry = entryStart; entry < entryEnd; ++entry) {
246  Long64_t nb = m_tree->GetEntry(entry);
247  if (nb < 0) {
248  std::cout << "Error reading entry " << entry << ", quitting" << std::endl;
249  throw std::runtime_error("Error reading entry for ZDCNLCalibration");
250  }
251 
252  int testMask = zdc_ZdcModuleMask >> side*4;
253  if ((testMask & 0xf) != 0xf) continue;
254 
255  double sumAmp = CalculateEnergy(zdc_ZdcModuleAmp[side], calib);
256 
257  if (triggerSame && sumAmp > heSumThresh && !excludeHE) {
258  AddToSums(sumsHE, sums2DHE, zdc_ZdcModuleAmp[side]);
259  numHEEvent++;
260  }
261 
262  // For 1N we require the opposite side trigger
263  //
264  if (triggerOpp && sumAmp > nNeutERange[0].first && sumAmp < nNeutERange[0].second) {
265  AddToSums(sumsNeutVec[0], sums2DNeutVec[0], zdc_ZdcModuleAmp[side]);
266  numNeutEvent[0]++;
267  }
268 
269  for (size_t iNNeut = 1; iNNeut < numNeutMax; iNNeut++) {
270  if (sumAmp > nNeutERange[iNNeut].first && sumAmp < nNeutERange[iNNeut].second) {
271  AddToSums(sumsNeutVec[iNNeut], sums2DNeutVec[iNNeut], zdc_ZdcModuleAmp[side]);
272  numNeutEvent[iNNeut]++;
273  }
274  }
275  }
276  }
277 
278  if (m_debugLevel >= 0) {
279  std::cout << "Statistics for side " << side << ", lb range = " << LBLow << ", " << LBHigh
280  << " : # high energy events = " << numHEEvent;
281 
282  for (size_t iNNeut = 0; iNNeut < numNeutMax; iNNeut++) {
283  std::cout << ", # of " + std::to_string(iNNeut) + " N events = " << numNeutEvent[iNNeut];
284  }
285 
286  std::cout << std::endl;
287  }
288 
289 
290  // Normalize the averages
291  //
292  for (size_t iNNeut = 0; iNNeut < numNeutMax; iNNeut++) {
293  if (numNeutEvent[iNNeut] == 0) {
294  std::cout << "Found no " + std::to_string(iNNeut) + "N events, quitting" << std::endl;
295  return;
296  }
297 
298  std::for_each(sumsNeutVec[iNNeut].begin(), sumsNeutVec[iNNeut].end(), [&](double& sum) {sum /= numNeutEvent[iNNeut];});
299  std::for_each(sums2DNeutVec[iNNeut].begin(), sums2DNeutVec[iNNeut].end(), [&](double& sum) {sum /= numNeutEvent[iNNeut];});
300  }
301 
302  if (!excludeHE) {
303  if (numHEEvent == 0) {
304  std::cout << "Found no high-energy events, quitting" << std::endl;
305  return;
306  }
307 
308  std::for_each(sumsHE.begin(), sumsHE.end(), [&](double& sum) {sum /= numHEEvent;});
309  std::for_each(sums2DHE.begin(), sums2DHE.end(), [&](double& sum) {sum /= numHEEvent;});
310  }
311 
312 
313  // Now we have the averages, set up to do the minimization
314  //
315  size_t numFitParams = numWeights;
316  TMatrixD minimMatrix(numFitParams, numFitParams);
317  TVectorD minimVector(numFitParams);
318 
319  FillMinimizationData(minimMatrix, minimVector, maxPowerModule, HEDeweight,
320  sumsNeutVec, sumsHE, sums2DNeutVec, sums2DHE);
321 
322  if (m_debugLevel >= 2) {
323  std::cout << "Dumping matrix equation: " << std::endl;
324 
325  for (size_t index1 = 0; index1 < numFitParams; index1++) {
326  std::cout << "| ";
327 
328  for (size_t index2 = 0; index2 < numFitParams; index2++) {
329  std::cout << std::scientific << std::setw(9) << std::setprecision(2) << minimMatrix(index1, index2) << " ";
330  }
331 
332  std::cout << "| | w_" << index1 << " = | " << minimVector(index1) << std::endl;
333  }
334  std::cout << std::endl;
335  }
336 
337 
338  TDecompLU lu(minimMatrix);
339  lu.Decompose();
340 
341  bool status = lu.Solve(minimVector);
342 
343  if (!status) {
344  std::cout << "Minimization failed for side " << side << ", lb range "
345  << LBLow << ", " << LBHigh << endl;
346  }
347  else {
348  CalibData calibData;
349  calibData.LBStart = LBLow;
350  calibData.LBEnd = LBHigh;
351 
352  int weightIndex = 0;
353  for (size_t module : {0, 1, 2, 3}) {
354  calibData.weights[module].assign(maxPowerModule[module], 0);
355 
356  for (size_t power = 0; power < maxPowerModule[module]; power++) {
357  calibData.weights.at(module).at(power) = minimVector(weightIndex++);
358 
359  if (m_debugLevel >= 0) {
360  std::cout << "Weight for module,power " << module << ", " << power << " = " << calibData.weights.at(module).at(power)
361  << std::endl;
362  }
363  }
364  }
365 
366  // Now store this calibration
367  //
368  AddCalibration(side, calibOutput, calibData);
369  }
370 }

◆ FillLumiBlockEvtMap()

void ZDCNLCalibration::FillLumiBlockEvtMap ( )
private

Definition at line 17 of file ZDCNLCalibration.cxx.

18 {
19  Long64_t numEntries = m_tree->GetEntriesFast();
20 
21  unsigned int countPassed = 0;
22 
23  int lastLumiBlock = -1;
24  int LBStartEntry = -1;
25 
26  for (Long64_t entry = 0; entry < numEntries; entry++) {
27  Long64_t numBytes = m_tree->GetEntry(entry);
28 
29  if (m_useGRL && passBits !=0) {
30  //std::cout << "lumiBlock = " << lumiBlock << " passBits = " << passBits << std::endl;
31  continue; // reject based on GRL or errors
32  }
33 
34  if (lumiBlock != lastLumiBlock) {
35  if (lastLumiBlock > -1) {
36  m_LumiBlockEvtMap.insert(LBEvtMap::value_type(lumiBlock, std::pair<unsigned int, unsigned int>(LBStartEntry, entry - 1)));
37  }
38 
39  LBStartEntry = entry;
40  lastLumiBlock = lumiBlock;
41  }
42 
43  countPassed++;
44  }
45 
46  if (m_debugLevel >= 0) {
47  std::cout << "From a total of " << numEntries << " entries, " << countPassed << " were selected for analysis" << std::endl;
48  }
49 }

◆ FillMinimizationData()

void ZDCNLCalibration::FillMinimizationData ( TMatrixD &  minimMatrix,
TVectorD &  minimVector,
std::array< int, 4 >  maxPowerModule,
float  HEDeweight,
const std::vector< std::vector< double > > &  sums1DVec,
const std::vector< double > &  sumsHE,
const std::vector< std::vector< double > > &  sums2DVec,
const std::vector< double > &  sumsHE2D 
)
private

Definition at line 484 of file ZDCNLCalibration.cxx.

488 {
489  size_t numWeights = minimVector.GetNrows();
490  size_t numColumns = minimMatrix.GetNcols();
491  size_t numRows = minimMatrix.GetNrows();
492  size_t numNeutMax = sums1DVec.size();
493 
494  if (numRows != numWeights || numColumns != numWeights){
495  throw std::runtime_error("ZDCNLCalibration::FillMinimizationData; incompatible row/col numbers");
496  }
497 
498  // if (m_debugLevel > 0) std::cout << "Filling single neutorn entries in matrix" << std::endl;
499 
500  double sumFracSqr = 0;
501  for (size_t module : {0, 1, 2, 3}) {
502  sumFracSqr += m_HEFraction.at(module)*m_HEFraction.at(module);
503  }
504 
505  // Start with the single neutron contributions
506  //
507  {
508  size_t index1 = 0;
509 
510  for (size_t module1 : {0, 1, 2, 3}) {
511  for (size_t power1 = 0; power1 < maxPowerModule[module1]; power1++) {
512 
513  size_t index2 = 0;
514  size_t sumIndexOffset = (module1*m_maxNLPower + power1)*4*m_maxNLPower;
515 
516  for (size_t module2 : {0, 1, 2, 3}) {
517  for (size_t power2 = 0; power2 < maxPowerModule[module2]; power2++) {
518 
519  // The N neutron contributions
520  //
521  unsigned int sum2DIndex = sumIndexOffset + power2;
522  double element = 0;
523  std::for_each(sums2DVec.begin(), sums2DVec.end(), [&] (const std::vector<double>& sums) {element += 2*sums.at(sum2DIndex);});
524 
525  // The TN contribution
526  //
527  // element += 2*sumsTN2D.at(sumIndexOffset + power2);
528 
529  // The HE fraction contribution
530  //
531  double sum = sumsHE2D.at(sumIndexOffset + power2);
532  sum *= HEDeweight;
533 
534  if (module2 == module1) element += 2*sum;
535  element -= 2*m_HEFraction.at(module2)*sum;
536  element -= 2*m_HEFraction.at(module1)*sum;
537  element += 2*sumFracSqr*sum;
538 
539  minimMatrix(index1, index2) += element;
540 
541  if (m_debugLevel > 2) {
542  std::cout << "Filling module 1, power 1 = " << module1 << ", " << power1 << ", module 2, power 2 = "
543  << module2 << ", " << power2 << " with sum at index " << sumIndexOffset + power2
544  << ", value = " << minimMatrix(index1, index2) << std::endl;
545  }
546 
547  index2++;
548  }
549 
550  sumIndexOffset += m_maxNLPower;
551  }
552 
553  size_t sum1DIndex = module1*m_maxNLPower + power1;
554 
555  minimVector[index1] = 0;
556  for (unsigned int iNNeut = 0; iNNeut < numNeutMax; iNNeut++) {
557  minimVector[index1] += 2*((double) iNNeut + 1)*m_SNEnergy*sums1DVec[iNNeut].at(sum1DIndex);
558  }
559 
560  if (m_debugLevel > 2) {
561  std::cout << "Filling vector module , power = " << module1 << ", " << power1 << " with sum at index " << sum1DIndex
562  << ", value = " << minimVector[index1] << std::endl;
563  }
564 
565  // minimVector[index1] = 2*m_SNEnergy*sumsSN.at(sum1DIndex) + 2*2*m_SNEnergy*sumsTN.at(sum1DIndex);
566  index1++;
567  }
568  }
569  }
570 
571  // if (m_debugLevel > 0) std::cout << "Filling high energy entries in matrix" << std::endl;
572 
573  // {
574 
575  // size_t sumIndexOffset = 0;
576  // size_t index1 = 0;
577 
578  // for (size_t module1 : {0, 1, 2, 3}) {
579  // for (size_t power1 = 0; power1 < maxPowerModule[module1]; power1++) {
580  // size_t index2 = 0;
581 
582  // for (size_t module2 : {0, 1, 2, 3}) {
583  // for (size_t power2 = 0; power2 < maxPowerModule[module2]; power2++) {
584 
585  // double element = 0;
586  // double sum = sumsHE2D.at(sumIndexOffset + power2);
587 
588  // // Now assemble to terms contributing to this matrix element
589  // //
590  // if (module2 == module1) element += 2*sum;
591  // element -= 2*m_HEFraction.at(module2)*sum;
592  // element -= 2*m_HEFraction.at(module1)*sum;
593  // element += 2*sumFracSqr*sum;
594 
595  // minimMatrix(index1, index2) += m_HEDeweight*element;
596  // index2++;
597  // }
598 
599  // sumIndexOffset += m_maxNLPower;
600  // }
601 
602  // index1++;
603  // }
604  // }
605  // }
606 
607 }

◆ FindSNPeaks()

std::array<float, 4> ZDCNLCalibration::FindSNPeaks ( size_t  LBLow,
size_t  LBHigh,
size_t  side 
)

◆ FindSNRange()

std::pair< float, float > ZDCNLCalibration::FindSNRange ( size_t  LBLow,
size_t  LBHigh,
size_t  side 
)

Definition at line 72 of file ZDCNLCalibration.cxx.

73 {
74  static TH1D* snHist = new TH1D("ZDCNLCalibSNHist", "", 200, 0, 2000);
75  snHist->Reset();
76 
77  LBEvtMap::const_iterator LBStart = m_LumiBlockEvtMap.lower_bound(LBLow);
78  LBEvtMap::const_iterator LBEnd = m_LumiBlockEvtMap.upper_bound(LBHigh);
79 
80  // Select which trigger we're going to look at based on the side
81  //
82  const bool& trigger = (side == 1 ? L1_ZDC_C : L1_ZDC_A);
83 
84  for ( LBEvtMap::const_iterator iter = LBStart; iter != LBEnd; ++iter) {
85  unsigned int entryStart = iter->second.first;
86  unsigned int entryEnd = iter->second.second;
87 
88  for (unsigned int entry = entryStart; entry < entryEnd; entry++) {
89  Long64_t nb = m_tree->GetEntry(entry);
90  if (nb < 0) {
91  std::cout << "Error reading entry " << entry << ", quitting" << std::endl;
92  throw std::runtime_error("Error reading entry for ZDCNLCalibration");
93  }
94 
95  if (!trigger) continue;
96 
97  int testMask = zdc_ZdcModuleMask >> side*4;
98  if ((testMask & 0xf) != 0xf) continue;
99 
100  float sumAmp = zdc_ZdcAmp[side];
101  snHist->Fill(sumAmp);
102  }
103  }
104 
105  int ibin = snHist->GetMaximumBin();
106  double center = snHist->GetBinCenter(ibin);
107  double maxCounts = snHist->GetBinContent(ibin);
108 
109  static TF1* gaussFunc = new TF1("gaussFunc", "gaus", 100, 2000);
110 
111  gaussFunc->SetParameters(maxCounts, center, 100);
112  snHist->Fit("gaussFunc", "", "", center*0.7, center*1.3);
113 
114  double gaussMean = gaussFunc->GetParameter(1);
115  double gaussWidth = gaussFunc->GetParameter(2);
116 
117  double cutLow = gaussMean -2.*gaussWidth;
118  double cutHigh = gaussMean + 2*gaussWidth;
119 
120  if (m_debugLevel >= 0) {
121  std::cout << "SN cut range = " << cutLow << " to " << cutHigh << std::endl;
122  }
123 
124  return std::pair<float, float>(cutLow, cutHigh);
125 }

◆ FindSNTwoNRanges()

std::pair< std::pair< float, float >, std::pair< float, float > > ZDCNLCalibration::FindSNTwoNRanges ( size_t  LBLow,
size_t  LBHigh,
size_t  side 
)

Definition at line 128 of file ZDCNLCalibration.cxx.

129 {
130  static TH1D* sntnHist = new TH1D("ZDCNLCalibSNTNHist", "", 200, 0, 2000);
131  sntnHist->Reset();
132 
133  LBEvtMap::const_iterator LBStart = m_LumiBlockEvtMap.lower_bound(LBLow);
134  LBEvtMap::const_iterator LBEnd = m_LumiBlockEvtMap.upper_bound(LBHigh);
135 
136  // Select which trigger we're going to look at based on the side
137  //
138  const bool& trigger = (side == 1 ? L1_ZDC_C : L1_ZDC_A);
139 
140  for ( LBEvtMap::const_iterator iter = LBStart; iter != LBEnd; ++iter) {
141  unsigned int entryStart = iter->second.first;
142  unsigned int entryEnd = iter->second.second;
143 
144  for (unsigned int entry = entryStart; entry < entryEnd; ++entry) {
145  Long64_t nb = m_tree->GetEntry(entry);
146  if (nb < 0) {
147  std::cout << "Error reading entry " << entry << ", quitting" << std::endl;
148  throw std::runtime_error("Error reading entry for ZDCNLCalibration");
149  }
150 
151  if (!trigger) continue;
152 
153  int testMask = zdc_ZdcModuleMask >> side*4;
154  if ((testMask & 0xf) != 0xf) continue;
155 
156  float sumAmp = zdc_ZdcAmp[side];
157  sntnHist->Fill(sumAmp);
158  }
159  }
160 
161  int ibin = sntnHist->GetMaximumBin();
162  double center = sntnHist->GetBinCenter(ibin);
163  double maxCounts = sntnHist->GetBinContent(ibin);
164 
165  static TF1* twogaussFunc = new TF1("twogaussFunc", "[0]*exp(-pow((x-[1])/[2],2)/2.) + [3]*exp(-pow((x-[4])/[5],2)/2.)", 100, 2000);
166 
167  twogaussFunc->SetParameters(maxCounts, center, 100, maxCounts/2, 2*center, 200);
168  sntnHist->Fit("twogaussFunc", "", "", center*0.7, 2*center*1.3);
169 
170  double gaussMean1 = twogaussFunc->GetParameter(1);
171  double gaussMean2 = twogaussFunc->GetParameter(4);
172 
173  double gaussWidth1 = std::abs(twogaussFunc->GetParameter(2));
174  double gaussWidth2 = std::abs(twogaussFunc->GetParameter(5));
175 
176  double cutLow1 = gaussMean1 -2.*gaussWidth1;
177  double cutHigh1 = gaussMean1 + 2*gaussWidth1
178 ;
179  double cutLow2 = gaussMean2 -2.*gaussWidth2;
180  double cutHigh2 = gaussMean2 + 2*gaussWidth2;
181 
182  if (cutHigh1 > cutLow2) {
183  double middle = 0.5*(cutHigh1 + cutLow2);
184  cutHigh1 = middle;
185  cutLow2 = middle;
186  }
187 
188  if (m_debugLevel >= 0) {
189  std::cout << "SN cut range = " << cutLow1 << " to " << cutHigh1 << std::endl;
190  std::cout << "TN cut range = " << cutLow2 << " to " << cutHigh2 << std::endl;
191  }
192 
193  return std::pair<std::pair<float, float>, std::pair<float, float> >(std::pair<float, float>(cutLow1, cutHigh1),
194  std::pair<float, float>(cutLow2, cutHigh2));
195 }

◆ GetCalibration()

CalibData ZDCNLCalibration::GetCalibration ( size_t  side,
const std::string &  tag 
)
private

Definition at line 61 of file ZDCNLCalibration.cxx.

62 {
63  CalibData null;
64 
65  auto iter = m_calibrations.at(side).find(tag);
66  if (iter != m_calibrations.at(side).end()) {
67  return iter->second;
68  }
69  else return null;
70 }

◆ GetTestFracHist()

TH1* ZDCNLCalibration::GetTestFracHist ( size_t  module)
inline

Definition at line 142 of file ZDCNLCalibration.h.

142 {return (m_haveTest ? m_testCalibHEFracHist.at(module) : 0);}

◆ GetTestSNHist()

TH1* ZDCNLCalibration::GetTestSNHist ( )
inline

Definition at line 141 of file ZDCNLCalibration.h.

141 {return (m_haveTest ? m_testCalibSNHist : 0);}

◆ GetTestTree()

TTree* ZDCNLCalibration::GetTestTree ( )
inline

Definition at line 144 of file ZDCNLCalibration.h.

144 {return m_testTree;}

◆ SetDefaultCalibration()

void ZDCNLCalibration::SetDefaultCalibration ( size_t  side,
const CalibData calib 
)
inline

Definition at line 130 of file ZDCNLCalibration.h.

130  {
131  AddCalibration(side, "default", calib);
132  }

◆ TestCalibration()

void ZDCNLCalibration::TestCalibration ( int  side,
const std::string &  calibName 
)

Definition at line 372 of file ZDCNLCalibration.cxx.

373 {
374  std::string calibName = (name != "" ? name : "default");
375  CalibData calib = GetCalibration(side, calibName);
376 
377  if (m_debugLevel >= 0) {
378  std::cout << "Testing calibration for side " << side << "name: " << calibName << std::endl;
379  }
380 
381  m_testCalibSNHist = new TH1D("testCalibSNHist" , "", 4000, 0, 200000);
382  for (int module : {0, 1, 2, 3}) {
383  std::string fracHistName = "testCalibHEFracHist_" + std::to_string(module);
384  std::string energyHistName = "testCalibEnergyHist_" + std::to_string(module);
385 
386  m_testCalibHEFracHist[module] = new TH1D(fracHistName.c_str(), "", 200, 0., 1.);
387  m_testCalibEnergyHist[module] = new TH1D(energyHistName.c_str(), "", 1000, 0., 100000);
388  }
389 
390  static std::array<float, 4> treeRawAmp;
391  static std::array<float, 4> treeCalibAmp;
392  static float treeRawSum, treeCalibSum, treeRawSumOpp;
393  static unsigned int treeEventNumber;
394  static bool passedSN, passedTN, passedHE;
395 
396  if (m_testTree != 0) delete m_testTree;
397 
398  m_testTree = new TTree("ZDCNLCalibTest", "ZDCNLCalibTest");
399  m_testTree->SetDirectory(0);
400 
401  m_testTree->Branch("EventNumber", &treeEventNumber, "eventNumber/I");
402  m_testTree->Branch("lumiBlock", &lumiBlock, "lumiBlock/I");
403  m_testTree->Branch("L1_ZDC_A", &L1_ZDC_A, "L1_ZDC_A/b");
404  m_testTree->Branch("L1_ZDC_C", &L1_ZDC_C, "L1_ZDC_C/b");
405  m_testTree->Branch("L1_ZDC_A_C", &L1_ZDC_A_C, "L1_ZDC_A_C/b");
406 
407  m_testTree->Branch("RawSum", &treeRawSum, "RawSum/f");
408  m_testTree->Branch("CalibSum", &treeCalibSum, "CalibSum/f");
409  m_testTree->Branch("RawAmp", &treeRawAmp[0], "RawAmp[4]/f");
410  m_testTree->Branch("CalibAmp", &treeCalibAmp[0], "CalibAmp[4]/f");
411  m_testTree->Branch("RawSumOpp", &treeRawSumOpp, "RawSumOpp/f");
412 
413  LBEvtMap::const_iterator LBStart = m_LumiBlockEvtMap.lower_bound(calib.LBStart);
414  LBEvtMap::const_iterator LBEnd = m_LumiBlockEvtMap.upper_bound(calib.LBEnd);
415 
416  int count = 0;
417  for ( LBEvtMap::const_iterator iter = LBStart; iter != LBEnd; ++iter) {
418  unsigned int entryStart = iter->second.first;
419  unsigned int entryEnd = iter->second.second;
420 
421  for (unsigned int entry = entryStart; entry < entryEnd; ++entry) {
422  Long64_t nb = m_tree->GetEntry(entry);
423  if (nb < 0) {
424  std::cout << "Error reading entry " << entry << ", quitting" << std::endl;
425  throw std::runtime_error("Error reading entry for ZDCNLCalibration");
426  }
427 
428  int testMask = zdc_ZdcModuleMask >> side*4;
429  if ((testMask & 0xf) != 0xf) continue;
430 
431  double sumAmp = zdc_ZdcAmp[side];
432 
433  // Calculate the calibrated energy
434  //
435  float sumCalibAmp = 0;
436  float weightIndex = 0;
437 
438  std::array<double, 4> calibAmpModule = {0, 0, 0, 0};
439 
440  if (m_debugLevel > 2 && count < 20) {
441  std::cout << "TestCalibration:: dumping calibration results for event " << count << std::endl;
442  }
443 
444  for (size_t module : {0, 1, 2, 3}) {
445  float amp = zdc_ZdcModuleAmp[side][module];
446  float ampPow = amp;
447 
448  for (size_t power = 0; power < calib.weights[module].size(); power++) {
449  calibAmpModule[module] += calib.weights[module][power]*ampPow;
450  ampPow *= amp;
451  }
452 
453  double energy = calibAmpModule[module];
454  sumCalibAmp += energy;
456 
457  treeRawAmp[module] = amp;
458  treeCalibAmp[module] = energy;
459 
460  if (m_debugLevel > 2 && count < 20) {
461  std::cout << "Module " << module << " : raw, calibrated energy = " << std::fixed << amp << ", " << energy << std::endl;
462  }
463  }
464 
465  if (m_debugLevel > 2 && count < 20) {
466  std::cout << "Total: raw, calibrated energy = " << std::fixed << sumAmp << ", " << sumCalibAmp << std::endl;
467  }
468 
469  treeRawSum = sumAmp;
470  treeRawSumOpp = (side == 0 ? zdc_ZdcAmp[1] : zdc_ZdcAmp[0]);
471  treeCalibSum = sumCalibAmp;
472 
473  m_testTree->Fill();
474 
475  m_testCalibSNHist->Fill(sumCalibAmp);
476 
477  count++;
478  }
479  }
480 
481  m_haveTest = true;
482 }

Member Data Documentation

◆ b_bcid

TBranch* ZDCNLCalibration::b_bcid
private

Definition at line 90 of file ZDCNLCalibration.h.

◆ b_eventNumber

TBranch* ZDCNLCalibration::b_eventNumber
private

Definition at line 88 of file ZDCNLCalibration.h.

◆ b_L1_ZDC_A

TBranch* ZDCNLCalibration::b_L1_ZDC_A
private

Definition at line 99 of file ZDCNLCalibration.h.

◆ b_L1_ZDC_A_C

TBranch* ZDCNLCalibration::b_L1_ZDC_A_C
private

Definition at line 102 of file ZDCNLCalibration.h.

◆ b_L1_ZDC_AND

TBranch* ZDCNLCalibration::b_L1_ZDC_AND
private

Definition at line 101 of file ZDCNLCalibration.h.

◆ b_L1_ZDC_C

TBranch* ZDCNLCalibration::b_L1_ZDC_C
private

Definition at line 100 of file ZDCNLCalibration.h.

◆ b_lumiBlock

TBranch* ZDCNLCalibration::b_lumiBlock
private

Definition at line 89 of file ZDCNLCalibration.h.

◆ b_passBits

TBranch* ZDCNLCalibration::b_passBits
private

Definition at line 92 of file ZDCNLCalibration.h.

◆ b_runNumber

TBranch* ZDCNLCalibration::b_runNumber
private

Definition at line 87 of file ZDCNLCalibration.h.

◆ b_zdc_ZdcAmp

TBranch* ZDCNLCalibration::b_zdc_ZdcAmp
private

Definition at line 94 of file ZDCNLCalibration.h.

◆ b_zdc_ZdcModuleAmp

TBranch* ZDCNLCalibration::b_zdc_ZdcModuleAmp
private

Definition at line 97 of file ZDCNLCalibration.h.

◆ b_zdc_ZdcModuleMask

TBranch* ZDCNLCalibration::b_zdc_ZdcModuleMask
private

Definition at line 96 of file ZDCNLCalibration.h.

◆ bcid

UInt_t ZDCNLCalibration::bcid
private

Definition at line 73 of file ZDCNLCalibration.h.

◆ eventNumber

UInt_t ZDCNLCalibration::eventNumber
private

Definition at line 71 of file ZDCNLCalibration.h.

◆ L1_ZDC_A

Bool_t ZDCNLCalibration::L1_ZDC_A
private

Definition at line 81 of file ZDCNLCalibration.h.

◆ L1_ZDC_A_C

Bool_t ZDCNLCalibration::L1_ZDC_A_C
private

Definition at line 84 of file ZDCNLCalibration.h.

◆ L1_ZDC_AND

Bool_t ZDCNLCalibration::L1_ZDC_AND
private

Definition at line 83 of file ZDCNLCalibration.h.

◆ L1_ZDC_C

Bool_t ZDCNLCalibration::L1_ZDC_C
private

Definition at line 82 of file ZDCNLCalibration.h.

◆ lumiBlock

UInt_t ZDCNLCalibration::lumiBlock
private

Definition at line 72 of file ZDCNLCalibration.h.

◆ m_calibrations

std::array<std::map<std::string, CalibData>, 2> ZDCNLCalibration::m_calibrations

Definition at line 112 of file ZDCNLCalibration.h.

◆ m_debugLevel

int ZDCNLCalibration::m_debugLevel
private

Definition at line 62 of file ZDCNLCalibration.h.

◆ m_haveTest

bool ZDCNLCalibration::m_haveTest {}

Definition at line 114 of file ZDCNLCalibration.h.

◆ m_HEFraction

const std::vector<float> ZDCNLCalibration::m_HEFraction
private

Definition at line 67 of file ZDCNLCalibration.h.

◆ m_LumiBlockEvtMap

LBEvtMap ZDCNLCalibration::m_LumiBlockEvtMap
private

Definition at line 108 of file ZDCNLCalibration.h.

◆ m_maxNLPower

size_t ZDCNLCalibration::m_maxNLPower
private

Definition at line 60 of file ZDCNLCalibration.h.

◆ m_SNEnergy

const float ZDCNLCalibration::m_SNEnergy
private

Definition at line 66 of file ZDCNLCalibration.h.

◆ m_testCalibEnergyHist

std::array<TH1D*, 4> ZDCNLCalibration::m_testCalibEnergyHist

Definition at line 117 of file ZDCNLCalibration.h.

◆ m_testCalibHEFracHist

std::array<TH1D*, 4> ZDCNLCalibration::m_testCalibHEFracHist

Definition at line 116 of file ZDCNLCalibration.h.

◆ m_testCalibSNHist

TH1D* ZDCNLCalibration::m_testCalibSNHist {}

Definition at line 115 of file ZDCNLCalibration.h.

◆ m_testTree

TTree* ZDCNLCalibration::m_testTree {}

Definition at line 118 of file ZDCNLCalibration.h.

◆ m_TFile

TFile* ZDCNLCalibration::m_TFile
private

Definition at line 57 of file ZDCNLCalibration.h.

◆ m_tree

TTree* ZDCNLCalibration::m_tree
private

Definition at line 58 of file ZDCNLCalibration.h.

◆ m_useGRL

bool ZDCNLCalibration::m_useGRL
private

Definition at line 61 of file ZDCNLCalibration.h.

◆ passBits

UInt_t ZDCNLCalibration::passBits
private

Definition at line 74 of file ZDCNLCalibration.h.

◆ runNumber

UInt_t ZDCNLCalibration::runNumber
private

Definition at line 70 of file ZDCNLCalibration.h.

◆ zdc_ZdcAmp

Float_t ZDCNLCalibration::zdc_ZdcAmp[2]
private

Definition at line 76 of file ZDCNLCalibration.h.

◆ zdc_ZdcModuleAmp

Float_t ZDCNLCalibration::zdc_ZdcModuleAmp[2][4]
private

Definition at line 79 of file ZDCNLCalibration.h.

◆ zdc_ZdcModuleMask

UInt_t ZDCNLCalibration::zdc_ZdcModuleMask
private

Definition at line 78 of file ZDCNLCalibration.h.


The documentation for this class was generated from the following files:
ZDCNLCalibration::m_debugLevel
int m_debugLevel
Definition: ZDCNLCalibration.h:62
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
python.SystemOfUnits.second
int second
Definition: SystemOfUnits.py:120
ZDCNLCalibration::L1_ZDC_A
Bool_t L1_ZDC_A
Definition: ZDCNLCalibration.h:81
ZDCNLCalibration::CalculateEnergy
static double CalculateEnergy(const float *moduleAmps, const CalibData &calib)
Definition: ZDCNLCalibration.h:187
ZDCNLCalibration::m_calibrations
std::array< std::map< std::string, CalibData >, 2 > m_calibrations
Definition: ZDCNLCalibration.h:112
ZDCNLCalibration::m_testCalibSNHist
TH1D * m_testCalibSNHist
Definition: ZDCNLCalibration.h:115
CalibData::LBEnd
unsigned int LBEnd
Definition: ZDCNLCalibration.h:38
ZDCNLCalibration::m_testTree
TTree * m_testTree
Definition: ZDCNLCalibration.h:118
ZDCNLCalibration::m_testCalibHEFracHist
std::array< TH1D *, 4 > m_testCalibHEFracHist
Definition: ZDCNLCalibration.h:116
accumulate
bool accumulate(AccumulateMap &map, std::vector< module_t > const &modules, FPGATrackSimMatrixAccumulator const &acc)
Accumulates an accumulator (e.g.
Definition: FPGATrackSimMatrixAccumulator.cxx:22
PlotCalibFromCool.begin
begin
Definition: PlotCalibFromCool.py:94
ZDCNLCalibration::m_LumiBlockEvtMap
LBEvtMap m_LumiBlockEvtMap
Definition: ZDCNLCalibration.h:108
ZDCNLCalibration::zdc_ZdcAmp
Float_t zdc_ZdcAmp[2]
Definition: ZDCNLCalibration.h:76
PURW_create_mc_default_profile.LBStart
LBStart
Definition: PURW_create_mc_default_profile.py:58
ZDCNLCalibration::GetCalibration
CalibData GetCalibration(size_t side, const std::string &tag)
Definition: ZDCNLCalibration.cxx:61
XMLtoHeader.count
count
Definition: XMLtoHeader.py:85
mergePhysValFiles.end
end
Definition: DataQuality/DataQualityUtils/scripts/mergePhysValFiles.py:93
TRT::Hit::side
@ side
Definition: HitInfo.h:83
ZDCNLCalibration::AddCalibration
void AddCalibration(size_t side, const std::string &tag, const CalibData &calib)
Definition: ZDCNLCalibration.cxx:51
python.PyAthena.module
module
Definition: PyAthena.py:131
CalibData::weights
std::array< std::vector< float >, 4 > weights
Definition: ZDCNLCalibration.h:36
ZDCNLCalibration::m_tree
TTree * m_tree
Definition: ZDCNLCalibration.h:58
convertTimingResiduals.sum
sum
Definition: convertTimingResiduals.py:55
Trk::index1
@ index1
Definition: BoundarySurfaceFace.h:48
ParticleGun_FastCalo_ChargeFlip_Config.energy
energy
Definition: ParticleGun_FastCalo_ChargeFlip_Config.py:78
CalibData
Definition: ZDCNLCalibration.h:35
ZDCNLCalibration::m_maxNLPower
size_t m_maxNLPower
Definition: ZDCNLCalibration.h:60
ZDCNLCalibration::m_haveTest
bool m_haveTest
Definition: ZDCNLCalibration.h:114
xAOD::double
double
Definition: CompositeParticle_v1.cxx:159
ZDCNLCalibration::m_useGRL
bool m_useGRL
Definition: ZDCNLCalibration.h:61
GetAllXsec.entry
list entry
Definition: GetAllXsec.py:132
Trk::index2
@ index2
Definition: BoundarySurfaceFace.h:49
ZDCNLCalibration::lumiBlock
UInt_t lumiBlock
Definition: ZDCNLCalibration.h:72
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:228
ActsTrk::to_string
std::string to_string(const DetectorType &type)
Definition: GeometryDefs.h:34
CalibData::LBStart
unsigned int LBStart
Definition: ZDCNLCalibration.h:37
PlotSFuncertainty.calib
calib
Definition: PlotSFuncertainty.py:110
ZDCNLCalibration::FillMinimizationData
void FillMinimizationData(TMatrixD &minimMatrix, TVectorD &minimVector, std::array< int, 4 > maxPowerModule, float HEDeweight, const std::vector< std::vector< double > > &sums1DVec, const std::vector< double > &sumsHE, const std::vector< std::vector< double > > &sums2DVec, const std::vector< double > &sumsHE2D)
Definition: ZDCNLCalibration.cxx:484
ZDCNLCalibration::AddToSums
void AddToSums(std::vector< double > &sums1D, std::vector< double > &sums2D, float *amps)
Definition: ZDCNLCalibration.h:161
DeMoScan.first
bool first
Definition: DeMoScan.py:536
merge.status
status
Definition: merge.py:17
ZDCNLCalibration::zdc_ZdcModuleAmp
Float_t zdc_ZdcModuleAmp[2][4]
Definition: ZDCNLCalibration.h:79
ZDCNLCalibration::passBits
UInt_t passBits
Definition: ZDCNLCalibration.h:74
CaloCondBlobAlgs_fillNoiseFromASCII.tag
string tag
Definition: CaloCondBlobAlgs_fillNoiseFromASCII.py:24
hotSpotInTAG.nb
nb
Definition: hotSpotInTAG.py:164
ZDCNLCalibration::m_testCalibEnergyHist
std::array< TH1D *, 4 > m_testCalibEnergyHist
Definition: ZDCNLCalibration.h:117
ZDCNLCalibration::m_HEFraction
const std::vector< float > m_HEFraction
Definition: ZDCNLCalibration.h:67
ZDCNLCalibration::L1_ZDC_A_C
Bool_t L1_ZDC_A_C
Definition: ZDCNLCalibration.h:84
ZDCNLCalibration::L1_ZDC_C
Bool_t L1_ZDC_C
Definition: ZDCNLCalibration.h:82
ZDCNLCalibration::m_SNEnergy
const float m_SNEnergy
Definition: ZDCNLCalibration.h:66
ZDCNLCalibration::zdc_ZdcModuleMask
UInt_t zdc_ZdcModuleMask
Definition: ZDCNLCalibration.h:78