|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
|
#include "CalibTracker/SiStripAPVAnalysis/interface/TT6NoiseCalculator.h"
#include "FWCore/MessageLogger/interface/MessageLogger.h"
#include <cmath>
#include <numeric>
#include <algorithm>
using namespace std;
//
// Constructors
//
TT6NoiseCalculator::TT6NoiseCalculator() : numberOfEvents(0), alreadyUsedEvent(false) {
if (false)
cout << "Constructing TT6NoiseCalculator " << endl;
init();
}
//
TT6NoiseCalculator::TT6NoiseCalculator(int evnt_ini, int evnt_iter, float sig_cut)
: numberOfEvents(0), alreadyUsedEvent(false) {
if (false)
cout << "Constructing TT6NoiseCalculator " << endl;
eventsRequiredToCalibrate_ = evnt_ini;
eventsRequiredToUpdate_ = evnt_iter;
cutToAvoidSignal_ = sig_cut;
init();
}
//
// Initialization.
//
void TT6NoiseCalculator::init() {
theCMPSubtractedSignal.clear();
theNoise.clear();
theNoiseSum.clear();
theNoiseSqSum.clear();
theEventPerStrip.clear();
theStatus.setCalibrating();
}
//
// Destructor
//
TT6NoiseCalculator::~TT6NoiseCalculator() {
if (false)
cout << "Destructing TT6NoiseCalculator " << endl;
}
//
// Update the Status of Noise Calculation
//
void TT6NoiseCalculator::updateStatus() {
if (theStatus.isCalibrating() && numberOfEvents >= eventsRequiredToCalibrate_) {
theStatus.setUpdating();
}
}
//
// Calculate and update (when needed) Noise Values
//
void TT6NoiseCalculator::updateNoise(ApvAnalysis::PedestalType& in) {
if (alreadyUsedEvent == false) {
alreadyUsedEvent = true;
numberOfEvents++;
if (numberOfEvents == 1 && theNoise.size() != in.size()) {
edm::LogError("TT6NoiseCalculator:updateNoise")
<< " You did not initialize the Noise correctly prior to noise calibration.";
}
// Initialize sums used for estimating noise.
if ((theStatus.isCalibrating() && numberOfEvents == 1) ||
(theStatus.isUpdating() && (numberOfEvents - eventsRequiredToCalibrate_) % eventsRequiredToUpdate_ == 1)) {
theNoiseSum.clear();
theNoiseSqSum.clear();
theEventPerStrip.clear();
theNoiseSum.resize(in.size(), 0.0);
theNoiseSqSum.resize(in.size(), 0.0);
theEventPerStrip.resize(in.size(), 0);
}
unsigned int i;
// Update sums used for estimating noise.
for (i = 0; i < in.size(); i++) {
if (fabs(in[i]) < cutToAvoidSignal_ * theNoise[i]) {
theNoiseSum[i] += in[i];
theNoiseSqSum[i] += in[i] * in[i];
theEventPerStrip[i]++;
}
}
// Calculate noise.
if ((theStatus.isCalibrating() && numberOfEvents == eventsRequiredToCalibrate_) ||
(theStatus.isUpdating() && (numberOfEvents - eventsRequiredToCalibrate_) % eventsRequiredToUpdate_ == 0)) {
theCMPSubtractedSignal.clear();
theNoise.clear();
for (i = 0; i < in.size(); i++) {
double avVal = (theEventPerStrip[i]) ? theNoiseSum[i] / (theEventPerStrip[i]) : 0.0;
double sqAvVal = (theEventPerStrip[i]) ? theNoiseSqSum[i] / (theEventPerStrip[i]) : 0.0;
double corr_fac = (theEventPerStrip[i] > 1) ? (theEventPerStrip[i] / (theEventPerStrip[i] - 1)) : 1.0;
double rmsVal = (sqAvVal - avVal * avVal > 0.0) ? sqrt(corr_fac * (sqAvVal - avVal * avVal)) : 0.0;
theCMPSubtractedSignal.push_back(static_cast<float>(avVal));
theNoise.push_back(static_cast<float>(rmsVal));
if (false)
cout << " TT6NoiseCalculator::updateNoise " << theNoiseSum[i] << " " << theNoiseSqSum[i] << " "
<< theEventPerStrip[i] << " " << avVal << " " << sqAvVal << " " << (sqAvVal - avVal * avVal) << " "
<< rmsVal << endl;
}
}
updateStatus();
}
}
//
// Define New Event
//
void TT6NoiseCalculator::newEvent() { alreadyUsedEvent = false; }
|