1 |
flucke |
1.1.2.1 |
/// \class SiStripBackplaneCalibration
|
2 |
|
|
///
|
3 |
flucke |
1.1.2.2 |
/// Calibration of back plane corrections for the strip tracker, integrated in the
|
4 |
|
|
/// alignment algorithms. Note that not all algorithms support this...
|
5 |
flucke |
1.1.2.1 |
///
|
6 |
flucke |
1.1.2.2 |
/// Usally use one instance for deco mode data since peak mode should give the
|
7 |
|
|
// reference - but the strip local reco foresees also calibration parameters
|
8 |
|
|
// for peak data like for the Lorentz angle.
|
9 |
flucke |
1.1.2.1 |
///
|
10 |
|
|
/// \author : Gero Flucke
|
11 |
|
|
/// date : November 2012
|
12 |
flucke |
1.1.2.3 |
/// $Revision: 1.1.2.2 $
|
13 |
|
|
/// $Date: 2013/04/17 09:50:03 $
|
14 |
flucke |
1.1.2.1 |
/// (last update by $Author: flucke $)
|
15 |
|
|
|
16 |
|
|
#include "Alignment/CommonAlignmentAlgorithm/interface/IntegratedCalibrationBase.h"
|
17 |
flucke |
1.1.2.2 |
#include "Alignment/CommonAlignmentAlgorithm/plugins/SiStripReadoutModeEnums.h"
|
18 |
flucke |
1.1.2.1 |
|
19 |
flucke |
1.1.2.2 |
#include "CondCore/DBOutputService/interface/PoolDBOutputService.h"
|
20 |
flucke |
1.1.2.1 |
#include "CondFormats/DataRecord/interface/SiStripCondDataRecords.h"
|
21 |
|
|
#include "CondFormats/SiStripObjects/interface/SiStripLatency.h"
|
22 |
|
|
#include "CondFormats/SiStripObjects/interface/SiStripLorentzAngle.h"
|
23 |
flucke |
1.1.2.2 |
// #include "CalibTracker/Records/interface/SiStripDependentRecords.h"
|
24 |
|
|
#include "CondFormats/SiStripObjects/interface/SiStripBackPlaneCorrection.h"
|
25 |
flucke |
1.1.2.1 |
|
26 |
|
|
#include "DataFormats/DetId/interface/DetId.h"
|
27 |
|
|
#include "DataFormats/SiStripDetId/interface/SiStripDetId.h"
|
28 |
|
|
|
29 |
|
|
#include "FWCore/Framework/interface/ESHandle.h"
|
30 |
|
|
#include "FWCore/Framework/interface/ESWatcher.h"
|
31 |
|
|
#include "FWCore/MessageLogger/interface/MessageLogger.h"
|
32 |
flucke |
1.1.2.2 |
#include "FWCore/ServiceRegistry/interface/Service.h"
|
33 |
flucke |
1.1.2.1 |
#include "FWCore/ParameterSet/interface/ParameterSet.h"
|
34 |
|
|
|
35 |
|
|
#include "Geometry/CommonDetUnit/interface/GeomDet.h"
|
36 |
|
|
#include "MagneticField/Engine/interface/MagneticField.h"
|
37 |
|
|
#include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
|
38 |
|
|
#include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHit.h"
|
39 |
|
|
|
40 |
flucke |
1.1.2.2 |
#include "TTree.h"
|
41 |
|
|
#include "TFile.h"
|
42 |
|
|
#include "TString.h"
|
43 |
flucke |
1.1.2.1 |
|
44 |
|
|
// #include <iostream>
|
45 |
|
|
#include <vector>
|
46 |
|
|
#include <map>
|
47 |
flucke |
1.1.2.2 |
#include <sstream>
|
48 |
|
|
#include <cstdio>
|
49 |
flucke |
1.1.2.3 |
#include <functional>
|
50 |
flucke |
1.1.2.1 |
|
51 |
|
|
class SiStripBackplaneCalibration : public IntegratedCalibrationBase
|
52 |
|
|
{
|
53 |
|
|
public:
|
54 |
|
|
/// Constructor
|
55 |
|
|
explicit SiStripBackplaneCalibration(const edm::ParameterSet &cfg);
|
56 |
|
|
|
57 |
|
|
/// Destructor
|
58 |
|
|
virtual ~SiStripBackplaneCalibration();
|
59 |
|
|
|
60 |
|
|
/// How many parameters does this calibration define?
|
61 |
|
|
virtual unsigned int numParameters() const;
|
62 |
|
|
|
63 |
|
|
// /// Return all derivatives,
|
64 |
|
|
// /// default implementation uses other derivatives(..) method,
|
65 |
|
|
// /// but can be overwritten in derived class for efficiency.
|
66 |
|
|
// virtual std::vector<double> derivatives(const TransientTrackingRecHit &hit,
|
67 |
|
|
// const TrajectoryStateOnSurface &tsos,
|
68 |
|
|
// const edm::EventSetup &setup,
|
69 |
|
|
// const EventInfo &eventInfo) const;
|
70 |
|
|
|
71 |
|
|
/// Return non-zero derivatives for x- and y-measurements with their indices by reference.
|
72 |
|
|
/// Return value is their number.
|
73 |
|
|
virtual unsigned int derivatives(std::vector<ValuesIndexPair> &outDerivInds,
|
74 |
|
|
const TransientTrackingRecHit &hit,
|
75 |
|
|
const TrajectoryStateOnSurface &tsos,
|
76 |
|
|
const edm::EventSetup &setup,
|
77 |
|
|
const EventInfo &eventInfo) const;
|
78 |
|
|
|
79 |
|
|
/// Setting the determined parameter identified by index,
|
80 |
|
|
/// returns false if out-of-bounds, true otherwise.
|
81 |
|
|
virtual bool setParameter(unsigned int index, double value);
|
82 |
|
|
|
83 |
|
|
/// Setting the determined parameter uncertainty identified by index,
|
84 |
|
|
/// returns false if out-of-bounds, true otherwise.
|
85 |
|
|
virtual bool setParameterError(unsigned int index, double error);
|
86 |
|
|
|
87 |
|
|
/// Return current value of parameter identified by index.
|
88 |
|
|
/// Returns 0. if index out-of-bounds.
|
89 |
|
|
virtual double getParameter(unsigned int index) const;
|
90 |
|
|
|
91 |
|
|
/// Return current value of parameter identified by index.
|
92 |
|
|
/// Returns 0. if index out-of-bounds or if errors undetermined.
|
93 |
|
|
virtual double getParameterError(unsigned int index) const;
|
94 |
|
|
|
95 |
|
|
// /// Call at beginning of job:
|
96 |
|
|
// virtual void beginOfJob(const AlignableTracker *tracker,
|
97 |
|
|
// const AlignableMuon *muon,
|
98 |
|
|
// const AlignableExtras *extras);
|
99 |
|
|
|
100 |
|
|
/// Called at end of a the job of the AlignmentProducer.
|
101 |
|
|
/// Write out determined parameters.
|
102 |
|
|
virtual void endOfJob();
|
103 |
|
|
|
104 |
|
|
private:
|
105 |
flucke |
1.1.2.2 |
/// If called the first time, fill 'siStripBackPlaneCorrInput_',
|
106 |
|
|
/// later check that Backplane has not changed.
|
107 |
|
|
bool checkBackPlaneCorrectionInput(const edm::EventSetup &setup, const EventInfo &eventInfo);
|
108 |
|
|
/// Input BackPlaneCorrection values:
|
109 |
flucke |
1.1.2.1 |
/// - either from EventSetup of first call to derivatives(..)
|
110 |
flucke |
1.1.2.2 |
/// - or created from files of passed by configuration (i.e. from parallel processing)
|
111 |
|
|
const SiStripBackPlaneCorrection* getBackPlaneCorrectionInput();
|
112 |
flucke |
1.1.2.1 |
|
113 |
|
|
/// Determined parameter value for this detId (detId not treated => 0.)
|
114 |
|
|
/// and the given run.
|
115 |
|
|
double getParameterForDetId(unsigned int detId, edm::RunNumber_t run) const;
|
116 |
|
|
/// Index of parameter for given detId (detId not treated => < 0)
|
117 |
|
|
/// and the given run.
|
118 |
|
|
int getParameterIndexFromDetId(unsigned int detId, edm::RunNumber_t run) const;
|
119 |
|
|
/// Total number of IOVs.
|
120 |
|
|
unsigned int numIovs() const;
|
121 |
|
|
/// First run of iov (0 if iovNum not treated).
|
122 |
|
|
edm::RunNumber_t firstRunOfIOV(unsigned int iovNum) const;
|
123 |
|
|
|
124 |
flucke |
1.1.2.2 |
void writeTree(const SiStripBackPlaneCorrection *backPlaneCorr,
|
125 |
|
|
const std::map<unsigned int,float> &errors, const char *treeName) const;
|
126 |
|
|
SiStripBackPlaneCorrection* createFromTree(const char *fileName, const char *treeName) const;
|
127 |
flucke |
1.1.2.1 |
|
128 |
flucke |
1.1.2.2 |
const std::string readoutModeName_;
|
129 |
|
|
int16_t readoutMode_;
|
130 |
flucke |
1.1.2.1 |
const bool saveToDB_;
|
131 |
|
|
const std::string recordNameDBwrite_;
|
132 |
|
|
const std::string outFileName_;
|
133 |
|
|
const std::vector<std::string> mergeFileNames_;
|
134 |
|
|
|
135 |
flucke |
1.1.2.2 |
edm::ESWatcher<SiStripBackPlaneCorrectionRcd> watchBackPlaneCorrRcd_;
|
136 |
flucke |
1.1.2.1 |
|
137 |
|
|
// const AlignableTracker *alignableTracker_;
|
138 |
flucke |
1.1.2.2 |
SiStripBackPlaneCorrection *siStripBackPlaneCorrInput_;
|
139 |
flucke |
1.1.2.1 |
|
140 |
|
|
std::vector<double> parameters_;
|
141 |
|
|
std::vector<double> paramUncertainties_;
|
142 |
|
|
};
|
143 |
|
|
|
144 |
|
|
//======================================================================
|
145 |
|
|
//======================================================================
|
146 |
|
|
//======================================================================
|
147 |
|
|
|
148 |
|
|
SiStripBackplaneCalibration::SiStripBackplaneCalibration(const edm::ParameterSet &cfg)
|
149 |
|
|
: IntegratedCalibrationBase(cfg),
|
150 |
flucke |
1.1.2.2 |
readoutModeName_(cfg.getParameter<std::string>("readoutMode")),
|
151 |
flucke |
1.1.2.1 |
saveToDB_(cfg.getParameter<bool>("saveToDB")),
|
152 |
|
|
recordNameDBwrite_(cfg.getParameter<std::string>("recordNameDBwrite")),
|
153 |
|
|
outFileName_(cfg.getParameter<std::string>("treeFile")),
|
154 |
|
|
mergeFileNames_(cfg.getParameter<std::vector<std::string> >("mergeTreeFiles")),
|
155 |
|
|
// alignableTracker_(0),
|
156 |
flucke |
1.1.2.2 |
siStripBackPlaneCorrInput_(0)
|
157 |
flucke |
1.1.2.1 |
{
|
158 |
flucke |
1.1.2.2 |
// SiStripLatency::singleReadOutMode() returns
|
159 |
|
|
// 1: all in peak, 0: all in deco, -1: mixed state
|
160 |
|
|
// (in principle one could treat even mixed state APV by APV...)
|
161 |
|
|
if (readoutModeName_ == "peak") {
|
162 |
|
|
readoutMode_ = kPeakMode;
|
163 |
|
|
} else if (readoutModeName_ == "deconvolution") {
|
164 |
|
|
readoutMode_ = kDeconvolutionMode;
|
165 |
|
|
} else {
|
166 |
|
|
throw cms::Exception("BadConfig")
|
167 |
|
|
<< "SiStripBackplaneCalibration:\n" << "Unknown mode '"
|
168 |
|
|
<< readoutModeName_ << "', should be 'peak' or 'deconvolution' .\n";
|
169 |
|
|
}
|
170 |
|
|
|
171 |
flucke |
1.1.2.1 |
// FIXME: Which granularity, leading to how many parameters?
|
172 |
flucke |
1.1.2.2 |
parameters_.resize(2, 0.); // currently two parameters (TIB, TOB), start value 0.
|
173 |
|
|
paramUncertainties_.resize(2, 0.); // dito for errors
|
174 |
flucke |
1.1.2.1 |
|
175 |
|
|
edm::LogInfo("Alignment") << "@SUB=SiStripBackplaneCalibration" << "Created with name "
|
176 |
flucke |
1.1.2.2 |
<< this->name() << " for readout mode '" << readoutModeName_
|
177 |
flucke |
1.1.2.1 |
<< "',\n" << this->numParameters() << " parameters to be determined."
|
178 |
|
|
<< "\nsaveToDB = " << saveToDB_
|
179 |
|
|
<< "\n outFileName = " << outFileName_
|
180 |
|
|
<< "\n N(merge files) = " << mergeFileNames_.size();
|
181 |
|
|
if (mergeFileNames_.size()) {
|
182 |
|
|
edm::LogInfo("Alignment") << "@SUB=SiStripBackplaneCalibration"
|
183 |
|
|
<< "First file to merge: " << mergeFileNames_[0];
|
184 |
|
|
}
|
185 |
|
|
}
|
186 |
|
|
|
187 |
|
|
//======================================================================
|
188 |
|
|
SiStripBackplaneCalibration::~SiStripBackplaneCalibration()
|
189 |
|
|
{
|
190 |
flucke |
1.1.2.2 |
// std::cout << "Destroy SiStripBackplaneCalibration named " << this->name() << std::endl;
|
191 |
|
|
delete siStripBackPlaneCorrInput_;
|
192 |
flucke |
1.1.2.1 |
}
|
193 |
|
|
|
194 |
|
|
//======================================================================
|
195 |
|
|
unsigned int SiStripBackplaneCalibration::numParameters() const
|
196 |
|
|
{
|
197 |
|
|
return parameters_.size();
|
198 |
|
|
}
|
199 |
|
|
|
200 |
|
|
//======================================================================
|
201 |
|
|
unsigned int
|
202 |
|
|
SiStripBackplaneCalibration::derivatives(std::vector<ValuesIndexPair> &outDerivInds,
|
203 |
|
|
const TransientTrackingRecHit &hit,
|
204 |
|
|
const TrajectoryStateOnSurface &tsos,
|
205 |
|
|
const edm::EventSetup &setup,
|
206 |
|
|
const EventInfo &eventInfo) const
|
207 |
|
|
{
|
208 |
|
|
// ugly const-cast:
|
209 |
|
|
// But it is either only first initialisation or throwing an exception...
|
210 |
flucke |
1.1.2.2 |
const_cast<SiStripBackplaneCalibration*>(this)->checkBackPlaneCorrectionInput(setup, eventInfo);
|
211 |
flucke |
1.1.2.1 |
|
212 |
|
|
outDerivInds.clear();
|
213 |
|
|
|
214 |
|
|
edm::ESHandle<SiStripLatency> latency;
|
215 |
|
|
setup.get<SiStripLatencyRcd>().get(latency);
|
216 |
|
|
const int16_t mode = latency->singleReadOutMode();
|
217 |
flucke |
1.1.2.3 |
// std::cout << "SiStripBackplaneCalibration in mode '" << readoutModeName_ << "' finds mode "
|
218 |
|
|
// << mode << std::endl;
|
219 |
|
|
//
|
220 |
|
|
// std::vector<SiStripLatency::Latency> latencies = const_cast<SiStripLatency*>(latency.product())->allUniqueLatencyAndModes();
|
221 |
|
|
// for (auto i = latencies.begin(); i != latencies.end(); ++i) {
|
222 |
|
|
// std::cout << static_cast<int>(i->latency) << ", mode "
|
223 |
|
|
// << static_cast<int>(i->mode) << ", id_apv " << i->detIdAndApv
|
224 |
|
|
// << std::endl;
|
225 |
|
|
// }
|
226 |
|
|
|
227 |
flucke |
1.1.2.2 |
if(mode == readoutMode_) {
|
228 |
flucke |
1.1.2.1 |
if (hit.det()) { // otherwise 'constraint hit' or whatever
|
229 |
flucke |
1.1.2.2 |
|
230 |
flucke |
1.1.2.1 |
const int index = this->getParameterIndexFromDetId(hit.det()->geographicalId(),
|
231 |
|
|
eventInfo.eventId_.run());
|
232 |
|
|
if (index >= 0) { // otherwise not treated
|
233 |
|
|
edm::ESHandle<MagneticField> magneticField;
|
234 |
|
|
setup.get<IdealMagneticFieldRecord>().get(magneticField);
|
235 |
|
|
const GlobalVector bField(magneticField->inTesla(hit.det()->surface().position()));
|
236 |
|
|
const LocalVector bFieldLocal(hit.det()->surface().toLocal(bField));
|
237 |
flucke |
1.1.2.2 |
//std::cout << "SiStripBackplaneCalibration derivatives " << readoutModeName_ << std::endl;
|
238 |
flucke |
1.1.2.1 |
const double dZ = hit.det()->surface().bounds().thickness(); // it's a float only...
|
239 |
|
|
const double tanPsi = tsos.localParameters().mixedFormatVector()[1]; //float...
|
240 |
|
|
|
241 |
|
|
edm::ESHandle<SiStripLorentzAngle> lorentzAngleHandle;
|
242 |
flucke |
1.1.2.2 |
setup.get<SiStripLorentzAngleRcd>().get(readoutModeName_, lorentzAngleHandle);
|
243 |
flucke |
1.1.2.1 |
// Yes, mobility (= LA/By) stored in object called LA...
|
244 |
|
|
const double mobility = lorentzAngleHandle->getLorentzAngle(hit.det()->geographicalId());
|
245 |
|
|
// shift due to dead back plane has two parts:
|
246 |
|
|
// 1) Lorentz Angle correction formula gets reduced thickness dz*(1-bp)
|
247 |
|
|
// 2) 'Direct' effect is shift of effective module position in local z by bp*dz/2
|
248 |
|
|
// (see GF's presentation in alignment meeting 25.10.2012,
|
249 |
|
|
// https://indico.cern.ch/conferenceDisplay.py?confId=174266#2012-10-25)
|
250 |
|
|
const double xDerivative = 0.5 * dZ * (mobility * bFieldLocal.y() - tanPsi);
|
251 |
|
|
// std::cout << "derivative is " << xDerivative << " for index " << index
|
252 |
|
|
// << std::endl;
|
253 |
|
|
const Values derivs(xDerivative, 0.); // yDerivative = 0.
|
254 |
|
|
outDerivInds.push_back(ValuesIndexPair(derivs, index));
|
255 |
flucke |
1.1.2.2 |
}
|
256 |
flucke |
1.1.2.1 |
} else {
|
257 |
flucke |
1.1.2.2 |
edm::LogWarning("Alignment") << "@SUB=SiStripBackplaneCalibration::derivatives1"
|
258 |
flucke |
1.1.2.1 |
<< "Hit without GeomDet, skip!";
|
259 |
|
|
}
|
260 |
flucke |
1.1.2.2 |
} else if (mode != kDeconvolutionMode && mode != kPeakMode) {
|
261 |
|
|
// warn only if unknown/mixed mode
|
262 |
flucke |
1.1.2.1 |
edm::LogWarning("Alignment") << "@SUB=SiStripBackplaneCalibration::derivatives2"
|
263 |
flucke |
1.1.2.2 |
<< "Readout mode is " << mode << ", but looking for "
|
264 |
|
|
<< readoutMode_ << " (" << readoutModeName_ << ").";
|
265 |
flucke |
1.1.2.1 |
}
|
266 |
|
|
|
267 |
|
|
return outDerivInds.size();
|
268 |
|
|
}
|
269 |
|
|
|
270 |
|
|
//======================================================================
|
271 |
|
|
bool SiStripBackplaneCalibration::setParameter(unsigned int index, double value)
|
272 |
|
|
{
|
273 |
|
|
if (index >= parameters_.size()) {
|
274 |
|
|
return false;
|
275 |
|
|
} else {
|
276 |
|
|
parameters_[index] = value;
|
277 |
|
|
return true;
|
278 |
|
|
}
|
279 |
|
|
}
|
280 |
|
|
|
281 |
|
|
//======================================================================
|
282 |
|
|
bool SiStripBackplaneCalibration::setParameterError(unsigned int index, double error)
|
283 |
|
|
{
|
284 |
|
|
if (index >= paramUncertainties_.size()) {
|
285 |
|
|
return false;
|
286 |
|
|
} else {
|
287 |
|
|
paramUncertainties_[index] = error;
|
288 |
|
|
return true;
|
289 |
|
|
}
|
290 |
|
|
}
|
291 |
|
|
|
292 |
|
|
//======================================================================
|
293 |
|
|
double SiStripBackplaneCalibration::getParameter(unsigned int index) const
|
294 |
|
|
{
|
295 |
|
|
// if (index >= parameters_.size()) {
|
296 |
|
|
// return 0.;
|
297 |
|
|
// } else {
|
298 |
|
|
// return parameters_[index];
|
299 |
|
|
// }
|
300 |
|
|
return (index >= parameters_.size() ? 0. : parameters_[index]);
|
301 |
|
|
}
|
302 |
|
|
|
303 |
|
|
//======================================================================
|
304 |
|
|
double SiStripBackplaneCalibration::getParameterError(unsigned int index) const
|
305 |
|
|
{
|
306 |
|
|
// if (index >= paramUncertainties_.size()) {
|
307 |
|
|
// return 0.;
|
308 |
|
|
// } else {
|
309 |
|
|
// return paramUncertainties_[index];
|
310 |
|
|
// }
|
311 |
|
|
return (index >= paramUncertainties_.size() ? 0. : paramUncertainties_[index]);
|
312 |
|
|
}
|
313 |
|
|
|
314 |
|
|
// //======================================================================
|
315 |
|
|
// void SiStripBackplaneCalibration::beginOfJob(const AlignableTracker *tracker,
|
316 |
flucke |
1.1.2.2 |
// const AlignableMuon */*muon*/,
|
317 |
|
|
// const AlignableExtras */*extras*/)
|
318 |
flucke |
1.1.2.1 |
// {
|
319 |
|
|
// alignableTracker_ = tracker;
|
320 |
|
|
// }
|
321 |
|
|
|
322 |
|
|
|
323 |
|
|
//======================================================================
|
324 |
|
|
void SiStripBackplaneCalibration::endOfJob()
|
325 |
|
|
{
|
326 |
|
|
// loginfo output
|
327 |
|
|
std::ostringstream out;
|
328 |
flucke |
1.1.2.2 |
out << "Parameter results for readout mode '" << readoutModeName_ << "'\n";
|
329 |
flucke |
1.1.2.1 |
for (unsigned int iPar = 0; iPar < parameters_.size(); ++iPar) {
|
330 |
|
|
out << iPar << ": " << parameters_[iPar] << " +- " << paramUncertainties_[iPar] << "\n";
|
331 |
|
|
}
|
332 |
|
|
edm::LogInfo("Alignment") << "@SUB=SiStripBackplaneCalibration::endOfJob" << out.str();
|
333 |
|
|
|
334 |
flucke |
1.1.2.2 |
std::map<unsigned int, float> errors; // Array of errors for each detId
|
335 |
flucke |
1.1.2.1 |
|
336 |
flucke |
1.1.2.2 |
// now write 'input' tree
|
337 |
|
|
const SiStripBackPlaneCorrection *input = this->getBackPlaneCorrectionInput(); // never NULL
|
338 |
|
|
const std::string treeName(this->name() + '_' + readoutModeName_ + '_');
|
339 |
|
|
this->writeTree(input, errors, (treeName + "input").c_str()); // empty errors for input...
|
340 |
|
|
|
341 |
|
|
if (input->getBackPlaneCorrections().empty()) {
|
342 |
|
|
edm::LogError("Alignment") << "@SUB=SiStripBackplaneCalibration::endOfJob"
|
343 |
|
|
<< "Input back plane correction map is empty ('"
|
344 |
|
|
<< readoutModeName_ << "' mode), skip writing output!";
|
345 |
|
|
return;
|
346 |
flucke |
1.1.2.1 |
}
|
347 |
|
|
|
348 |
flucke |
1.1.2.3 |
const unsigned int nonZeroParamsOrErrors = // Any determined value?
|
349 |
|
|
count_if (parameters_.begin(), parameters_.end(), std::bind2nd(std::not_equal_to<double>(),0.))
|
350 |
|
|
+ count_if(paramUncertainties_.begin(), paramUncertainties_.end(),
|
351 |
|
|
std::bind2nd(std::not_equal_to<double>(), 0.));
|
352 |
|
|
|
353 |
flucke |
1.1.2.2 |
for (unsigned int iIOV = 0; iIOV < this->numIovs(); ++iIOV) {
|
354 |
|
|
cond::Time_t firstRunOfIOV = this->firstRunOfIOV(iIOV);
|
355 |
|
|
SiStripBackPlaneCorrection *output = new SiStripBackPlaneCorrection;
|
356 |
|
|
// Loop on map of values from input and add (possible) parameter results
|
357 |
|
|
for (auto iterIdValue = input->getBackPlaneCorrections().begin();
|
358 |
|
|
iterIdValue != input->getBackPlaneCorrections().end(); ++iterIdValue) {
|
359 |
|
|
// type of (*iterIdValue) is pair<unsigned int, float>
|
360 |
|
|
const unsigned int detId = iterIdValue->first; // key of map is DetId
|
361 |
|
|
const float value = iterIdValue->second + this->getParameterForDetId(detId, firstRunOfIOV);
|
362 |
|
|
output->putBackPlaneCorrection(detId, value); // put result in output
|
363 |
|
|
int parameterIndex = this->getParameterIndexFromDetId(detId, firstRunOfIOV);
|
364 |
|
|
errors[detId] = this->getParameterError(parameterIndex);
|
365 |
|
|
}
|
366 |
flucke |
1.1.2.1 |
|
367 |
flucke |
1.1.2.3 |
if (saveToDB_ || nonZeroParamsOrErrors != 0) { // Skip writing mille jobs...
|
368 |
|
|
this->writeTree(output, errors, (treeName + Form("result_%lld", firstRunOfIOV)).c_str());
|
369 |
|
|
}
|
370 |
flucke |
1.1.2.1 |
|
371 |
flucke |
1.1.2.2 |
if (saveToDB_) { // If requested, write out to DB
|
372 |
|
|
edm::Service<cond::service::PoolDBOutputService> dbService;
|
373 |
|
|
if (dbService.isAvailable()) {
|
374 |
|
|
dbService->writeOne(output, firstRunOfIOV, recordNameDBwrite_.c_str());
|
375 |
|
|
// no 'delete output;': writeOne(..) took over ownership
|
376 |
|
|
} else {
|
377 |
|
|
delete output;
|
378 |
|
|
edm::LogError("BadConfig") << "@SUB=SiStripBackplaneCalibration::endOfJob"
|
379 |
|
|
<< "No PoolDBOutputService available, but saveToDB true!";
|
380 |
|
|
}
|
381 |
|
|
} else {
|
382 |
|
|
delete output;
|
383 |
|
|
}
|
384 |
|
|
} // end loop on IOVs
|
385 |
flucke |
1.1.2.1 |
}
|
386 |
|
|
|
387 |
|
|
//======================================================================
|
388 |
flucke |
1.1.2.2 |
bool SiStripBackplaneCalibration::checkBackPlaneCorrectionInput(const edm::EventSetup &setup,
|
389 |
|
|
const EventInfo &eventInfo)
|
390 |
flucke |
1.1.2.1 |
{
|
391 |
flucke |
1.1.2.2 |
edm::ESHandle<SiStripBackPlaneCorrection> backPlaneCorrHandle;
|
392 |
|
|
if (!siStripBackPlaneCorrInput_) {
|
393 |
|
|
setup.get<SiStripBackPlaneCorrectionRcd>().get(readoutModeName_, backPlaneCorrHandle);
|
394 |
|
|
siStripBackPlaneCorrInput_ = new SiStripBackPlaneCorrection(*backPlaneCorrHandle);
|
395 |
|
|
// FIXME: Should we call 'watchBackPlaneCorrRcd_.check(setup)' as well?
|
396 |
|
|
// Otherwise could be that next check has to check via following 'else', though
|
397 |
|
|
// no new IOV has started... (to be checked)
|
398 |
|
|
} else {
|
399 |
|
|
if (watchBackPlaneCorrRcd_.check(setup)) { // new IOV of input - but how to check peak vs deco?
|
400 |
|
|
setup.get<SiStripBackPlaneCorrectionRcd>().get(readoutModeName_, backPlaneCorrHandle);
|
401 |
|
|
if (backPlaneCorrHandle->getBackPlaneCorrections() // but only bad if non-identical values
|
402 |
|
|
!= siStripBackPlaneCorrInput_->getBackPlaneCorrections()) { // (comparing maps)
|
403 |
|
|
// Maps are containers sorted by key, but comparison problems may arise from
|
404 |
|
|
// 'floating point comparison' problems (FIXME?)
|
405 |
|
|
throw cms::Exception("BadInput")
|
406 |
|
|
<< "SiStripBackplaneCalibration::checkBackPlaneCorrectionInput:\n"
|
407 |
|
|
<< "Content of SiStripBackPlaneCorrection changed at run " << eventInfo.eventId_.run()
|
408 |
|
|
<< ", but algorithm expects constant input!\n";
|
409 |
|
|
return false; // not reached...
|
410 |
|
|
}
|
411 |
flucke |
1.1.2.1 |
}
|
412 |
|
|
}
|
413 |
|
|
|
414 |
|
|
return true;
|
415 |
|
|
}
|
416 |
|
|
|
417 |
|
|
//======================================================================
|
418 |
flucke |
1.1.2.2 |
const SiStripBackPlaneCorrection* SiStripBackplaneCalibration::getBackPlaneCorrectionInput()
|
419 |
flucke |
1.1.2.1 |
{
|
420 |
flucke |
1.1.2.2 |
// For parallel processing in Millepede II, create SiStripBackPlaneCorrection
|
421 |
|
|
// from info stored in files of parallel jobs and check that they are identical.
|
422 |
|
|
// If this job has run on events, still check that back plane corrections are
|
423 |
|
|
// identical to the ones from mergeFileNames_.
|
424 |
|
|
const std::string treeName(((this->name() + '_') += readoutModeName_) += "_input");
|
425 |
|
|
for (auto iFile = mergeFileNames_.begin(); iFile != mergeFileNames_.end(); ++iFile) {
|
426 |
|
|
SiStripBackPlaneCorrection* bpCorr = this->createFromTree(iFile->c_str(), treeName.c_str());
|
427 |
|
|
// siStripBackPlaneCorrInput_ could be non-null from previous file of this loop
|
428 |
|
|
// or from checkBackPlaneCorrectionInput(..) when running on data in this job as well
|
429 |
|
|
if (!siStripBackPlaneCorrInput_ || siStripBackPlaneCorrInput_->getBackPlaneCorrections().empty()) {
|
430 |
|
|
delete siStripBackPlaneCorrInput_; // NULL or empty
|
431 |
|
|
siStripBackPlaneCorrInput_ = bpCorr;
|
432 |
|
|
} else {
|
433 |
|
|
// about comparison of maps see comments in checkBackPlaneCorrectionInput
|
434 |
|
|
if (bpCorr && !bpCorr->getBackPlaneCorrections().empty() && // single job might not have got events
|
435 |
|
|
bpCorr->getBackPlaneCorrections() != siStripBackPlaneCorrInput_->getBackPlaneCorrections()) {
|
436 |
|
|
// Throw exception instead of error?
|
437 |
|
|
edm::LogError("NoInput") << "@SUB=SiStripBackplaneCalibration::getBackPlaneCorrectionInput"
|
438 |
|
|
<< "Different input values from tree " << treeName
|
439 |
|
|
<< " in file " << *iFile << ".";
|
440 |
flucke |
1.1.2.1 |
|
441 |
flucke |
1.1.2.2 |
}
|
442 |
|
|
delete bpCorr;
|
443 |
|
|
}
|
444 |
|
|
}
|
445 |
flucke |
1.1.2.1 |
|
446 |
flucke |
1.1.2.2 |
if (!siStripBackPlaneCorrInput_) { // no files nor ran on events
|
447 |
|
|
siStripBackPlaneCorrInput_ = new SiStripBackPlaneCorrection;
|
448 |
|
|
edm::LogError("NoInput") << "@SUB=SiStripBackplaneCalibration::getBackPlaneCorrectionInput"
|
449 |
|
|
<< "No input, create an empty one ('" << readoutModeName_ << "' mode)!";
|
450 |
|
|
} else if (siStripBackPlaneCorrInput_->getBackPlaneCorrections().empty()) {
|
451 |
|
|
edm::LogError("NoInput") << "@SUB=SiStripBackplaneCalibration::getBackPlaneCorrectionInput"
|
452 |
|
|
<< "Empty result ('" << readoutModeName_ << "' mode)!";
|
453 |
|
|
}
|
454 |
|
|
|
455 |
|
|
return siStripBackPlaneCorrInput_;
|
456 |
flucke |
1.1.2.1 |
}
|
457 |
|
|
|
458 |
|
|
//======================================================================
|
459 |
|
|
double SiStripBackplaneCalibration::getParameterForDetId(unsigned int detId,
|
460 |
|
|
edm::RunNumber_t run) const
|
461 |
|
|
{
|
462 |
|
|
const int index = this->getParameterIndexFromDetId(detId, run);
|
463 |
|
|
|
464 |
|
|
return (index < 0 ? 0. : parameters_[index]);
|
465 |
|
|
}
|
466 |
|
|
|
467 |
|
|
//======================================================================
|
468 |
|
|
int SiStripBackplaneCalibration::getParameterIndexFromDetId(unsigned int detId,
|
469 |
|
|
edm::RunNumber_t run) const
|
470 |
|
|
{
|
471 |
|
|
// Return the index of the parameter that is used for this DetId.
|
472 |
|
|
// If this DetId is not treated, return values < 0.
|
473 |
|
|
|
474 |
|
|
// FIXME: Extend to configurable granularity?
|
475 |
|
|
// Including treatment of run dependence?
|
476 |
|
|
const SiStripDetId id(detId);
|
477 |
|
|
if (id.det() == DetId::Tracker) {
|
478 |
flucke |
1.1.2.2 |
if (id.subDetector() == SiStripDetId::TIB) return 0;
|
479 |
|
|
else if (id.subDetector() == SiStripDetId::TOB) return 1;
|
480 |
flucke |
1.1.2.1 |
}
|
481 |
|
|
|
482 |
|
|
return -1;
|
483 |
|
|
}
|
484 |
|
|
|
485 |
|
|
//======================================================================
|
486 |
|
|
unsigned int SiStripBackplaneCalibration::numIovs() const
|
487 |
|
|
{
|
488 |
|
|
// FIXME: Needed to include treatment of run dependence!
|
489 |
|
|
return 1;
|
490 |
|
|
}
|
491 |
|
|
|
492 |
|
|
//======================================================================
|
493 |
|
|
edm::RunNumber_t SiStripBackplaneCalibration::firstRunOfIOV(unsigned int iovNum) const
|
494 |
|
|
{
|
495 |
|
|
// FIXME: Needed to include treatment of run dependence!
|
496 |
|
|
if (iovNum < this->numIovs()) return 1;
|
497 |
|
|
else return 0;
|
498 |
|
|
}
|
499 |
|
|
|
500 |
|
|
|
501 |
flucke |
1.1.2.2 |
//======================================================================
|
502 |
|
|
void SiStripBackplaneCalibration::writeTree(const SiStripBackPlaneCorrection *backPlaneCorrection,
|
503 |
|
|
const std::map<unsigned int, float> &errors,
|
504 |
|
|
const char *treeName) const
|
505 |
|
|
{
|
506 |
|
|
if (!backPlaneCorrection) return;
|
507 |
|
|
|
508 |
|
|
TFile* file = TFile::Open(outFileName_.c_str(), "UPDATE");
|
509 |
|
|
if (!file) {
|
510 |
|
|
edm::LogError("BadConfig") << "@SUB=SiStripBackplaneCalibration::writeTree"
|
511 |
|
|
<< "Could not open file '" << outFileName_ << "'.";
|
512 |
|
|
return;
|
513 |
|
|
}
|
514 |
flucke |
1.1.2.1 |
|
515 |
flucke |
1.1.2.2 |
TTree *tree = new TTree(treeName, treeName);
|
516 |
|
|
unsigned int id = 0;
|
517 |
|
|
float value = 0.;
|
518 |
|
|
float error = 0.;
|
519 |
|
|
tree->Branch("detId", &id, "detId/i");
|
520 |
|
|
tree->Branch("value", &value, "value/F");
|
521 |
|
|
tree->Branch("error", &error, "error/F");
|
522 |
|
|
|
523 |
|
|
for (auto iterIdValue = backPlaneCorrection->getBackPlaneCorrections().begin();
|
524 |
|
|
iterIdValue != backPlaneCorrection->getBackPlaneCorrections().end(); ++iterIdValue) {
|
525 |
|
|
// type of (*iterIdValue) is pair<unsigned int, float>
|
526 |
|
|
id = iterIdValue->first; // key of map is DetId
|
527 |
|
|
value = iterIdValue->second;
|
528 |
|
|
auto idErrPairIter = errors.find(id); // find error for this id - if none, fill 0. in tree
|
529 |
|
|
error = (idErrPairIter != errors.end() ? idErrPairIter->second : 0.f);
|
530 |
|
|
tree->Fill();
|
531 |
|
|
}
|
532 |
|
|
tree->Write();
|
533 |
|
|
delete file; // tree vanishes with the file... (?)
|
534 |
flucke |
1.1.2.1 |
|
535 |
flucke |
1.1.2.2 |
}
|
536 |
flucke |
1.1.2.1 |
|
537 |
flucke |
1.1.2.2 |
//======================================================================
|
538 |
|
|
SiStripBackPlaneCorrection*
|
539 |
|
|
SiStripBackplaneCalibration::createFromTree(const char *fileName, const char *treeName) const
|
540 |
|
|
{
|
541 |
|
|
// Check for file existence on your own to work around
|
542 |
|
|
// https://hypernews.cern.ch/HyperNews/CMS/get/swDevelopment/2715.html:
|
543 |
|
|
TFile* file = 0;
|
544 |
|
|
FILE* testFile = fopen(fileName,"r");
|
545 |
|
|
if (testFile) {
|
546 |
|
|
fclose(testFile);
|
547 |
|
|
file = TFile::Open(fileName, "READ");
|
548 |
|
|
} // else not existing, see error below
|
549 |
|
|
|
550 |
|
|
TTree *tree = 0;
|
551 |
|
|
if (file) file->GetObject(treeName, tree);
|
552 |
|
|
|
553 |
|
|
SiStripBackPlaneCorrection *result = 0;
|
554 |
|
|
if (tree) {
|
555 |
|
|
unsigned int id = 0;
|
556 |
|
|
float value = 0.;
|
557 |
|
|
tree->SetBranchAddress("detId", &id);
|
558 |
|
|
tree->SetBranchAddress("value", &value);
|
559 |
|
|
|
560 |
|
|
result = new SiStripBackPlaneCorrection;
|
561 |
|
|
const Long64_t nEntries = tree->GetEntries();
|
562 |
|
|
for (Long64_t iEntry = 0; iEntry < nEntries; ++iEntry) {
|
563 |
|
|
tree->GetEntry(iEntry);
|
564 |
|
|
result->putBackPlaneCorrection(id, value);
|
565 |
|
|
}
|
566 |
|
|
} else { // Warning only since could be parallel job on no events.
|
567 |
|
|
edm::LogWarning("Alignment") << "@SUB=SiStripBackplaneCalibration::createFromTree"
|
568 |
|
|
<< "Could not get TTree '" << treeName << "' from file '"
|
569 |
|
|
<< fileName << (file ? "'." : "' (file does not exist).");
|
570 |
|
|
}
|
571 |
flucke |
1.1.2.1 |
|
572 |
flucke |
1.1.2.2 |
delete file; // tree will vanish with file
|
573 |
|
|
return result;
|
574 |
|
|
}
|
575 |
flucke |
1.1.2.1 |
|
576 |
|
|
|
577 |
|
|
//======================================================================
|
578 |
|
|
//======================================================================
|
579 |
|
|
// Plugin definition
|
580 |
|
|
|
581 |
|
|
#include "Alignment/CommonAlignmentAlgorithm/interface/IntegratedCalibrationPluginFactory.h"
|
582 |
|
|
|
583 |
|
|
DEFINE_EDM_PLUGIN(IntegratedCalibrationPluginFactory,
|
584 |
|
|
SiStripBackplaneCalibration, "SiStripBackplaneCalibration");
|