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