ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/CMSSW/Alignment/CommonAlignmentAlgorithm/plugins/SiStripBackplaneCalibration.cc
Revision: 1.1.2.6
Committed: Wed May 15 15:20:35 2013 UTC (11 years, 11 months ago) by jbehr
Content type: text/plain
Branch: branch53X_calibration
Changes since 1.1.2.5: +36 -29 lines
Log Message:
Some code refactoring following GF comments.

File Contents

# User Rev Content
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.6 /// $Revision: 1.1.2.5 $
13     /// $Date: 2013/05/14 08:01:04 $
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 flucke 1.1.2.2 // SiStripLatency::singleReadOutMode() returns
163     // 1: all in peak, 0: all in deco, -1: mixed state
164     // (in principle one could treat even mixed state APV by APV...)
165     if (readoutModeName_ == "peak") {
166     readoutMode_ = kPeakMode;
167     } else if (readoutModeName_ == "deconvolution") {
168     readoutMode_ = kDeconvolutionMode;
169     } else {
170     throw cms::Exception("BadConfig")
171     << "SiStripBackplaneCalibration:\n" << "Unknown mode '"
172     << readoutModeName_ << "', should be 'peak' or 'deconvolution' .\n";
173     }
174    
175 flucke 1.1.2.1 // FIXME: Which granularity, leading to how many parameters?
176 flucke 1.1.2.2 parameters_.resize(2, 0.); // currently two parameters (TIB, TOB), start value 0.
177     paramUncertainties_.resize(2, 0.); // dito for errors
178 flucke 1.1.2.1
179 jbehr 1.1.2.6
180 flucke 1.1.2.1 }
181    
182     //======================================================================
183     SiStripBackplaneCalibration::~SiStripBackplaneCalibration()
184     {
185 flucke 1.1.2.2 // std::cout << "Destroy SiStripBackplaneCalibration named " << this->name() << std::endl;
186     delete siStripBackPlaneCorrInput_;
187 flucke 1.1.2.1 }
188    
189     //======================================================================
190     unsigned int SiStripBackplaneCalibration::numParameters() const
191     {
192     return parameters_.size();
193     }
194    
195     //======================================================================
196     unsigned int
197     SiStripBackplaneCalibration::derivatives(std::vector<ValuesIndexPair> &outDerivInds,
198     const TransientTrackingRecHit &hit,
199     const TrajectoryStateOnSurface &tsos,
200     const edm::EventSetup &setup,
201     const EventInfo &eventInfo) const
202     {
203     // ugly const-cast:
204     // But it is either only first initialisation or throwing an exception...
205 flucke 1.1.2.2 const_cast<SiStripBackplaneCalibration*>(this)->checkBackPlaneCorrectionInput(setup, eventInfo);
206 flucke 1.1.2.1
207     outDerivInds.clear();
208    
209     edm::ESHandle<SiStripLatency> latency;
210     setup.get<SiStripLatencyRcd>().get(latency);
211     const int16_t mode = latency->singleReadOutMode();
212 flucke 1.1.2.3 // std::cout << "SiStripBackplaneCalibration in mode '" << readoutModeName_ << "' finds mode "
213     // << mode << std::endl;
214     //
215     // std::vector<SiStripLatency::Latency> latencies = const_cast<SiStripLatency*>(latency.product())->allUniqueLatencyAndModes();
216     // for (auto i = latencies.begin(); i != latencies.end(); ++i) {
217     // std::cout << static_cast<int>(i->latency) << ", mode "
218     // << static_cast<int>(i->mode) << ", id_apv " << i->detIdAndApv
219     // << std::endl;
220     // }
221    
222 flucke 1.1.2.2 if(mode == readoutMode_) {
223 flucke 1.1.2.1 if (hit.det()) { // otherwise 'constraint hit' or whatever
224 flucke 1.1.2.2
225 jbehr 1.1.2.5 const int index = moduleGroupSelector_.getParameterIndexFromDetId(hit.det()->geographicalId(),
226 jbehr 1.1.2.4 eventInfo.eventId_.run());
227 flucke 1.1.2.1 if (index >= 0) { // otherwise not treated
228     edm::ESHandle<MagneticField> magneticField;
229     setup.get<IdealMagneticFieldRecord>().get(magneticField);
230     const GlobalVector bField(magneticField->inTesla(hit.det()->surface().position()));
231     const LocalVector bFieldLocal(hit.det()->surface().toLocal(bField));
232 flucke 1.1.2.2 //std::cout << "SiStripBackplaneCalibration derivatives " << readoutModeName_ << std::endl;
233 flucke 1.1.2.1 const double dZ = hit.det()->surface().bounds().thickness(); // it's a float only...
234     const double tanPsi = tsos.localParameters().mixedFormatVector()[1]; //float...
235    
236     edm::ESHandle<SiStripLorentzAngle> lorentzAngleHandle;
237 flucke 1.1.2.2 setup.get<SiStripLorentzAngleRcd>().get(readoutModeName_, lorentzAngleHandle);
238 flucke 1.1.2.1 // Yes, mobility (= LA/By) stored in object called LA...
239     const double mobility = lorentzAngleHandle->getLorentzAngle(hit.det()->geographicalId());
240     // shift due to dead back plane has two parts:
241     // 1) Lorentz Angle correction formula gets reduced thickness dz*(1-bp)
242     // 2) 'Direct' effect is shift of effective module position in local z by bp*dz/2
243     // (see GF's presentation in alignment meeting 25.10.2012,
244     // https://indico.cern.ch/conferenceDisplay.py?confId=174266#2012-10-25)
245     const double xDerivative = 0.5 * dZ * (mobility * bFieldLocal.y() - tanPsi);
246     // std::cout << "derivative is " << xDerivative << " for index " << index
247     // << std::endl;
248     const Values derivs(xDerivative, 0.); // yDerivative = 0.
249     outDerivInds.push_back(ValuesIndexPair(derivs, index));
250 flucke 1.1.2.2 }
251 flucke 1.1.2.1 } else {
252 flucke 1.1.2.2 edm::LogWarning("Alignment") << "@SUB=SiStripBackplaneCalibration::derivatives1"
253 flucke 1.1.2.1 << "Hit without GeomDet, skip!";
254     }
255 flucke 1.1.2.2 } else if (mode != kDeconvolutionMode && mode != kPeakMode) {
256     // warn only if unknown/mixed mode
257 flucke 1.1.2.1 edm::LogWarning("Alignment") << "@SUB=SiStripBackplaneCalibration::derivatives2"
258 flucke 1.1.2.2 << "Readout mode is " << mode << ", but looking for "
259     << readoutMode_ << " (" << readoutModeName_ << ").";
260 flucke 1.1.2.1 }
261    
262     return outDerivInds.size();
263     }
264    
265     //======================================================================
266     bool SiStripBackplaneCalibration::setParameter(unsigned int index, double value)
267     {
268     if (index >= parameters_.size()) {
269     return false;
270     } else {
271 jbehr 1.1.2.6 parameters_[index] = value;
272 flucke 1.1.2.1 return true;
273     }
274     }
275    
276     //======================================================================
277     bool SiStripBackplaneCalibration::setParameterError(unsigned int index, double error)
278     {
279     if (index >= paramUncertainties_.size()) {
280     return false;
281     } else {
282 jbehr 1.1.2.6 paramUncertainties_[index] = error;
283 flucke 1.1.2.1 return true;
284     }
285     }
286    
287     //======================================================================
288     double SiStripBackplaneCalibration::getParameter(unsigned int index) const
289     {
290     // if (index >= parameters_.size()) {
291     // return 0.;
292     // } else {
293 jbehr 1.1.2.6 // return parameters_[index];
294 flucke 1.1.2.1 // }
295 jbehr 1.1.2.6 return (index >= parameters_.size() ? 0. : parameters_[index]);
296 flucke 1.1.2.1 }
297    
298     //======================================================================
299     double SiStripBackplaneCalibration::getParameterError(unsigned int index) const
300     {
301     // if (index >= paramUncertainties_.size()) {
302     // return 0.;
303     // } else {
304     // return paramUncertainties_[index];
305     // }
306 jbehr 1.1.2.6 return (index >= paramUncertainties_.size() ? 0. : paramUncertainties_[index]);
307 flucke 1.1.2.1 }
308    
309 jbehr 1.1.2.6 //======================================================================
310     void SiStripBackplaneCalibration::beginOfJob(AlignableTracker *aliTracker,
311     AlignableMuon *aliMuon,
312     AlignableExtras *aliExtras)
313     {
314     moduleGroupSelector_.createModuleGroups(aliTracker,aliMuon,aliExtras);
315    
316     parameters_.resize(moduleGroupSelector_.getNumberOfParameters(), 0.);
317     paramUncertainties_.resize(moduleGroupSelector_.getNumberOfParameters(), 0.);
318    
319     edm::LogInfo("Alignment") << "@SUB=SiStripBackplaneCalibration" << "Created with name "
320     << this->name() << " for readout mode '" << readoutModeName_
321     << "',\n" << this->numParameters() << " parameters to be determined."
322     << "\nsaveToDB = " << saveToDB_
323     << "\n outFileName = " << outFileName_
324     << "\n N(merge files) = " << mergeFileNames_.size()
325     << "\n number of IOVs = " << moduleGroupSelector_.numIovs();
326     if (mergeFileNames_.size()) {
327     edm::LogInfo("Alignment") << "@SUB=SiStripBackplaneCalibration"
328     << "First file to merge: " << mergeFileNames_[0];
329     }
330    
331     }
332 flucke 1.1.2.1
333    
334     //======================================================================
335     void SiStripBackplaneCalibration::endOfJob()
336     {
337     // loginfo output
338     std::ostringstream out;
339 flucke 1.1.2.2 out << "Parameter results for readout mode '" << readoutModeName_ << "'\n";
340 flucke 1.1.2.1 for (unsigned int iPar = 0; iPar < parameters_.size(); ++iPar) {
341     out << iPar << ": " << parameters_[iPar] << " +- " << paramUncertainties_[iPar] << "\n";
342     }
343     edm::LogInfo("Alignment") << "@SUB=SiStripBackplaneCalibration::endOfJob" << out.str();
344    
345 flucke 1.1.2.2 std::map<unsigned int, float> errors; // Array of errors for each detId
346 flucke 1.1.2.1
347 flucke 1.1.2.2 // now write 'input' tree
348     const SiStripBackPlaneCorrection *input = this->getBackPlaneCorrectionInput(); // never NULL
349     const std::string treeName(this->name() + '_' + readoutModeName_ + '_');
350     this->writeTree(input, errors, (treeName + "input").c_str()); // empty errors for input...
351    
352     if (input->getBackPlaneCorrections().empty()) {
353     edm::LogError("Alignment") << "@SUB=SiStripBackplaneCalibration::endOfJob"
354     << "Input back plane correction map is empty ('"
355     << readoutModeName_ << "' mode), skip writing output!";
356     return;
357 flucke 1.1.2.1 }
358    
359 flucke 1.1.2.3 const unsigned int nonZeroParamsOrErrors = // Any determined value?
360     count_if (parameters_.begin(), parameters_.end(), std::bind2nd(std::not_equal_to<double>(),0.))
361     + count_if(paramUncertainties_.begin(), paramUncertainties_.end(),
362     std::bind2nd(std::not_equal_to<double>(), 0.));
363    
364 jbehr 1.1.2.5 for (unsigned int iIOV = 0; iIOV < moduleGroupSelector_.numIovs(); ++iIOV) {
365     cond::Time_t firstRunOfIOV = moduleGroupSelector_.firstRunOfIOV(iIOV);
366 flucke 1.1.2.2 SiStripBackPlaneCorrection *output = new SiStripBackPlaneCorrection;
367     // Loop on map of values from input and add (possible) parameter results
368     for (auto iterIdValue = input->getBackPlaneCorrections().begin();
369     iterIdValue != input->getBackPlaneCorrections().end(); ++iterIdValue) {
370     // type of (*iterIdValue) is pair<unsigned int, float>
371     const unsigned int detId = iterIdValue->first; // key of map is DetId
372     const float value = iterIdValue->second + this->getParameterForDetId(detId, firstRunOfIOV);
373     output->putBackPlaneCorrection(detId, value); // put result in output
374 jbehr 1.1.2.5 int parameterIndex = moduleGroupSelector_.getParameterIndexFromDetId(detId, firstRunOfIOV);
375 flucke 1.1.2.2 errors[detId] = this->getParameterError(parameterIndex);
376     }
377 flucke 1.1.2.1
378 flucke 1.1.2.3 if (saveToDB_ || nonZeroParamsOrErrors != 0) { // Skip writing mille jobs...
379     this->writeTree(output, errors, (treeName + Form("result_%lld", firstRunOfIOV)).c_str());
380     }
381 flucke 1.1.2.1
382 flucke 1.1.2.2 if (saveToDB_) { // If requested, write out to DB
383     edm::Service<cond::service::PoolDBOutputService> dbService;
384     if (dbService.isAvailable()) {
385     dbService->writeOne(output, firstRunOfIOV, recordNameDBwrite_.c_str());
386     // no 'delete output;': writeOne(..) took over ownership
387     } else {
388     delete output;
389     edm::LogError("BadConfig") << "@SUB=SiStripBackplaneCalibration::endOfJob"
390     << "No PoolDBOutputService available, but saveToDB true!";
391     }
392     } else {
393     delete output;
394     }
395     } // end loop on IOVs
396 flucke 1.1.2.1 }
397    
398     //======================================================================
399 flucke 1.1.2.2 bool SiStripBackplaneCalibration::checkBackPlaneCorrectionInput(const edm::EventSetup &setup,
400     const EventInfo &eventInfo)
401 flucke 1.1.2.1 {
402 flucke 1.1.2.2 edm::ESHandle<SiStripBackPlaneCorrection> backPlaneCorrHandle;
403     if (!siStripBackPlaneCorrInput_) {
404     setup.get<SiStripBackPlaneCorrectionRcd>().get(readoutModeName_, backPlaneCorrHandle);
405     siStripBackPlaneCorrInput_ = new SiStripBackPlaneCorrection(*backPlaneCorrHandle);
406     // FIXME: Should we call 'watchBackPlaneCorrRcd_.check(setup)' as well?
407     // Otherwise could be that next check has to check via following 'else', though
408     // no new IOV has started... (to be checked)
409     } else {
410     if (watchBackPlaneCorrRcd_.check(setup)) { // new IOV of input - but how to check peak vs deco?
411     setup.get<SiStripBackPlaneCorrectionRcd>().get(readoutModeName_, backPlaneCorrHandle);
412     if (backPlaneCorrHandle->getBackPlaneCorrections() // but only bad if non-identical values
413     != siStripBackPlaneCorrInput_->getBackPlaneCorrections()) { // (comparing maps)
414     // Maps are containers sorted by key, but comparison problems may arise from
415     // 'floating point comparison' problems (FIXME?)
416     throw cms::Exception("BadInput")
417     << "SiStripBackplaneCalibration::checkBackPlaneCorrectionInput:\n"
418     << "Content of SiStripBackPlaneCorrection changed at run " << eventInfo.eventId_.run()
419     << ", but algorithm expects constant input!\n";
420     return false; // not reached...
421     }
422 flucke 1.1.2.1 }
423     }
424    
425     return true;
426     }
427    
428     //======================================================================
429 flucke 1.1.2.2 const SiStripBackPlaneCorrection* SiStripBackplaneCalibration::getBackPlaneCorrectionInput()
430 flucke 1.1.2.1 {
431 flucke 1.1.2.2 // For parallel processing in Millepede II, create SiStripBackPlaneCorrection
432     // from info stored in files of parallel jobs and check that they are identical.
433     // If this job has run on events, still check that back plane corrections are
434     // identical to the ones from mergeFileNames_.
435     const std::string treeName(((this->name() + '_') += readoutModeName_) += "_input");
436     for (auto iFile = mergeFileNames_.begin(); iFile != mergeFileNames_.end(); ++iFile) {
437     SiStripBackPlaneCorrection* bpCorr = this->createFromTree(iFile->c_str(), treeName.c_str());
438     // siStripBackPlaneCorrInput_ could be non-null from previous file of this loop
439     // or from checkBackPlaneCorrectionInput(..) when running on data in this job as well
440     if (!siStripBackPlaneCorrInput_ || siStripBackPlaneCorrInput_->getBackPlaneCorrections().empty()) {
441     delete siStripBackPlaneCorrInput_; // NULL or empty
442     siStripBackPlaneCorrInput_ = bpCorr;
443     } else {
444     // about comparison of maps see comments in checkBackPlaneCorrectionInput
445     if (bpCorr && !bpCorr->getBackPlaneCorrections().empty() && // single job might not have got events
446     bpCorr->getBackPlaneCorrections() != siStripBackPlaneCorrInput_->getBackPlaneCorrections()) {
447     // Throw exception instead of error?
448     edm::LogError("NoInput") << "@SUB=SiStripBackplaneCalibration::getBackPlaneCorrectionInput"
449     << "Different input values from tree " << treeName
450     << " in file " << *iFile << ".";
451 flucke 1.1.2.1
452 flucke 1.1.2.2 }
453     delete bpCorr;
454     }
455     }
456 flucke 1.1.2.1
457 flucke 1.1.2.2 if (!siStripBackPlaneCorrInput_) { // no files nor ran on events
458     siStripBackPlaneCorrInput_ = new SiStripBackPlaneCorrection;
459     edm::LogError("NoInput") << "@SUB=SiStripBackplaneCalibration::getBackPlaneCorrectionInput"
460     << "No input, create an empty one ('" << readoutModeName_ << "' mode)!";
461     } else if (siStripBackPlaneCorrInput_->getBackPlaneCorrections().empty()) {
462     edm::LogError("NoInput") << "@SUB=SiStripBackplaneCalibration::getBackPlaneCorrectionInput"
463     << "Empty result ('" << readoutModeName_ << "' mode)!";
464     }
465    
466     return siStripBackPlaneCorrInput_;
467 flucke 1.1.2.1 }
468    
469     //======================================================================
470     double SiStripBackplaneCalibration::getParameterForDetId(unsigned int detId,
471     edm::RunNumber_t run) const
472     {
473 jbehr 1.1.2.5 const int index = moduleGroupSelector_.getParameterIndexFromDetId(detId, run);
474 flucke 1.1.2.1
475 jbehr 1.1.2.6 return (index < 0 ? 0. : parameters_[index]);
476 flucke 1.1.2.1 }
477    
478 flucke 1.1.2.2 //======================================================================
479     void SiStripBackplaneCalibration::writeTree(const SiStripBackPlaneCorrection *backPlaneCorrection,
480     const std::map<unsigned int, float> &errors,
481     const char *treeName) const
482     {
483     if (!backPlaneCorrection) return;
484    
485     TFile* file = TFile::Open(outFileName_.c_str(), "UPDATE");
486     if (!file) {
487     edm::LogError("BadConfig") << "@SUB=SiStripBackplaneCalibration::writeTree"
488     << "Could not open file '" << outFileName_ << "'.";
489     return;
490     }
491 flucke 1.1.2.1
492 flucke 1.1.2.2 TTree *tree = new TTree(treeName, treeName);
493     unsigned int id = 0;
494     float value = 0.;
495     float error = 0.;
496     tree->Branch("detId", &id, "detId/i");
497     tree->Branch("value", &value, "value/F");
498     tree->Branch("error", &error, "error/F");
499    
500     for (auto iterIdValue = backPlaneCorrection->getBackPlaneCorrections().begin();
501     iterIdValue != backPlaneCorrection->getBackPlaneCorrections().end(); ++iterIdValue) {
502     // type of (*iterIdValue) is pair<unsigned int, float>
503     id = iterIdValue->first; // key of map is DetId
504     value = iterIdValue->second;
505     auto idErrPairIter = errors.find(id); // find error for this id - if none, fill 0. in tree
506     error = (idErrPairIter != errors.end() ? idErrPairIter->second : 0.f);
507     tree->Fill();
508     }
509     tree->Write();
510     delete file; // tree vanishes with the file... (?)
511 flucke 1.1.2.1
512 flucke 1.1.2.2 }
513 flucke 1.1.2.1
514 flucke 1.1.2.2 //======================================================================
515     SiStripBackPlaneCorrection*
516     SiStripBackplaneCalibration::createFromTree(const char *fileName, const char *treeName) const
517     {
518     // Check for file existence on your own to work around
519     // https://hypernews.cern.ch/HyperNews/CMS/get/swDevelopment/2715.html:
520     TFile* file = 0;
521     FILE* testFile = fopen(fileName,"r");
522     if (testFile) {
523     fclose(testFile);
524     file = TFile::Open(fileName, "READ");
525     } // else not existing, see error below
526    
527     TTree *tree = 0;
528     if (file) file->GetObject(treeName, tree);
529    
530     SiStripBackPlaneCorrection *result = 0;
531     if (tree) {
532     unsigned int id = 0;
533     float value = 0.;
534     tree->SetBranchAddress("detId", &id);
535     tree->SetBranchAddress("value", &value);
536    
537     result = new SiStripBackPlaneCorrection;
538     const Long64_t nEntries = tree->GetEntries();
539     for (Long64_t iEntry = 0; iEntry < nEntries; ++iEntry) {
540     tree->GetEntry(iEntry);
541     result->putBackPlaneCorrection(id, value);
542     }
543     } else { // Warning only since could be parallel job on no events.
544     edm::LogWarning("Alignment") << "@SUB=SiStripBackplaneCalibration::createFromTree"
545     << "Could not get TTree '" << treeName << "' from file '"
546     << fileName << (file ? "'." : "' (file does not exist).");
547     }
548 flucke 1.1.2.1
549 flucke 1.1.2.2 delete file; // tree will vanish with file
550     return result;
551     }
552 flucke 1.1.2.1
553    
554     //======================================================================
555     //======================================================================
556     // Plugin definition
557    
558     #include "Alignment/CommonAlignmentAlgorithm/interface/IntegratedCalibrationPluginFactory.h"
559    
560     DEFINE_EDM_PLUGIN(IntegratedCalibrationPluginFactory,
561     SiStripBackplaneCalibration, "SiStripBackplaneCalibration");