00001 #include "CalibTracker/SiStripAPVAnalysis/interface/SimpleNoiseCalculator.h" 00002 #include "FWCore/MessageLogger/interface/MessageLogger.h" 00003 #include <cmath> 00004 #include <numeric> 00005 #include <algorithm> 00006 using namespace std; 00007 // 00008 // Constructors 00009 // 00010 SimpleNoiseCalculator::SimpleNoiseCalculator() : 00011 numberOfEvents(0) , 00012 alreadyUsedEvent(false) 00013 { 00014 if (0) cout << "Constructing SimpleNoiseCalculator " << endl; 00015 init(); 00016 } 00017 // 00018 SimpleNoiseCalculator::SimpleNoiseCalculator(int evnt_ini, bool use_DB) : 00019 numberOfEvents(0) , 00020 alreadyUsedEvent(false) 00021 { 00022 if (0) cout << "Constructing SimpleNoiseCalculator " << endl; 00023 useDB_ = use_DB; 00024 eventsRequiredToCalibrate_ = evnt_ini; 00025 // eventsRequiredToUpdate_ = evnt_iter; 00026 // cutToAvoidSignal_ = sig_cut; 00027 init(); 00028 } 00029 // 00030 // Initialization. 00031 // 00032 void SimpleNoiseCalculator::init() { 00033 theCMPSubtractedSignal.clear(); 00034 theNoise.clear(); 00035 theNoiseSum.clear(); 00036 theNoiseSqSum.clear(); 00037 theEventPerStrip.clear(); 00038 // theStatus.setCalibrating(); 00039 } 00040 // 00041 // Destructor 00042 // 00043 SimpleNoiseCalculator::~SimpleNoiseCalculator() { 00044 if (0) cout << "Destructing SimpleNoiseCalculator " << endl; 00045 } 00046 // 00047 // Update the Status of Noise Calculation 00048 // 00049 void SimpleNoiseCalculator::updateStatus(){ 00050 if ( (theStatus.isCalibrating() && numberOfEvents >= eventsRequiredToCalibrate_) || (useDB_==true && numberOfEvents ==1) ) { 00051 theStatus.setUpdating(); 00052 } 00053 } 00054 // 00055 // Calculate and update (when needed) Noise Values 00056 // 00057 void SimpleNoiseCalculator::updateNoise(ApvAnalysis::PedestalType& in){ 00058 if (alreadyUsedEvent == false) { 00059 alreadyUsedEvent = true; 00060 numberOfEvents++; 00061 00062 if (numberOfEvents == 1 && theNoise.size() != in.size()) { 00063 edm::LogError("SimpleNoiseCalculator:updateNoise") << " You did not initialize the Noise correctly prior to noise calibration."; 00064 } 00065 00066 // Initialize sums used for estimating noise. 00067 if (numberOfEvents == 1) 00068 { 00069 theNoiseSum.clear(); 00070 theNoiseSqSum.clear(); 00071 theEventPerStrip.clear(); 00072 00073 theNoiseSum.resize(in.size(),0.0); 00074 theNoiseSqSum.resize(in.size(),0.0); 00075 theEventPerStrip.resize(in.size(),0); 00076 } 00077 00078 unsigned int i; 00079 00080 // At every event Update sums used for estimating noise. 00081 for (i = 0; i < in.size(); i++) { 00082 00083 theNoiseSum[i] += in[i]; 00084 theNoiseSqSum[i] += in[i]*in[i]; 00085 theEventPerStrip[i]++; 00086 } 00087 00088 // Calculate noise. 00089 if ((theStatus.isCalibrating() && numberOfEvents == eventsRequiredToCalibrate_) || theStatus.isUpdating() ) 00090 { 00091 theCMPSubtractedSignal.clear(); 00092 theNoise.clear(); 00093 00094 for (i = 0; i < in.size(); i++) { 00095 double avVal = (theEventPerStrip[i]) 00096 ? theNoiseSum[i]/(theEventPerStrip[i]):0.0; 00097 double sqAvVal = (theEventPerStrip[i]) 00098 ? theNoiseSqSum[i]/(theEventPerStrip[i]):0.0; 00099 double corr_fac = (theEventPerStrip[i] > 1) 00100 ? (theEventPerStrip[i]/(theEventPerStrip[i]-1)) : 1.0; 00101 double rmsVal = (sqAvVal - avVal*avVal > 0.0) 00102 ? sqrt(corr_fac * (sqAvVal - avVal*avVal)) : 0.0; 00103 00104 theCMPSubtractedSignal.push_back(static_cast<float>(avVal)); 00105 00106 theNoise.push_back(static_cast<float>(rmsVal)); 00107 00108 if (0) cout << " SimpleNoiseCalculator::updateNoise " 00109 << theNoiseSum[i] << " " 00110 << theNoiseSqSum[i] << " " 00111 << theEventPerStrip[i] << " " 00112 << avVal << " " 00113 << sqAvVal << " " 00114 << (sqAvVal - avVal*avVal) << " " 00115 << rmsVal << endl; 00116 } 00117 } 00118 updateStatus(); 00119 } 00120 } 00121 // 00122 // Define New Event 00123 // 00124 void SimpleNoiseCalculator::newEvent() { 00125 alreadyUsedEvent = false; 00126 } 00127 00128