CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
SiStripGainFromCalibTree.cc
Go to the documentation of this file.
1 // Original Author: Loic QUERTENMONT
2 // Created: Mon Nov 16 08:55:18 CET 2009
3 
4 #include <memory>
5 #include <iostream>
6 
7 #include "FWCore/Framework/interface/Frameworkfwd.h"
8 #include "FWCore/Framework/interface/EDAnalyzer.h"
9 #include "FWCore/Framework/interface/Event.h"
10 #include "FWCore/Framework/interface/MakerMacros.h"
11 #include "FWCore/Framework/interface/ESHandle.h"
12 #include "FWCore/ParameterSet/interface/ParameterSet.h"
13 #include "FWCore/MessageLogger/interface/MessageLogger.h"
14 #include "FWCore/Utilities/interface/Exception.h"
15 #include "FWCore/Utilities/interface/EDGetToken.h"
16 
17 #include "CondFormats/RunInfo/interface/RunInfo.h"
18 
19 #include "Geometry/CommonDetUnit/interface/GeomDetUnit.h"
20 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
21 #include "Geometry/CommonTopologies/interface/StripTopology.h"
22 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
23 #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
24 
25 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
26 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
27 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
28 #include "Geometry/TrackerGeometryBuilder/interface/PixelGeomDetUnit.h"
29 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
30 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
31 
32 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
33 
34 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
35 #include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
36 
37 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
38 #include "DataFormats/TrackReco/interface/Track.h"
39 #include "DataFormats/TrackReco/interface/TrackFwd.h"
40 #include "DataFormats/TrackerRecHit2D/interface/SiPixelRecHit.h"
41 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
42 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
43 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
44 #include "DataFormats/SiStripDetId/interface/SiStripSubStructure.h"
45 #include "DataFormats/DetId/interface/DetId.h"
46 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
47 #include "DataFormats/TrackReco/interface/DeDxHit.h"
48 #include "DataFormats/TrackReco/interface/TrackDeDxHits.h"
49 
50 #include "CommonTools/ConditionDBWriter/interface/ConditionDBWriter.h"
51 #include "CondFormats/SiStripObjects/interface/SiStripApvGain.h"
52 
53 #include "TrackingTools/PatternTools/interface/Trajectory.h"
54 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
55 
56 #include "DQMServices/Core/interface/DQMStore.h"
57 #include "DQMServices/Core/interface/MonitorElement.h"
58 #include "FWCore/ServiceRegistry/interface/Service.h"
59 #include "CommonTools/UtilAlgos/interface/TFileService.h"
60 
61 #include "CalibFormats/SiStripObjects/interface/SiStripGain.h"
62 #include "CalibTracker/Records/interface/SiStripGainRcd.h"
63 
64 #include "CalibFormats/SiStripObjects/interface/SiStripQuality.h"
65 #include "CalibTracker/Records/interface/SiStripQualityRcd.h"
66 
67 
68 #include "TFile.h"
69 #include "TObjString.h"
70 #include "TString.h"
71 #include "TH1F.h"
72 #include "TH2S.h"
73 #include "TProfile.h"
74 #include "TF1.h"
75 #include "TROOT.h"
76 #include "TTree.h"
77 #include "TChain.h"
78 
79 #include <ext/hash_map>
80 
81 
82 
83 using namespace edm;
84 using namespace reco;
85 using namespace std;
86 using namespace __gnu_cxx;
87 using __gnu_cxx::hash_map;
88 using __gnu_cxx::hash;
89 
90 struct stAPVGain{
91  unsigned int Index;
92  int Bin;
93  unsigned int DetId;
94  unsigned int APVId;
95  unsigned int SubDet;
96  float x;
97  float y;
98  float z;
99  float Eta;
100  float R;
101  float Phi;
102  float Thickness;
103  double FitMPV;
104  double FitMPVErr;
105  double FitWidth;
106  double FitWidthErr;
107  double FitChi2;
108  double FitNorm;
109  double Gain;
110  double CalibGain;
111  double PreviousGain;
112  double PreviousGainTick;
113  double NEntries;
114  TH1F* HCharge;
115  TH1F* HChargeN;
116  bool isMasked;
117 };
118 
119 class SiStripGainFromCalibTree : public ConditionDBWriter<SiStripApvGain> {
120 public:
121  explicit SiStripGainFromCalibTree(const edm::ParameterSet&);
122  ~SiStripGainFromCalibTree();
123 
124 
125 private:
126 
127 
128  virtual void algoBeginRun(const edm::Run& run, const edm::EventSetup& iSetup) override;
129  virtual void algoEndRun (const edm::Run& run, const edm::EventSetup& iSetup) override;
130  virtual void algoBeginJob (const edm::EventSetup& iSetup) override;
131  virtual void algoEndJob () override;
132  virtual void algoAnalyze (const edm::Event &, const edm::EventSetup &) override;
133 
134  int statCollectionFromMode(const char* tag) const;
135  void bookDQMHistos(const char* dqm_dir, const char* tag);
136 
137  bool isBFieldConsistentWithMode( const edm::EventSetup& iSetup) const;
138  void swapBFieldMode(void);
139 
140  void merge(TH2* A, TH2* B); //needed to add histograms with different number of bins
141  void algoAnalyzeTheTree();
142  void algoComputeMPVandGain();
143  void processEvent(); //what really does the job
144 
145  void getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange=50, double HighRange=5400);
146  bool IsGoodLandauFit(double* FitResults);
147  void storeOnTree(TFileService* tfs);
148  void MakeCalibrationMap();
149  bool produceTagFilter();
150 
151  template<typename T>
152  inline edm::Handle<T> connect(const T* &ptr, edm::EDGetTokenT<T> token, const edm::Event &evt) {
153  edm::Handle<T> handle;
154  evt.getByToken(token, handle);
155  ptr = handle.product();
156  return handle; //return handle to keep alive pointer (safety first)
157  }
158 
159  SiStripApvGain* getNewObject() override;
160 
161  TFileService *tfs;
162  DQMStore* dbe;
163  double MagFieldCurrentTh;
164  double MinNrEntries;
165  double MaxMPVError;
166  double MaxChi2OverNDF;
167  double MinTrackMomentum;
168  double MaxTrackMomentum;
169  double MinTrackEta;
170  double MaxTrackEta;
171  unsigned int MaxNrStrips;
172  unsigned int MinTrackHits;
173  double MaxTrackChiOverNdf;
174  int MaxTrackingIteration;
175  bool AllowSaturation;
176  bool FirstSetOfConstants;
177  bool Validation;
178  bool OldGainRemoving;
179  int CalibrationLevel;
180 
181  bool saveSummary;
182  bool useCalibration;
183  bool m_harvestingMode;
184  bool m_splitDQMstat;
185  string m_calibrationMode;
186  string m_calibrationPath;
187  string m_DQMdir;
188 
189  double tagCondition_NClusters;
190  double tagCondition_GoodFrac;
191 
192  std::string AlgoMode;
193  std::string OutputGains;
194  vector<string> VInputFiles;
195 
196  enum statistic_type {None=-1, StdBunch, StdBunch0T, FaABunch, FaABunch0T, IsoBunch, IsoBunch0T, Harvest};
197 
198  std::vector<string> dqm_tag_;
199 
200  std::vector<MonitorElement*> Charge_Vs_Index;
201  //std::vector<MonitorElement*> Charge_Vs_Index_Absolute;
202  std::vector<MonitorElement*> Charge_Vs_PathlengthTIB;
203  std::vector<MonitorElement*> Charge_Vs_PathlengthTOB;
204  std::vector<MonitorElement*> Charge_Vs_PathlengthTIDP;
205  std::vector<MonitorElement*> Charge_Vs_PathlengthTIDM;
206  std::vector<MonitorElement*> Charge_Vs_PathlengthTECP1;
207  std::vector<MonitorElement*> Charge_Vs_PathlengthTECP2;
208  std::vector<MonitorElement*> Charge_Vs_PathlengthTECM1;
209  std::vector<MonitorElement*> Charge_Vs_PathlengthTECM2;
210 
211  unsigned int NEvent;
212  unsigned int NTrack;
213  unsigned int NClusterStrip;
214  unsigned int NClusterPixel;
215  int NStripAPVs;
216  int NPixelDets;
217  unsigned int SRun;
218  unsigned int ERun;
219  unsigned int GOOD;
220  unsigned int BAD;
221  unsigned int MASKED;
222 
223  //Data members for processing
224 
225  //Event data
226  unsigned int eventnumber =0;
227  unsigned int runnumber =0;
228  const std::vector<bool>* TrigTech =0; edm::EDGetTokenT<std::vector<bool> > TrigTech_token_;
229 
230  // Track data
231  const std::vector<double>* trackchi2ndof =0; edm::EDGetTokenT<std::vector<double> > trackchi2ndof_token_;
232  const std::vector<float>* trackp =0; edm::EDGetTokenT<std::vector<float> > trackp_token_;
233  const std::vector<float>* trackpt =0; edm::EDGetTokenT<std::vector<float> > trackpt_token_;
234  const std::vector<double>* tracketa =0; edm::EDGetTokenT<std::vector<double> > tracketa_token_;
235  const std::vector<double>* trackphi =0; edm::EDGetTokenT<std::vector<double> > trackphi_token_;
236  const std::vector<unsigned int>* trackhitsvalid =0; edm::EDGetTokenT<std::vector<unsigned int> > trackhitsvalid_token_;
237  const std::vector<int>* trackalgo =0; edm::EDGetTokenT<std::vector<int> > trackalgo_token_;
238 
239  // CalibTree data
240  const std::vector<int>* trackindex =0; edm::EDGetTokenT<std::vector<int> > trackindex_token_;
241  const std::vector<unsigned int>* rawid =0; edm::EDGetTokenT<std::vector<unsigned int> > rawid_token_;
242  const std::vector<double>* localdirx =0; edm::EDGetTokenT<std::vector<double> > localdirx_token_;
243  const std::vector<double>* localdiry =0; edm::EDGetTokenT<std::vector<double> > localdiry_token_;
244  const std::vector<double>* localdirz =0; edm::EDGetTokenT<std::vector<double> > localdirz_token_;
245  const std::vector<unsigned short>* firststrip =0; edm::EDGetTokenT<std::vector<unsigned short> > firststrip_token_;
246  const std::vector<unsigned short>* nstrips =0; edm::EDGetTokenT<std::vector<unsigned short> > nstrips_token_;
247  const std::vector<bool>* saturation =0; edm::EDGetTokenT<std::vector<bool> > saturation_token_;
248  const std::vector<bool>* overlapping =0; edm::EDGetTokenT<std::vector<bool> > overlapping_token_;
249  const std::vector<bool>* farfromedge =0; edm::EDGetTokenT<std::vector<bool> > farfromedge_token_;
250  const std::vector<unsigned int>* charge =0; edm::EDGetTokenT<std::vector<unsigned int> > charge_token_;
251  const std::vector<double>* path =0; edm::EDGetTokenT<std::vector<double> > path_token_;
252  const std::vector<double>* chargeoverpath =0; edm::EDGetTokenT<std::vector<double> > chargeoverpath_token_;
253  const std::vector<unsigned char>* amplitude =0; edm::EDGetTokenT<std::vector<unsigned char> > amplitude_token_;
254  const std::vector<double>* gainused =0; edm::EDGetTokenT<std::vector<double> > gainused_token_;
255 
256  string EventPrefix_; //("");
257  string EventSuffix_; //("");
258  string TrackPrefix_; //("track");
259  string TrackSuffix_; //("");
260  string CalibPrefix_; //("GainCalibration");
261  string CalibSuffix_; //("");
262 
263 private :
264  class isEqual{
265  public:
266  template <class T> bool operator () (const T& PseudoDetId1, const T& PseudoDetId2) { return PseudoDetId1==PseudoDetId2; }
267  };
268 
269  std::vector<stAPVGain*> APVsCollOrdered;
270  __gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual > APVsColl;
271 };
272 
273 inline int
274 SiStripGainFromCalibTree::statCollectionFromMode(const char* tag) const
275 {
276  std::vector<string>::const_iterator it=dqm_tag_.begin();
277  while(it!=dqm_tag_.end()) {
278  if(*it==std::string(tag)) return it-dqm_tag_.begin();
279  it++;
280  }
281 
282  if (std::string(tag)=="") return 0; // return StdBunch calibration mode for backward compatibility
283 
284  return None;
285 }
286 
287 void SiStripGainFromCalibTree::merge(TH2* A, TH2* B){
288  if(A->GetNbinsX() == B->GetNbinsX()){
289  A->Add(B);
290  }else{
291  for(int x=0;x<=B->GetNbinsX()+1; x++){
292  for(int y=0;y<=B->GetNbinsY()+1; y++){
293  A->SetBinContent(x,y,A->GetBinContent(x,y)+B->GetBinContent(x,y));
294  }}
295  }
296 }
297 
298 
299 
300 SiStripGainFromCalibTree::SiStripGainFromCalibTree(const edm::ParameterSet& iConfig) : ConditionDBWriter<SiStripApvGain>(iConfig)
301 {
302  OutputGains = iConfig.getParameter<std::string>("OutputGains");
303 
304  AlgoMode = iConfig.getUntrackedParameter<std::string>("AlgoMode", "CalibTree");
305  MagFieldCurrentTh = iConfig.getUntrackedParameter<double> ("MagFieldCurrentTh" , 2000.);
306  MinNrEntries = iConfig.getUntrackedParameter<double> ("minNrEntries" , 20);
307  MaxMPVError = iConfig.getUntrackedParameter<double> ("maxMPVError" , 500.0);
308  MaxChi2OverNDF = iConfig.getUntrackedParameter<double> ("maxChi2OverNDF" , 5.0);
309  MinTrackMomentum = iConfig.getUntrackedParameter<double> ("minTrackMomentum" , 3.0);
310  MaxTrackMomentum = iConfig.getUntrackedParameter<double> ("maxTrackMomentum" , 99999.0);
311  MinTrackEta = iConfig.getUntrackedParameter<double> ("minTrackEta" , -5.0);
312  MaxTrackEta = iConfig.getUntrackedParameter<double> ("maxTrackEta" , 5.0);
313  MaxNrStrips = iConfig.getUntrackedParameter<unsigned>("maxNrStrips" , 2);
314  MinTrackHits = iConfig.getUntrackedParameter<unsigned>("MinTrackHits" , 8);
315  MaxTrackChiOverNdf = iConfig.getUntrackedParameter<double> ("MaxTrackChiOverNdf" , 3);
316  MaxTrackingIteration = iConfig.getUntrackedParameter<int> ("MaxTrackingIteration" , 7);
317  AllowSaturation = iConfig.getUntrackedParameter<bool> ("AllowSaturation" , false);
318  FirstSetOfConstants = iConfig.getUntrackedParameter<bool> ("FirstSetOfConstants" , true);
319  Validation = iConfig.getUntrackedParameter<bool> ("Validation" , false);
320  OldGainRemoving = iConfig.getUntrackedParameter<bool> ("OldGainRemoving" , false);
321 
322  CalibrationLevel = iConfig.getUntrackedParameter<int> ("CalibrationLevel" , 0);
323  VInputFiles = iConfig.getUntrackedParameter<vector<string> > ("InputFiles");
324 
325 
326  useCalibration = iConfig.getUntrackedParameter<bool> ("UseCalibration" , false);
327  m_harvestingMode = iConfig.getUntrackedParameter<bool> ("harvestingMode" , false);
328  m_splitDQMstat = iConfig.getUntrackedParameter<bool> ("splitDQMstat" , false);
329  m_calibrationMode = iConfig.getUntrackedParameter<string> ("calibrationMode" , "StdBunch");
330  m_calibrationPath = iConfig.getUntrackedParameter<string> ("calibrationPath");
331  m_DQMdir = iConfig.getUntrackedParameter<string> ("DQMdir" , "AlCaReco/SiStripGains");
332 
333  tagCondition_NClusters = iConfig.getUntrackedParameter<double> ("NClustersForTagProd" , 2E8);
334  tagCondition_GoodFrac = iConfig.getUntrackedParameter<double> ("GoodFracForTagProd" , 0.95);
335 
336  saveSummary = iConfig.getUntrackedParameter<bool> ("saveSummary" , false);
337 
338 
339  // Gather DQM Service
340  dbe = edm::Service<DQMStore>().operator->();
341  dbe->setVerbose(10);
342 
343  //Set the monitoring element tag and store
344  dqm_tag_.reserve(7);
345  dqm_tag_.clear();
346  dqm_tag_.push_back( "StdBunch" ); // statistic collection from Standard Collision Bunch @ 3.8 T
347  dqm_tag_.push_back( "StdBunch0T" ); // statistic collection from Standard Collision Bunch @ 0 T
348  dqm_tag_.push_back( "AagBunch" ); // statistic collection from First Collision After Abort Gap @ 3.8 T
349  dqm_tag_.push_back( "AagBunch0T" ); // statistic collection from First Collision After Abort Gap @ 0 T
350  dqm_tag_.push_back( "IsoMuon" ); // statistic collection from Isolated Muon @ 3.8 T
351  dqm_tag_.push_back( "IsoMuon0T" ); // statistic collection from Isolated Muon @ 0 T
352  dqm_tag_.push_back( "Harvest" ); // statistic collection: Harvest
353 
354  Charge_Vs_Index.insert( Charge_Vs_Index.begin(), dqm_tag_.size(), 0);
355  //Charge_Vs_Index_Absolute.insert( Charge_Vs_Index_Absolute.begin(), dqm_tag_.size(), 0);
356  Charge_Vs_PathlengthTIB.insert( Charge_Vs_PathlengthTIB.begin(), dqm_tag_.size(), 0);
357  Charge_Vs_PathlengthTOB.insert( Charge_Vs_PathlengthTOB.begin(), dqm_tag_.size(), 0);
358  Charge_Vs_PathlengthTIDP.insert( Charge_Vs_PathlengthTIDP.begin(), dqm_tag_.size(), 0);
359  Charge_Vs_PathlengthTIDM.insert( Charge_Vs_PathlengthTIDM.begin(), dqm_tag_.size(), 0);
360  Charge_Vs_PathlengthTECP1.insert( Charge_Vs_PathlengthTECP1.begin(), dqm_tag_.size(), 0);
361  Charge_Vs_PathlengthTECP2.insert( Charge_Vs_PathlengthTECP2.begin(), dqm_tag_.size(), 0);
362  Charge_Vs_PathlengthTECM1.insert( Charge_Vs_PathlengthTECM1.begin(), dqm_tag_.size(), 0);
363  Charge_Vs_PathlengthTECM2.insert( Charge_Vs_PathlengthTECM2.begin(), dqm_tag_.size(), 0);
364 
365 
366 
367 
368  // configure token for gathering the ntuple variables
369  edm::ParameterSet swhallowgain_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("gain");
370 
371  string label = swhallowgain_pset.getUntrackedParameter<string>("label");
372  CalibPrefix_ = swhallowgain_pset.getUntrackedParameter<string>("prefix");
373  CalibSuffix_ = swhallowgain_pset.getUntrackedParameter<string>("suffix");
374 
375  trackindex_token_ = consumes<std::vector<int> >(edm::InputTag(label, CalibPrefix_ + "trackindex" + CalibSuffix_));
376  rawid_token_ = consumes<std::vector<unsigned int> >(edm::InputTag(label, CalibPrefix_ + "rawid" + CalibSuffix_));
377  localdirx_token_ = consumes<std::vector<double> >(edm::InputTag(label, CalibPrefix_ + "localdirx" + CalibSuffix_));
378  localdiry_token_ = consumes<std::vector<double> >(edm::InputTag(label, CalibPrefix_ + "localdiry" + CalibSuffix_));
379  localdirz_token_ = consumes<std::vector<double> >(edm::InputTag(label, CalibPrefix_ + "localdirz" + CalibSuffix_));
380  firststrip_token_ = consumes<std::vector<unsigned short> >(edm::InputTag(label, CalibPrefix_ + "firststrip" + CalibSuffix_));
381  nstrips_token_ = consumes<std::vector<unsigned short> >(edm::InputTag(label, CalibPrefix_ + "nstrips" + CalibSuffix_));
382  saturation_token_ = consumes<std::vector<bool> >(edm::InputTag(label, CalibPrefix_ + "saturation" + CalibSuffix_));
383  overlapping_token_ = consumes<std::vector<bool> >(edm::InputTag(label, CalibPrefix_ + "overlapping" + CalibSuffix_));
384  farfromedge_token_ = consumes<std::vector<bool> >(edm::InputTag(label, CalibPrefix_ + "farfromedge" + CalibSuffix_));
385  charge_token_ = consumes<std::vector<unsigned int> >(edm::InputTag(label, CalibPrefix_ + "charge" + CalibSuffix_));
386  path_token_ = consumes<std::vector<double> >(edm::InputTag(label, CalibPrefix_ + "path" + CalibSuffix_));
387  chargeoverpath_token_ = consumes<std::vector<double> >(edm::InputTag(label, CalibPrefix_ + "chargeoverpath"+ CalibSuffix_));
388  amplitude_token_ = consumes<std::vector<unsigned char> >(edm::InputTag(label, CalibPrefix_ + "amplitude" + CalibSuffix_));
389  gainused_token_ = consumes<std::vector<double> >(edm::InputTag(label, CalibPrefix_ + "gainused" + CalibSuffix_));
390 
391  edm::ParameterSet evtinfo_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("evtinfo");
392  label = evtinfo_pset.getUntrackedParameter<string>("label");
393  EventPrefix_ = evtinfo_pset.getUntrackedParameter<string>("prefix");
394  EventSuffix_ = evtinfo_pset.getUntrackedParameter<string>("suffix");
395  TrigTech_token_ = consumes<std::vector<bool> >(edm::InputTag(label, EventPrefix_ + "TrigTech" + EventSuffix_));
396 
397  edm::ParameterSet track_pset = iConfig.getUntrackedParameter<edm::ParameterSet>("tracks");
398  label = track_pset.getUntrackedParameter<string>("label");
399  TrackPrefix_ = track_pset.getUntrackedParameter<string>("prefix");
400  TrackSuffix_ = track_pset.getUntrackedParameter<string>("suffix");
401 
402  trackchi2ndof_token_ = consumes<std::vector<double> >(edm::InputTag(label, TrackPrefix_ + "chi2ndof" + TrackSuffix_));
403  trackp_token_ = consumes<std::vector<float> >(edm::InputTag(label, TrackPrefix_ + "momentum" + TrackSuffix_));
404  trackpt_token_ = consumes<std::vector<float> >(edm::InputTag(label, TrackPrefix_ + "pt" + TrackSuffix_));
405  tracketa_token_ = consumes<std::vector<double> >(edm::InputTag(label, TrackPrefix_ + "eta" + TrackSuffix_));
406  trackphi_token_ = consumes<std::vector<double> >(edm::InputTag(label, TrackPrefix_ + "phi" + TrackSuffix_));
407  trackhitsvalid_token_ = consumes<std::vector<unsigned int> >(edm::InputTag(label, TrackPrefix_ + "hitsvalid" + TrackSuffix_));
408  trackalgo_token_ = consumes<std::vector<int> >(edm::InputTag(label, TrackPrefix_ + "algo" + TrackSuffix_));
409 }
410 
411 void SiStripGainFromCalibTree::bookDQMHistos(const char* dqm_dir, const char* tag)
412 {
413  static std::string booked_dir = "";
414 
415  edm::LogInfo("SiStripGainFromCalibTree") << "Setting " << dqm_dir << "in DQM and booking histograms for tag "
416  << tag << std::endl;
417 
418  if ( strcmp(booked_dir.c_str(),dqm_dir)!=0 ) {
419  booked_dir = dqm_dir;
420  dbe->setCurrentFolder(dqm_dir);
421  }
422 
423  std::string stag(tag);
424  if(stag.size()!=0 && stag[0]!='_') stag.insert(0,1,'_');
425 
426  std::string cvi = std::string("Charge_Vs_Index") + stag;
427  //std::string cviA = std::string("Charge_Vs_Index_Absolute") + stag;
428  std::string cvpTIB = std::string("Charge_Vs_PathlengthTIB") + stag;
429  std::string cvpTOB = std::string("Charge_Vs_PathlengthTOB") + stag;
430  std::string cvpTIDP = std::string("Charge_Vs_PathlengthTIDP") + stag;
431  std::string cvpTIDM = std::string("Charge_Vs_PathlengthTIDM") + stag;
432  std::string cvpTECP1 = std::string("Charge_Vs_PathlengthTECP1") + stag;
433  std::string cvpTECP2 = std::string("Charge_Vs_PathlengthTECP2") + stag;
434  std::string cvpTECM1 = std::string("Charge_Vs_PathlengthTECM1") + stag;
435  std::string cvpTECM2 = std::string("Charge_Vs_PathlengthTECM2") + stag;
436 
437  int elepos = (m_harvestingMode && AlgoMode=="PCL")? Harvest : statCollectionFromMode(tag);
438 
439  Charge_Vs_Index[elepos] = dbe->book2S(cvi.c_str() , cvi.c_str() , 88625, 0 , 88624,2000,0,4000);
440  //Charge_Vs_Index_Absolute[elepos] = dbe->book2S(cviA.c_str() , cviA.c_str() , 88625, 0 , 88624,1000,0,4000);
441  Charge_Vs_PathlengthTIB[elepos] = dbe->book2S(cvpTIB.c_str() , cvpTIB.c_str() , 20 , 0.3 , 1.3 , 250,0,2000);
442  Charge_Vs_PathlengthTOB[elepos] = dbe->book2S(cvpTOB.c_str() , cvpTOB.c_str() , 20 , 0.3 , 1.3 , 250,0,2000);
443  Charge_Vs_PathlengthTIDP[elepos] = dbe->book2S(cvpTIDP.c_str() , cvpTIDP.c_str() , 20 , 0.3 , 1.3 , 250,0,2000);
444  Charge_Vs_PathlengthTIDM[elepos] = dbe->book2S(cvpTIDM.c_str() , cvpTIDM.c_str() , 20 , 0.3 , 1.3 , 250,0,2000);
445  Charge_Vs_PathlengthTECP1[elepos] = dbe->book2S(cvpTECP1.c_str(), cvpTECP1.c_str(), 20 , 0.3 , 1.3 , 250,0,2000);
446  Charge_Vs_PathlengthTECP2[elepos] = dbe->book2S(cvpTECP2.c_str(), cvpTECP2.c_str(), 20 , 0.3 , 1.3 , 250,0,2000);
447  Charge_Vs_PathlengthTECM1[elepos] = dbe->book2S(cvpTECM1.c_str(), cvpTECM1.c_str(), 20 , 0.3 , 1.3 , 250,0,2000);
448  Charge_Vs_PathlengthTECM2[elepos] = dbe->book2S(cvpTECM2.c_str(), cvpTECM2.c_str(), 20 , 0.3 , 1.3 , 250,0,2000);
449 }
450 
451 void SiStripGainFromCalibTree::algoBeginJob(const edm::EventSetup& iSetup)
452 {
453  edm::LogInfo("SiStripGainFromCalibTree") << "AlgoMode : " << AlgoMode << "\n"
454  << "CalibrationMode : " << m_calibrationMode << "\n"
455  << "HarvestingMode : " << m_harvestingMode << std::endl;
456  //Setup DQM histograms
457  if(AlgoMode != "PCL" or m_harvestingMode) {
458  const char * dqm_dir = "AlCaReco/SiStripGainsHarvesting/";
459  this->bookDQMHistos( dqm_dir, dqm_tag_[statCollectionFromMode(m_calibrationMode.c_str())].c_str() );
460  } else {
461  //Check consistency of calibration Mode and BField only for the ALCAPROMPT in the PCL workflow
462  if (!isBFieldConsistentWithMode(iSetup)) {
463  string prevMode = m_calibrationMode;
464  swapBFieldMode();
465  edm::LogInfo("SiStripGainFromCalibTree") << "Switching calibration mode for endorsing BField status: "
466  << prevMode << " ==> " << m_calibrationMode << std::endl;
467  }
468  std::string dqm_dir = m_DQMdir + ((m_splitDQMstat)? m_calibrationMode:"") + "/";
469  int elem = statCollectionFromMode(m_calibrationMode.c_str());
470  this->bookDQMHistos( dqm_dir.c_str(), dqm_tag_[elem].c_str() );
471  }
472 
473 
474  edm::ESHandle<TrackerGeometry> tkGeom;
475  iSetup.get<TrackerDigiGeometryRecord>().get( tkGeom );
476  auto const & Det = tkGeom->dets();
477 
478  NPixelDets = 0;
479  NStripAPVs = 0;
480  unsigned int Index=0;
481  for(unsigned int i=0;i<Det.size();i++){
482  DetId Detid = Det[i]->geographicalId();
483  int SubDet = Detid.subdetId();
484 
485  if( SubDet == StripSubdetector::TIB || SubDet == StripSubdetector::TID ||
486  SubDet == StripSubdetector::TOB || SubDet == StripSubdetector::TEC ){
487 
488  auto DetUnit = dynamic_cast<const StripGeomDetUnit*> (Det[i]);
489  if(!DetUnit)continue;
490 
491  const StripTopology& Topo = DetUnit->specificTopology();
492  unsigned int NAPV = Topo.nstrips()/128;
493 
494  for(unsigned int j=0;j<NAPV;j++){
495  stAPVGain* APV = new stAPVGain;
496  APV->Index = Index;
497  APV->Bin = -1;
498  APV->DetId = Detid.rawId();
499  APV->APVId = j;
500  APV->SubDet = SubDet;
501  APV->FitMPV = -1;
502  APV->FitMPVErr = -1;
503  APV->FitWidth = -1;
504  APV->FitWidthErr = -1;
505  APV->FitChi2 = -1;
506  APV->FitNorm = -1;
507  APV->Gain = -1;
508  APV->PreviousGain = 1;
509  APV->PreviousGainTick = 1;
510  APV->x = DetUnit->position().basicVector().x();
511  APV->y = DetUnit->position().basicVector().y();
512  APV->z = DetUnit->position().basicVector().z();
513  APV->Eta = DetUnit->position().basicVector().eta();
514  APV->Phi = DetUnit->position().basicVector().phi();
515  APV->R = DetUnit->position().basicVector().transverse();
516  APV->Thickness = DetUnit->surface().bounds().thickness();
517  APV->NEntries = 0;
518  APV->isMasked = false;
519 
520  APVsCollOrdered.push_back(APV);
521  APVsColl[(APV->DetId<<4) | APV->APVId] = APV;
522  Index++;
523  NStripAPVs++;
524  }
525  }
526  }
527 
528  for(unsigned int i=0;i<Det.size();i++){ //Make two loop such that the Pixel information is added at the end --> make transition simpler
529  DetId Detid = Det[i]->geographicalId();
530  int SubDet = Detid.subdetId();
531  if( SubDet == PixelSubdetector::PixelBarrel || PixelSubdetector::PixelEndcap ){
532  auto DetUnit = dynamic_cast<const PixelGeomDetUnit*> (Det[i]);
533  if(!DetUnit) continue;
534 
535  const PixelTopology& Topo = DetUnit->specificTopology();
536  unsigned int NROCRow = Topo.nrows()/(80.);
537  unsigned int NROCCol = Topo.ncolumns()/(52.);
538 
539  for(unsigned int j=0;j<NROCRow;j++){
540  for(unsigned int i=0;i<NROCCol;i++){
541 
542  stAPVGain* APV = new stAPVGain;
543  APV->Index = Index;
544  APV->Bin = -1;
545  APV->DetId = Detid.rawId();
546  APV->APVId = (j<<3 | i);
547  APV->SubDet = SubDet;
548  APV->FitMPV = -1;
549  APV->FitMPVErr = -1;
550  APV->FitWidth = -1;
551  APV->FitWidthErr = -1;
552  APV->FitChi2 = -1;
553  APV->Gain = -1;
554  APV->PreviousGain = 1;
555  APV->x = DetUnit->position().basicVector().x();
556  APV->y = DetUnit->position().basicVector().y();
557  APV->z = DetUnit->position().basicVector().z();
558  APV->Eta = DetUnit->position().basicVector().eta();
559  APV->Phi = DetUnit->position().basicVector().phi();
560  APV->R = DetUnit->position().basicVector().transverse();
561  APV->Thickness = DetUnit->surface().bounds().thickness();
562  APV->isMasked = false; //SiPixelQuality_->IsModuleBad(Detid.rawId());
563  APV->NEntries = 0;
564 
565  APVsCollOrdered.push_back(APV);
566  APVsColl[(APV->DetId<<4) | APV->APVId] = APV;
567  Index++;
568  NPixelDets++;
569  }}
570  }
571  }
572 
573 
574  MakeCalibrationMap();
575 
576  NEvent = 0;
577  NTrack = 0;
578  NClusterStrip = 0;
579  NClusterPixel = 0;
580  SRun = 1<<31;
581  ERun = 0;
582  GOOD = 0;
583  BAD = 0;
584  MASKED = 0;
585 }
586 
587 bool SiStripGainFromCalibTree::isBFieldConsistentWithMode( const edm::EventSetup& iSetup) const {
588  edm::ESHandle<RunInfo> runInfo;
589  iSetup.get<RunInfoRcd>().get(runInfo);
590 
591  double average_current = runInfo.product()->m_avg_current;
592  bool isOn = (average_current > MagFieldCurrentTh);
593  bool is0T = (m_calibrationMode.substr( m_calibrationMode.length()-2 )=="0T");
594 
595  return ( (isOn && !is0T) || (!isOn && is0T) );
596 }
597 
598 void SiStripGainFromCalibTree::swapBFieldMode() {
599  if ( m_calibrationMode.substr( m_calibrationMode.length()-2 )=="0T" ) {
600  m_calibrationMode.erase( m_calibrationMode.length()-2,2);
601  } else {
602  m_calibrationMode.append( "0T" );
603  }
604 }
605 
606 void SiStripGainFromCalibTree::algoBeginRun(const edm::Run& run, const edm::EventSetup& iSetup)
607 {
608  if( !m_harvestingMode && AlgoMode=="PCL") {
609  //Check consistency of calibration Mode and BField only for the ALCAPROMPT in the PCL workflow
610  if (!isBFieldConsistentWithMode(iSetup)) {
611  string prevMode = m_calibrationMode;
612  swapBFieldMode();
613  edm::LogInfo("SiStripGainFromCalibTree") << "Switching calibration mode for endorsing BField status: "
614  << prevMode << " ==> " << m_calibrationMode << std::endl;
615  }
616  }
617 
618  edm::ESHandle<SiStripGain> gainHandle;
619  iSetup.get<SiStripGainRcd>().get(gainHandle);
620  if(!gainHandle.isValid()){edm::LogError("SiStripGainFromCalibTree")<< "gainHandle is not valid\n"; exit(0);}
621 
622  edm::ESHandle<SiStripQuality> SiStripQuality_;
623  iSetup.get<SiStripQualityRcd>().get(SiStripQuality_);
624 
625  for(unsigned int a=0;a<APVsCollOrdered.size();a++){
626  stAPVGain* APV = APVsCollOrdered[a];
627 
628  APV->isMasked = SiStripQuality_->IsApvBad(APV->DetId,APV->APVId);
629 // if(!FirstSetOfConstants){
630  if(gainHandle->getNumberOfTags()!=2){edm::LogError("SiStripGainFromCalibTree")<< "NUMBER OF GAIN TAG IS EXPECTED TO BE 2\n";fflush(stdout);exit(0);};
631  float newPreviousGain = gainHandle->getApvGain(APV->APVId,gainHandle->getRange(APV->DetId, 1),1);
632  if(APV->PreviousGain!=1 and newPreviousGain!=APV->PreviousGain)edm::LogWarning("SiStripGainFromCalibTree")<< "WARNING: ParticleGain in the global tag changed\n";
633  APV->PreviousGain = newPreviousGain;
634 
635  float newPreviousGainTick = gainHandle->getApvGain(APV->APVId,gainHandle->getRange(APV->DetId, 0),1);
636  if(APV->PreviousGainTick!=1 and newPreviousGainTick!=APV->PreviousGainTick)edm::LogWarning("SiStripGainFromCalibTree")<< "WARNING: TickMarkGain in the global tag changed\n";
637  APV->PreviousGainTick = newPreviousGainTick;
638 
639 
640  //printf("DETID = %7i APVID=%1i Previous Gain=%8.4f (G1) x %8.4f (G2)\n",APV->DetId,APV->APVId,APV->PreviousGainTick, APV->PreviousGain);
641 
642 
643 
644 // }
645  }
646 }
647 
648 void SiStripGainFromCalibTree::algoEndRun(const edm::Run& run, const edm::EventSetup& iSetup){
649  if(AlgoMode == "PCL" && !m_harvestingMode) return;//nothing to do in that case
650 
651  if(AlgoMode == "PCL" and m_harvestingMode){
652  // Load the 2D histograms from the DQM objects
653  // When running in AlCaHarvesting mode the histos are already booked and should be just retrieved from
654  // DQMStore so that they can be used in the fit
655 
656  edm::LogInfo("SiStripGainFromCalibTree") << "Starting harvesting statistics" << std::endl;
657 
658  // check the required tag before adding histograms
659  int elepos = statCollectionFromMode(m_calibrationMode.c_str());
660  if (elepos!=Harvest) {
661 
662  //collect statistics from DQM into the related monitored elements
663  std::string stag = m_calibrationMode;
664  if(stag.size()!=0 && stag[0]!='_') stag.insert(0,1,'_');
665 
666 
667  if (elepos==-1) {
668  //implememt backward compatibility
669  elepos=0;
670  stag="";
671  }
672 
673  std::string DQM_dir = m_DQMdir + ((m_splitDQMstat)? m_calibrationMode:"");
674 
675  std::string cvi = DQM_dir + std::string("/Charge_Vs_Index") + stag;
676  //std::string cviA = DQM_dir + std::string("/Charge_Vs_Index_Absolute") + stag;
677  std::string cvpTIB = DQM_dir + std::string("/Charge_Vs_PathlengthTIB") + stag;
678  std::string cvpTOB = DQM_dir + std::string("/Charge_Vs_PathlengthTOB") + stag;
679  std::string cvpTIDP = DQM_dir + std::string("/Charge_Vs_PathlengthTIDP") + stag;
680  std::string cvpTIDM = DQM_dir + std::string("/Charge_Vs_PathlengthTIDM") + stag;
681  std::string cvpTECP1 = DQM_dir + std::string("/Charge_Vs_PathlengthTECP1") + stag;
682  std::string cvpTECP2 = DQM_dir + std::string("/Charge_Vs_PathlengthTECP2") + stag;
683  std::string cvpTECM1 = DQM_dir + std::string("/Charge_Vs_PathlengthTECM1") + stag;
684  std::string cvpTECM2 = DQM_dir + std::string("/Charge_Vs_PathlengthTECM2") + stag;
685 
686  Charge_Vs_Index[elepos] = dbe->get(cvi.c_str());
687  //Charge_Vs_Index_Absolute[elepos] = dbe->get(cviA.c_str());
688  Charge_Vs_PathlengthTIB[elepos] = dbe->get(cvpTIB.c_str());
689  Charge_Vs_PathlengthTOB[elepos] = dbe->get(cvpTOB.c_str());
690  Charge_Vs_PathlengthTIDP[elepos] = dbe->get(cvpTIDP.c_str());
691  Charge_Vs_PathlengthTIDM[elepos] = dbe->get(cvpTIDM.c_str());
692  Charge_Vs_PathlengthTECP1[elepos] = dbe->get(cvpTECP1.c_str());
693  Charge_Vs_PathlengthTECP2[elepos] = dbe->get(cvpTECP2.c_str());
694  Charge_Vs_PathlengthTECM1[elepos] = dbe->get(cvpTECM1.c_str());
695  Charge_Vs_PathlengthTECM2[elepos] = dbe->get(cvpTECM2.c_str());
696 
697  if (Charge_Vs_Index[elepos]==0) {
698  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvi.c_str()
699  << ", statistics will not be summed!" << std::endl;
700  } else {
701  merge( (Charge_Vs_Index[Harvest])->getTH2S(), (Charge_Vs_Index[elepos])->getTH2S() );
702  edm::LogInfo("SiStripGainFromCalibTree") << "Harvesting "
703  << (Charge_Vs_Index[elepos])->getTH2S()->GetEntries() << " more clusters" << std::endl;
704  }
705 
706  //if (Charge_Vs_Index_Absolute[elepos]==0) {
707  // edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cviA.c_str()
708  // << ", statistics will not be summed!" << std::endl;
709  //} else merge( (Charge_Vs_Index_Absolute[Harvest])->getTH2S(), (Charge_Vs_Index_Absolute[elepos])->getTH2S() );
710 
711  if (Charge_Vs_PathlengthTIB[elepos]==0) {
712  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTIB.c_str()
713  << ", statistics will not be summed!" << std::endl;
714  } else (Charge_Vs_PathlengthTIB[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTIB[elepos])->getTH2S());
715 
716  if (Charge_Vs_PathlengthTOB[elepos]==0) {
717  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTOB.c_str()
718  << ", statistics will not be summed!" << std::endl;
719  } else (Charge_Vs_PathlengthTOB[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTOB[elepos])->getTH2S());
720 
721  if (Charge_Vs_PathlengthTIDP[elepos]==0) {
722  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTIDP.c_str()
723  << ", statistics will not be summed!" << std::endl;
724  } else (Charge_Vs_PathlengthTIDP[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTIDP[elepos])->getTH2S());
725 
726  if (Charge_Vs_PathlengthTIDM[elepos]==0) {
727  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTIDM.c_str()
728  << ", statistics will not be summed!" << std::endl;
729  } else (Charge_Vs_PathlengthTIDM[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTIDM[elepos])->getTH2S());
730 
731  if (Charge_Vs_PathlengthTECP1[elepos]==0) {
732  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTECP1.c_str()
733  << ", statistics will not be summed!" << std::endl;
734  } else (Charge_Vs_PathlengthTECP1[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECP1[elepos])->getTH2S());
735 
736  if (Charge_Vs_PathlengthTECP2[elepos]==0) {
737  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTECP2.c_str()
738  << ", statistics will not be summed!" << std::endl;
739  } else (Charge_Vs_PathlengthTECP2[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECP2[elepos])->getTH2S());
740 
741  if (Charge_Vs_PathlengthTECM1[elepos]==0) {
742  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTECM1.c_str()
743  << ", statistics will not be summed!" << std::endl;
744  } else (Charge_Vs_PathlengthTECM1[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECM1[elepos])->getTH2S());
745 
746  if (Charge_Vs_PathlengthTECM2[elepos]==0) {
747  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not retrieve " << cvpTECM2.c_str()
748  << ", statistics will not be summed!" << std::endl;
749  } else (Charge_Vs_PathlengthTECM2[Harvest])->getTH2S()->Add((Charge_Vs_PathlengthTECM2[elepos])->getTH2S());
750  }
751  }
752 }
753 
754 void
755 SiStripGainFromCalibTree::algoEndJob() {
756  if(AlgoMode == "PCL" && !m_harvestingMode) return;//nothing to do in that case
757 
758  if(AlgoMode == "CalibTree"){
759  edm::LogInfo("SiStripGainFromCalibTree") << "Analyzing calibration tree" << std::endl;
760  // Loop on calibTrees to fill the 2D histograms
761  algoAnalyzeTheTree();
762  }else if(m_harvestingMode){
763  NClusterStrip = (Charge_Vs_Index[Harvest])->getTH2S()->Integral(0,NStripAPVs+1, 0, 99999 );
764  NClusterPixel = (Charge_Vs_Index[Harvest])->getTH2S()->Integral(NStripAPVs+2, NStripAPVs+NPixelDets+2, 0, 99999 );
765  }
766 
767  // Now that we have the full statistics we can extract the information of the 2D histograms
768  algoComputeMPVandGain();
769 
770  if(AlgoMode != "PCL" or saveSummary){
771  edm::LogInfo("SiStripGainFromCalibTree") << "Saving summary into root file" << std::endl;
772 
773  //also save the 2D monitor elements to this file as TH2D tfs
774  tfs = edm::Service<TFileService>().operator->();
775 
776  //save only the statistics for the calibrationTag
777  int elepos = statCollectionFromMode(m_calibrationMode.c_str());
778 
779  if( Charge_Vs_Index[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_Index[elepos])->getTH2S() );
780  //if( Charge_Vs_Index_Absolute[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_Index_Absolute[elepos])->getTH2S() );
781  if( Charge_Vs_PathlengthTIB[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTIB[elepos])->getTH2S() );
782  if( Charge_Vs_PathlengthTOB[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTOB[elepos])->getTH2S() );
783  if( Charge_Vs_PathlengthTIDP[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTIDP[elepos])->getTH2S() );
784  if( Charge_Vs_PathlengthTIDM[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTIDM[elepos])->getTH2S() );
785  if( Charge_Vs_PathlengthTECP1[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTECP1[elepos])->getTH2S() );
786  if( Charge_Vs_PathlengthTECP2[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTECP2[elepos])->getTH2S() );
787  if( Charge_Vs_PathlengthTECM1[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTECM1[elepos])->getTH2S() );
788  if( Charge_Vs_PathlengthTECM2[elepos]!=0 ) tfs->make<TH2S> ( *(Charge_Vs_PathlengthTECM2[elepos])->getTH2S() );
789 
790  storeOnTree(tfs);
791  }
792 }
793 
794 
795 void SiStripGainFromCalibTree::getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange, double HighRange)
796 {
797  FitResults[0] = -0.5; //MPV
798  FitResults[1] = 0; //MPV error
799  FitResults[2] = -0.5; //Width
800  FitResults[3] = 0; //Width error
801  FitResults[4] = -0.5; //Fit Chi2/NDF
802  FitResults[5] = 0; //Normalization
803 
804  if( InputHisto->GetEntries() < MinNrEntries)return;
805 
806  // perform fit with standard landau
807  TF1* MyLandau = new TF1("MyLandau","landau",LowRange, HighRange);
808  MyLandau->SetParameter(1,300);
809  InputHisto->Fit(MyLandau,"0QR WW");
810 
811  // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
812  FitResults[0] = MyLandau->GetParameter(1); //MPV
813  FitResults[1] = MyLandau->GetParError(1); //MPV error
814  FitResults[2] = MyLandau->GetParameter(2); //Width
815  FitResults[3] = MyLandau->GetParError(2); //Width error
816  FitResults[4] = MyLandau->GetChisquare() / MyLandau->GetNDF(); //Fit Chi2/NDF
817  FitResults[5] = MyLandau->GetParameter(0);
818 
819  delete MyLandau;
820 }
821 
822 bool SiStripGainFromCalibTree::IsGoodLandauFit(double* FitResults){
823  if(FitResults[0] <= 0 )return false;
824 // if(FitResults[1] > MaxMPVError )return false;
825 // if(FitResults[4] > MaxChi2OverNDF)return false;
826  return true;
827 }
828 
829 void SiStripGainFromCalibTree::processEvent() {
830 
831  edm::LogInfo("SiStripGainFromCalibTree") << "Processing run " << runnumber
832  << " and event " << eventnumber
833  << " for " << m_calibrationMode << " calibration." << std::endl;
834 
835  if(runnumber<SRun)SRun=runnumber;
836  if(runnumber>ERun)ERun=runnumber;
837 
838  NEvent++;
839  NTrack+=(*trackp).size();
840 
841  int elepos = statCollectionFromMode(m_calibrationMode.c_str());
842 
843  unsigned int FirstAmplitude=0;
844  for(unsigned int i=0;i<(*chargeoverpath).size();i++){
845  FirstAmplitude+=(*nstrips)[i];
846  int TI = (*trackindex)[i];
847 
848  //printf("%i - %i - %i %i %i\n", (int)(*rawid)[i], (int)(*firststrip)[i]/128, (int)(*farfromedge)[i], (int)(*overlapping)[i], (int)(*saturation )[i] );
849  if((*tracketa )[TI] < MinTrackEta )continue;
850  if((*tracketa )[TI] > MaxTrackEta )continue;
851  if((*trackp )[TI] < MinTrackMomentum )continue;
852  if((*trackp )[TI] > MaxTrackMomentum )continue;
853  if((*trackhitsvalid)[TI] < MinTrackHits )continue;
854  if((*trackchi2ndof )[TI] > MaxTrackChiOverNdf )continue;
855  if((*trackalgo )[TI] > MaxTrackingIteration )continue;
856 
857  stAPVGain* APV = APVsColl[((*rawid)[i]<<4) | ((*firststrip)[i]/128)]; //works for both strip and pixel thanks to firstStrip encoding for pixel in the calibTree
858 
859  if(APV->SubDet>2 && (*farfromedge)[i] == false )continue;
860  if(APV->SubDet>2 && (*overlapping)[i] == true )continue;
861  if(APV->SubDet>2 && (*saturation )[i] && !AllowSaturation)continue;
862  if(APV->SubDet>2 && (*nstrips )[i] > MaxNrStrips )continue;
863 
864 
865  //printf("detId=%7i run=%7i event=%9i charge=%5i cs=%3i\n",(*rawid)[i],runnumber,eventnumber,(*charge)[i],(*nstrips)[i]);
866 
867  //double trans = atan2((*localdiry)[i],(*localdirx)[i])*(180/3.14159265);
868  //double alpha = acos ((*localdirx)[i] / sqrt( pow((*localdirx)[i],2) + pow((*localdirz)[i],2) ) ) * (180/3.14159265);
869  //double beta = acos ((*localdiry)[i] / sqrt( pow((*localdirx)[i],2) + pow((*localdirz)[i],2) ) ) * (180/3.14159265);
870 
871  //printf("NStrip = %i : Charge = %i --> Path = %f --> ChargeOverPath=%f\n",(*nstrips)[i],(*charge)[i],(*path)[i],(*chargeoverpath)[i]);
872  //printf("Amplitudes: ");
873  //for(unsigned int a=0;a<(*nstrips)[i];a++){printf("%i ",(*amplitude)[FirstAmplitude+a]);}
874  //printf("\n");
875 
876  if(APV->SubDet>2){NClusterStrip++;}else{NClusterPixel++;}
877 
878  int Charge = 0;
879  if(APV->SubDet>2 && (useCalibration || !FirstSetOfConstants)){
880  bool Saturation = false;
881  for(unsigned int s=0;s<(*nstrips)[i];s++){
882  int StripCharge = (*amplitude)[FirstAmplitude-(*nstrips)[i]+s];
883  if(useCalibration && !FirstSetOfConstants){ StripCharge=(int)(StripCharge*(APV->PreviousGain/APV->CalibGain));
884  }else if(useCalibration){ StripCharge=(int)(StripCharge/APV->CalibGain);
885  }else if(!FirstSetOfConstants){ StripCharge=(int)(StripCharge*APV->PreviousGain);}
886  if(StripCharge>1024){
887  StripCharge = 255;
888  Saturation = true;
889  }else if(StripCharge>254){
890  StripCharge = 254;
891  Saturation = true;
892  }
893  Charge += StripCharge;
894  }
895  if(Saturation && !AllowSaturation)continue;
896  }else if(APV->SubDet>2){
897  Charge = (*charge)[i];
898  }else{
899  Charge = (*charge)[i]/265.0; //expected scale factor between pixel and strip charge
900  }
901 
902  //printf("ChargeDifference = %i Vs %i with Gain = %f\n",(*charge)[i],Charge,APV->CalibGain);
903 
904  double ClusterChargeOverPath = ( (double) Charge )/(*path)[i] ;
905  if(APV->SubDet>2){
906  if(Validation) {ClusterChargeOverPath/=(*gainused)[i];}
907  if(OldGainRemoving){ClusterChargeOverPath*=(*gainused)[i];}
908  }
909  //(Charge_Vs_Index_Absolute[elepos])->Fill(APV->Index,Charge);
910  (Charge_Vs_Index[elepos]) ->Fill(APV->Index,ClusterChargeOverPath);
911 
912  if(APV->SubDet==StripSubdetector::TIB){ (Charge_Vs_PathlengthTIB[elepos]) ->Fill((*path)[i],Charge);
913  }else if(APV->SubDet==StripSubdetector::TOB){ (Charge_Vs_PathlengthTOB[elepos]) ->Fill((*path)[i],Charge);
914  }else if(APV->SubDet==StripSubdetector::TID){
915  if(APV->Eta<0){ (Charge_Vs_PathlengthTIDM[elepos]) ->Fill((*path)[i],Charge);
916  }else if(APV->Eta>0){ (Charge_Vs_PathlengthTIDP[elepos]) ->Fill((*path)[i],Charge);
917  }
918  }else if(APV->SubDet==StripSubdetector::TEC){
919  if(APV->Eta<0){
920  if(APV->Thickness<0.04){ (Charge_Vs_PathlengthTECM1[elepos])->Fill((*path)[i],Charge);
921  }else if(APV->Thickness>0.04){ (Charge_Vs_PathlengthTECM2[elepos])->Fill((*path)[i],Charge);
922  }
923  }else if(APV->Eta>0){
924  if(APV->Thickness<0.04){ (Charge_Vs_PathlengthTECP1[elepos])->Fill((*path)[i],Charge);
925  }else if(APV->Thickness>0.04){ (Charge_Vs_PathlengthTECP2[elepos])->Fill((*path)[i],Charge);
926  }
927  }
928  }
929 
930  }// END OF ON-CLUSTER LOOP
931 }//END OF processEvent()
932 
933 void SiStripGainFromCalibTree::algoAnalyzeTheTree()
934 {
935  for(unsigned int i=0;i<VInputFiles.size();i++){
936  printf("Openning file %3i/%3i --> %s\n",i+1, (int)VInputFiles.size(), (char*)(VInputFiles[i].c_str())); fflush(stdout);
937  TFile *tfile = TFile::Open(VInputFiles[i].c_str());
938  TString tree_path = TString::Format("gainCalibrationTree%s/tree",m_calibrationMode.c_str());
939  TTree* tree = dynamic_cast<TTree*> (tfile->Get(tree_path.Data()));
940 
941  tree->SetBranchAddress((EventPrefix_ + "event" + EventSuffix_).c_str(), &eventnumber , NULL);
942  tree->SetBranchAddress((EventPrefix_ + "run" + EventSuffix_).c_str(), &runnumber , NULL);
943  tree->SetBranchAddress((EventPrefix_ + "TrigTech" + EventSuffix_).c_str(), &TrigTech , NULL);
944 
945  tree->SetBranchAddress((TrackPrefix_ + "chi2ndof" + TrackSuffix_).c_str(), &trackchi2ndof , NULL);
946  tree->SetBranchAddress((TrackPrefix_ + "momentum" + TrackSuffix_).c_str(), &trackp , NULL);
947  tree->SetBranchAddress((TrackPrefix_ + "pt" + TrackSuffix_).c_str(), &trackpt , NULL);
948  tree->SetBranchAddress((TrackPrefix_ + "eta" + TrackSuffix_).c_str(), &tracketa , NULL);
949  tree->SetBranchAddress((TrackPrefix_ + "phi" + TrackSuffix_).c_str(), &trackphi , NULL);
950  tree->SetBranchAddress((TrackPrefix_ + "hitsvalid" + TrackSuffix_).c_str(), &trackhitsvalid, NULL);
951  tree->SetBranchAddress((TrackPrefix_ + "algo" + TrackSuffix_).c_str(), &trackalgo , NULL);
952 
953  tree->SetBranchAddress((CalibPrefix_ + "trackindex" + CalibSuffix_).c_str(), &trackindex , NULL);
954  tree->SetBranchAddress((CalibPrefix_ + "rawid" + CalibSuffix_).c_str(), &rawid , NULL);
955  tree->SetBranchAddress((CalibPrefix_ + "localdirx" + CalibSuffix_).c_str(), &localdirx , NULL);
956  tree->SetBranchAddress((CalibPrefix_ + "localdiry" + CalibSuffix_).c_str(), &localdiry , NULL);
957  tree->SetBranchAddress((CalibPrefix_ + "localdirz" + CalibSuffix_).c_str(), &localdirz , NULL);
958  tree->SetBranchAddress((CalibPrefix_ + "firststrip" + CalibSuffix_).c_str(), &firststrip , NULL);
959  tree->SetBranchAddress((CalibPrefix_ + "nstrips" + CalibSuffix_).c_str(), &nstrips , NULL);
960  tree->SetBranchAddress((CalibPrefix_ + "saturation" + CalibSuffix_).c_str(), &saturation , NULL);
961  tree->SetBranchAddress((CalibPrefix_ + "overlapping" + CalibSuffix_).c_str(), &overlapping , NULL);
962  tree->SetBranchAddress((CalibPrefix_ + "farfromedge" + CalibSuffix_).c_str(), &farfromedge , NULL);
963  tree->SetBranchAddress((CalibPrefix_ + "charge" + CalibSuffix_).c_str(), &charge , NULL);
964  tree->SetBranchAddress((CalibPrefix_ + "path" + CalibSuffix_).c_str(), &path , NULL);
965  tree->SetBranchAddress((CalibPrefix_ + "chargeoverpath" + CalibSuffix_).c_str(), &chargeoverpath, NULL);
966  tree->SetBranchAddress((CalibPrefix_ + "amplitude" + CalibSuffix_).c_str(), &amplitude , NULL);
967  tree->SetBranchAddress((CalibPrefix_ + "gainused" + CalibSuffix_).c_str(), &gainused , NULL);
968 
969 
970  unsigned int nentries = tree->GetEntries();
971  printf("Number of Events = %i + %i = %i\n",NEvent,nentries,(NEvent+nentries));
972  printf("Progressing Bar :0%% 20%% 40%% 60%% 80%% 100%%\n");
973  printf("Looping on the Tree :");
974  int TreeStep = nentries/50;if(TreeStep<=1)TreeStep=1;
975  for (unsigned int ientry = 0; ientry < tree->GetEntries(); ientry++) {
976  if(ientry%TreeStep==0){printf(".");fflush(stdout);}
977  tree->GetEntry(ientry);
978  processEvent();
979  }printf("\n");// END OF EVENT LOOP
980  }
981 }
982 
983 
984 
985 void SiStripGainFromCalibTree::algoComputeMPVandGain() {
986  unsigned int I=0;
987  TH1F* Proj = NULL;
988  double FitResults[6];
989  double MPVmean = 300;
990 
991  int elepos = (AlgoMode == "PCL")? Harvest : statCollectionFromMode(m_calibrationMode.c_str());
992 
993  if ( Charge_Vs_Index[elepos]==0 ) {
994  edm::LogError("SiStripGainFromCalibTree") << "Harvesting: could not execute algoComputeMPVandGain method because "
995  << m_calibrationMode.c_str() << " statistics cannot be retrieved.\n"
996  << "Please check if input contains "
997  << m_calibrationMode.c_str() << " data." << std::endl;
998  return;
999  }
1000 
1001  TH2S *chvsidx = (Charge_Vs_Index[elepos])->getTH2S();
1002 
1003 
1004  printf("Progressing Bar :0%% 20%% 40%% 60%% 80%% 100%%\n");
1005  printf("Fitting Charge Distribution :");
1006  int TreeStep = APVsColl.size()/50;
1007  for(__gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual >::iterator it = APVsColl.begin();it!=APVsColl.end();it++,I++){
1008  if(I%TreeStep==0){printf(".");fflush(stdout);}
1009  stAPVGain* APV = it->second;
1010  if(APV->Bin<0) APV->Bin = chvsidx->GetXaxis()->FindBin(APV->Index);
1011 
1012  if(APV->isMasked){APV->Gain=APV->PreviousGain; MASKED++; continue;}
1013 
1014  Proj = (TH1F*)(chvsidx->ProjectionY("",chvsidx->GetXaxis()->FindBin(APV->Index),chvsidx->GetXaxis()->FindBin(APV->Index),"e"));
1015  if(!Proj)continue;
1016 
1017  if(CalibrationLevel==0){
1018  }else if(CalibrationLevel==1){
1019  int SecondAPVId = APV->APVId;
1020  if(SecondAPVId%2==0){ SecondAPVId = SecondAPVId+1; }else{ SecondAPVId = SecondAPVId-1; }
1021  stAPVGain* APV2 = APVsColl[(APV->DetId<<4) | SecondAPVId];
1022  if(APV2->Bin<0) APV2->Bin = chvsidx->GetXaxis()->FindBin(APV2->Index);
1023  TH1F* Proj2 = (TH1F*)(chvsidx->ProjectionY("",APV2->Bin,APV2->Bin,"e"));
1024  if(Proj2){Proj->Add(Proj2,1);delete Proj2;}
1025  }else if(CalibrationLevel==2){
1026  for(unsigned int i=0;i<16;i++){ //loop up to 6APV for Strip and up to 16 for Pixels
1027  __gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual >::iterator tmpit;
1028  tmpit = APVsColl.find((APV->DetId<<4) | i);
1029  if(tmpit==APVsColl.end())continue;
1030  stAPVGain* APV2 = tmpit->second;
1031  if(APV2->DetId != APV->DetId || APV2->APVId == APV->APVId)continue;
1032  if(APV2->Bin<0) APV2->Bin = chvsidx->GetXaxis()->FindBin(APV2->Index);
1033  TH1F* Proj2 = (TH1F*)(chvsidx->ProjectionY("",APV2->Bin,APV2->Bin,"e"));
1034  if(Proj2){Proj->Add(Proj2,1);delete Proj2;}
1035  }
1036  }else{
1037  CalibrationLevel = 0;
1038  printf("Unknown Calibration Level, will assume %i\n",CalibrationLevel);
1039  }
1040 
1041  getPeakOfLandau(Proj,FitResults);
1042  APV->FitMPV = FitResults[0];
1043  APV->FitMPVErr = FitResults[1];
1044  APV->FitWidth = FitResults[2];
1045  APV->FitWidthErr = FitResults[3];
1046  APV->FitChi2 = FitResults[4];
1047  APV->FitNorm = FitResults[5];
1048  APV->NEntries = Proj->GetEntries();
1049 
1050  if(IsGoodLandauFit(FitResults)){
1051  APV->Gain = APV->FitMPV / MPVmean;
1052  if(APV->SubDet>2)GOOD++;
1053  }else{
1054  APV->Gain = APV->PreviousGain;
1055  if(APV->SubDet>2)BAD++;
1056  }
1057  if(APV->Gain<=0) APV->Gain = 1;
1058 
1059  //printf("%5i/%5i: %6i - %1i %5E Entries --> MPV = %f +- %f\n",I,APVsColl.size(),APV->DetId, APV->APVId, Proj->GetEntries(), FitResults[0], FitResults[1]);fflush(stdout);
1060  delete Proj;
1061  }printf("\n");
1062 }
1063 
1064 
1065 void SiStripGainFromCalibTree::storeOnTree(TFileService* tfs)
1066 {
1067  unsigned int tree_Index;
1068  unsigned int tree_Bin;
1069  unsigned int tree_DetId;
1070  unsigned char tree_APVId;
1071  unsigned char tree_SubDet;
1072  float tree_x;
1073  float tree_y;
1074  float tree_z;
1075  float tree_Eta;
1076  float tree_R;
1077  float tree_Phi;
1078  float tree_Thickness;
1079  float tree_FitMPV;
1080  float tree_FitMPVErr;
1081  float tree_FitWidth;
1082  float tree_FitWidthErr;
1083  float tree_FitChi2NDF;
1084  float tree_FitNorm;
1085  double tree_Gain;
1086  double tree_PrevGain;
1087  double tree_PrevGainTick;
1088  double tree_NEntries;
1089  bool tree_isMasked;
1090 
1091  TTree* MyTree;
1092  MyTree = tfs->make<TTree> ("APVGain","APVGain");
1093  MyTree->Branch("Index" ,&tree_Index ,"Index/i");
1094  MyTree->Branch("Bin" ,&tree_Bin ,"Bin/i");
1095  MyTree->Branch("DetId" ,&tree_DetId ,"DetId/i");
1096  MyTree->Branch("APVId" ,&tree_APVId ,"APVId/b");
1097  MyTree->Branch("SubDet" ,&tree_SubDet ,"SubDet/b");
1098  MyTree->Branch("x" ,&tree_x ,"x/F");
1099  MyTree->Branch("y" ,&tree_y ,"y/F");
1100  MyTree->Branch("z" ,&tree_z ,"z/F");
1101  MyTree->Branch("Eta" ,&tree_Eta ,"Eta/F");
1102  MyTree->Branch("R" ,&tree_R ,"R/F");
1103  MyTree->Branch("Phi" ,&tree_Phi ,"Phi/F");
1104  MyTree->Branch("Thickness" ,&tree_Thickness ,"Thickness/F");
1105  MyTree->Branch("FitMPV" ,&tree_FitMPV ,"FitMPV/F");
1106  MyTree->Branch("FitMPVErr" ,&tree_FitMPVErr ,"FitMPVErr/F");
1107  MyTree->Branch("FitWidth" ,&tree_FitWidth ,"FitWidth/F");
1108  MyTree->Branch("FitWidthErr" ,&tree_FitWidthErr,"FitWidthErr/F");
1109  MyTree->Branch("FitChi2NDF" ,&tree_FitChi2NDF ,"FitChi2NDF/F");
1110  MyTree->Branch("FitNorm" ,&tree_FitNorm ,"FitNorm/F");
1111  MyTree->Branch("Gain" ,&tree_Gain ,"Gain/D");
1112  MyTree->Branch("PrevGain" ,&tree_PrevGain ,"PrevGain/D");
1113  MyTree->Branch("PrevGainTick" ,&tree_PrevGainTick,"PrevGainTick/D");
1114  MyTree->Branch("NEntries" ,&tree_NEntries ,"NEntries/D");
1115  MyTree->Branch("isMasked" ,&tree_isMasked ,"isMasked/O");
1116 
1117 
1118  FILE* Gains = stdout;
1119  fprintf(Gains,"NEvents = %i\n",NEvent);
1120  fprintf(Gains,"NTracks = %i\n",NTrack);
1121  fprintf(Gains,"NClustersPixel = %i\n",NClusterPixel);
1122  fprintf(Gains,"NClustersStrip = %i\n",NClusterStrip);
1123  fprintf(Gains,"Number of Pixel Dets = %lu\n",static_cast<unsigned long>(NPixelDets));
1124  fprintf(Gains,"Number of Strip APVs = %lu\n",static_cast<unsigned long>(NStripAPVs));
1125  fprintf(Gains,"GoodFits = %i BadFits = %i ratio = %f%% (MASKED=%i)\n",GOOD,BAD,(100.0*GOOD)/(GOOD+BAD), MASKED);
1126 
1127  Gains=fopen(OutputGains.c_str(),"w");
1128  fprintf(Gains,"NEvents = %i\n",NEvent);
1129  fprintf(Gains,"NTracks = %i\n",NTrack);
1130  fprintf(Gains,"NClustersPixel = %i\n",NClusterPixel);
1131  fprintf(Gains,"NClustersStrip = %i\n",NClusterStrip);
1132  fprintf(Gains,"Number of Strip APVs = %lu\n",static_cast<unsigned long>(NStripAPVs));
1133  fprintf(Gains,"Number of Pixel Dets = %lu\n",static_cast<unsigned long>(NPixelDets));
1134  fprintf(Gains,"GoodFits = %i BadFits = %i ratio = %f%% (MASKED=%i)\n",GOOD,BAD,(100.0*GOOD)/(GOOD+BAD), MASKED);
1135 
1136  int elepos = statCollectionFromMode(m_calibrationMode.c_str());
1137 
1138  for(unsigned int a=0;a<APVsCollOrdered.size();a++){
1139  stAPVGain* APV = APVsCollOrdered[a];
1140  if(APV==NULL)continue;
1141 // printf( "%i | %i | PreviousGain = %7.5f NewGain = %7.5f (#clusters=%8.0f)\n", APV->DetId,APV->APVId,APV->PreviousGain,APV->Gain, APV->NEntries);
1142  fprintf(Gains,"%i | %i | PreviousGain = %7.5f(tick) x %7.5f(particle) NewGain (particle) = %7.5f (#clusters=%8.0f)\n", APV->DetId,APV->APVId,APV->PreviousGainTick, APV->PreviousGain,APV->Gain, APV->NEntries);
1143 
1144  tree_Index = APV->Index;
1145  tree_Bin = (Charge_Vs_Index[elepos])->getTH2S()->GetXaxis()->FindBin(APV->Index);
1146  tree_DetId = APV->DetId;
1147  tree_APVId = APV->APVId;
1148  tree_SubDet = APV->SubDet;
1149  tree_x = APV->x;
1150  tree_y = APV->y;
1151  tree_z = APV->z;
1152  tree_Eta = APV->Eta;
1153  tree_R = APV->R;
1154  tree_Phi = APV->Phi;
1155  tree_Thickness = APV->Thickness;
1156  tree_FitMPV = APV->FitMPV;
1157  tree_FitMPVErr = APV->FitMPVErr;
1158  tree_FitWidth = APV->FitWidth;
1159  tree_FitWidthErr= APV->FitWidthErr;
1160  tree_FitChi2NDF = APV->FitChi2;
1161  tree_FitNorm = APV->FitNorm;
1162  tree_Gain = APV->Gain;
1163  tree_PrevGain = APV->PreviousGain;
1164  tree_PrevGainTick = APV->PreviousGainTick;
1165  tree_NEntries = APV->NEntries;
1166  tree_isMasked = APV->isMasked;
1167 
1168 
1169  if(tree_DetId==402673324){
1170  printf("%i | %i : %f --> %f (%f)\n", tree_DetId, tree_APVId, tree_PrevGain, tree_Gain, tree_NEntries);
1171  }
1172 
1173 
1174  MyTree->Fill();
1175  }
1176  if(Gains)fclose(Gains);
1177 
1178 
1179 }
1180 
1181 bool SiStripGainFromCalibTree::produceTagFilter(){
1182 
1183  // The goal of this function is to check wether or not there is enough statistics
1184  // to produce a meaningful tag for the DB
1185  int elepos = (AlgoMode == "PCL")? Harvest : statCollectionFromMode(m_calibrationMode.c_str());
1186  if( Charge_Vs_Index[elepos]==0 ) {
1187  edm::LogError("SiStripGainFromCalibTree") << "produceTagFilter -> Return false: could not retrieve the "
1188  << m_calibrationMode.c_str() << " statistics.\n"
1189  << "Please check if input contains "
1190  << m_calibrationMode.c_str() << " data." << std::endl;
1191  return false;
1192  }
1193 
1194 
1195  float integral = (Charge_Vs_Index[elepos])->getTH2S()->Integral();
1196  if( (Charge_Vs_Index[elepos])->getTH2S()->Integral(0,NStripAPVs+1, 0, 99999 ) < tagCondition_NClusters) {
1197  edm::LogWarning("SiStripGainFromCalibTree")
1198  << "calibrationMode -> " << m_calibrationMode << "\n"
1199  << "produceTagFilter -> Return false: Statistics is too low : " << integral << endl;
1200  return false;
1201  }
1202  if((1.0 * GOOD) / (GOOD+BAD) < tagCondition_GoodFrac) {
1203  edm::LogWarning("SiStripGainFromCalibTree")
1204  << "calibrationMode -> " << m_calibrationMode << "\n"
1205  << "produceTagFilter -> Return false: ratio of GOOD/TOTAL is too low: " << (1.0 * GOOD) / (GOOD+BAD) << endl;
1206  return false;
1207  }
1208  return true;
1209 }
1210 
1211 SiStripApvGain* SiStripGainFromCalibTree::getNewObject()
1212 {
1213  SiStripApvGain* obj = new SiStripApvGain();
1214  if(!m_harvestingMode) return obj;
1215 
1216  if(!produceTagFilter()){
1217  edm::LogWarning("SiStripGainFromCalibTree")<< "getNewObject -> will not produce a paylaod because produceTagFilter returned false " << endl;
1218  setDoStore(false);
1219  return obj;
1220  }
1221 
1222 
1223  std::vector<float>* theSiStripVector = NULL;
1224  unsigned int PreviousDetId = 0;
1225  for(unsigned int a=0;a<APVsCollOrdered.size();a++){
1226  stAPVGain* APV = APVsCollOrdered[a];
1227  if(APV==NULL){ printf("Bug\n"); continue; }
1228  if(APV->SubDet<=2)continue;
1229  if(APV->DetId != PreviousDetId){
1230  if(theSiStripVector!=NULL){
1231  SiStripApvGain::Range range(theSiStripVector->begin(),theSiStripVector->end());
1232  if ( !obj->put(PreviousDetId,range) ) printf("Bug to put detId = %i\n",PreviousDetId);
1233  }
1234  theSiStripVector = new std::vector<float>;
1235  PreviousDetId = APV->DetId;
1236  }
1237  theSiStripVector->push_back(APV->Gain);
1238  }
1239  if(theSiStripVector!=NULL){
1240  SiStripApvGain::Range range(theSiStripVector->begin(),theSiStripVector->end());
1241  if ( !obj->put(PreviousDetId,range) ) printf("Bug to put detId = %i\n",PreviousDetId);
1242  }
1243 
1244  return obj;
1245 }
1246 
1247 
1248 SiStripGainFromCalibTree::~SiStripGainFromCalibTree()
1249 {
1250 }
1251 
1252 void SiStripGainFromCalibTree::MakeCalibrationMap(){
1253  if(!useCalibration)return;
1254 
1255  TChain* t1 = new TChain("SiStripCalib/APVGain");
1256  t1->Add(m_calibrationPath.c_str());
1257 
1258  unsigned int tree_DetId;
1259  unsigned char tree_APVId;
1260  double tree_Gain;
1261 
1262  t1->SetBranchAddress("DetId" ,&tree_DetId );
1263  t1->SetBranchAddress("APVId" ,&tree_APVId );
1264  t1->SetBranchAddress("Gain" ,&tree_Gain );
1265 
1266  for (unsigned int ientry = 0; ientry < t1->GetEntries(); ientry++) {
1267  t1->GetEntry(ientry);
1268  stAPVGain* APV = APVsColl[(tree_DetId<<4) | (unsigned int)tree_APVId];
1269  APV->CalibGain = tree_Gain;
1270  }
1271 }
1272 
1273 void
1274 SiStripGainFromCalibTree::algoAnalyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
1275 {
1276  // in AlCaHarvesting mode we just need to run the logic in the endJob step
1277  if(m_harvestingMode) return;
1278 
1279  if(AlgoMode=="CalibTree")return;
1280 
1281  eventnumber = iEvent.id().event();
1282  runnumber = iEvent.id().run();
1283  auto handle01 = connect(TrigTech , TrigTech_token_ , iEvent);
1284  auto handle02 = connect(trackchi2ndof , trackchi2ndof_token_ , iEvent);
1285  auto handle03 = connect(trackp , trackp_token_ , iEvent);
1286  auto handle04 = connect(trackpt , trackpt_token_ , iEvent);
1287  auto handle05 = connect(tracketa , tracketa_token_ , iEvent);
1288  auto handle06 = connect(trackphi , trackphi_token_ , iEvent);
1289  auto handle07 = connect(trackhitsvalid, trackhitsvalid_token_, iEvent);
1290  auto handle08 = connect(trackindex , trackindex_token_ , iEvent);
1291  auto handle09 = connect(rawid , rawid_token_ , iEvent);
1292  auto handle11 = connect(localdirx , localdirx_token_ , iEvent);
1293  auto handle12 = connect(localdiry , localdiry_token_ , iEvent);
1294  auto handle13 = connect(localdirz , localdirz_token_ , iEvent);
1295  auto handle14 = connect(firststrip , firststrip_token_ , iEvent);
1296  auto handle15 = connect(nstrips , nstrips_token_ , iEvent);
1297  auto handle16 = connect(saturation , saturation_token_ , iEvent);
1298  auto handle17 = connect(overlapping , overlapping_token_ , iEvent);
1299  auto handle18 = connect(farfromedge , farfromedge_token_ , iEvent);
1300  auto handle19 = connect(charge , charge_token_ , iEvent);
1301  auto handle21 = connect(path , path_token_ , iEvent);
1302  auto handle22 = connect(chargeoverpath, chargeoverpath_token_, iEvent);
1303  auto handle23 = connect(amplitude , amplitude_token_ , iEvent);
1304  auto handle24 = connect(gainused , gainused_token_ , iEvent);
1305 
1306  auto handle25 = connect(trackalgo , trackalgo_token_ , iEvent);
1307 
1308  processEvent();
1309 }
1310 
1311 DEFINE_FWK_MODULE(SiStripGainFromCalibTree);
edm::EventID::run
RunNumber_t run() const
Definition: EventID.h:39
SiStripGainFromCalibTree::trackindex_token_
edm::EDGetTokenT< std::vector< int > > trackindex_token_
Definition: SiStripGainFromCalibTree.cc:240
SiStripGainFromCalibTree::Charge_Vs_Index
std::vector< MonitorElement * > Charge_Vs_Index
Definition: SiStripGainFromCalibTree.cc:200
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::EventID::event
EventNumber_t event() const
Definition: EventID.h:41
StripTopology::nstrips
virtual int nstrips() const =0
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
i
int i
Definition: DBlmapReader.cc:9
SiStripGainFromCalibTree::trackp
const std::vector< float > * trackp
Definition: SiStripGainFromCalibTree.cc:232
SiStripGainFromCalibTree::tracketa_token_
edm::EDGetTokenT< std::vector< double > > tracketa_token_
Definition: SiStripGainFromCalibTree.cc:234
SiStripGainFromCalibTree::amplitude_token_
edm::EDGetTokenT< std::vector< unsigned char > > amplitude_token_
Definition: SiStripGainFromCalibTree.cc:253
SiStripGainFromCalibTree::chargeoverpath
const std::vector< double > * chargeoverpath
Definition: SiStripGainFromCalibTree.cc:252
SiStripGainFromCalibTree::rawid_token_
edm::EDGetTokenT< std::vector< unsigned int > > rawid_token_
Definition: SiStripGainFromCalibTree.cc:241
SiStripGainFromCalibTree::amplitude
const std::vector< unsigned char > * amplitude
Definition: SiStripGainFromCalibTree.cc:253
SiStripGainFromCalibTree::isBFieldConsistentWithMode
bool isBFieldConsistentWithMode(const edm::EventSetup &iSetup) const
Definition: SiStripGainFromCalibTree.cc:587
or
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e< void, edm::EventIDconst &, edm::Timestampconst & > We also list in braces which AR_WATCH_USING_METHOD_ is used for those or
Definition: Activities.doc:12
SiStripGainFromCalibTree::getNewObject
SiStripApvGain * getNewObject() override
Definition: SiStripGainFromCalibTree.cc:1211
PixelTopology::ncolumns
virtual int ncolumns() const =0
SiStripGainFromCalibTree::localdirz
const std::vector< double > * localdirz
Definition: SiStripGainFromCalibTree.cc:244
edm::Event::getByToken
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:462
SiStripGainFromCalibTree::trackchi2ndof
const std::vector< double > * trackchi2ndof
Definition: SiStripGainFromCalibTree.cc:231
SiStripGainFromCalibTree::algoBeginRun
virtual void algoBeginRun(const edm::Run &run, const edm::EventSetup &iSetup) override
Definition: SiStripGainFromCalibTree.cc:606
SiStripGainFromCalibTree::localdirx
const std::vector< double > * localdirx
Definition: SiStripGainFromCalibTree.cc:242
AlcaSiStripGainsAfterAbortGapHarvester_cff.DQMStore
tuple DQMStore
Definition: AlcaSiStripGainsAfterAbortGapHarvester_cff.py:11
processEvent
dispatcher processEvent(e, inputTag, standby)
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
SiStripGainFromCalibTree::algoEndRun
virtual void algoEndRun(const edm::Run &run, const edm::EventSetup &iSetup) override
Definition: SiStripGainFromCalibTree.cc:648
SiStripGainFromCalibTree::rawid
const std::vector< unsigned int > * rawid
Definition: SiStripGainFromCalibTree.cc:241
SiStripGainFromCalibTree::overlapping_token_
edm::EDGetTokenT< std::vector< bool > > overlapping_token_
Definition: SiStripGainFromCalibTree.cc:248
SiStripGainFromCalibTree::SiStripGainFromCalibTree
SiStripGainFromCalibTree(const edm::ParameterSet &)
Definition: SiStripGainFromCalibTree.cc:300
SiStripGainFromCalibTree::connect
edm::Handle< T > connect(const T *&ptr, edm::EDGetTokenT< T > token, const edm::Event &evt)
Definition: SiStripGainFromCalibTree.cc:152
PixelTopology::nrows
virtual int nrows() const =0
SiStripGainFromCalibTree::trackalgo
const std::vector< int > * trackalgo
Definition: SiStripGainFromCalibTree.cc:237
SiStripGainFromCalibTree::Charge_Vs_PathlengthTOB
std::vector< MonitorElement * > Charge_Vs_PathlengthTOB
Definition: SiStripGainFromCalibTree.cc:203
NULL
#define NULL
Definition: scimark2.h:8
SiStripGainFromCalibTree::storeOnTree
void storeOnTree(TFileService *tfs)
Definition: SiStripGainFromCalibTree.cc:1065
TFileService::make
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
SiStripGainFromCalibTree::trackpt
const std::vector< float > * trackpt
Definition: SiStripGainFromCalibTree.cc:233
SiStripGainFromCalibTree::Charge_Vs_PathlengthTECM1
std::vector< MonitorElement * > Charge_Vs_PathlengthTECM1
Definition: SiStripGainFromCalibTree.cc:208
dt_dqm_sourceclient_common_cff.reco
tuple reco
Definition: dt_dqm_sourceclient_common_cff.py:107
SiStripGainFromCalibTree::trackhitsvalid
const std::vector< unsigned int > * trackhitsvalid
Definition: SiStripGainFromCalibTree.cc:236
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
SiStripGainFromCalibTree::Charge_Vs_PathlengthTECM2
std::vector< MonitorElement * > Charge_Vs_PathlengthTECM2
Definition: SiStripGainFromCalibTree.cc:209
SiStripGainFromCalibTree::farfromedge_token_
edm::EDGetTokenT< std::vector< bool > > farfromedge_token_
Definition: SiStripGainFromCalibTree.cc:249
SiStripApvGain::put
bool put(const uint32_t &detID, Range input)
Definition: SiStripApvGain.cc:8
double_binary
Definition: DDStreamer.cc:227
DetId::rawId
uint32_t rawId() const
get the raw id
Definition: DetId.h:43
SiStripGainFromCalibTree::trackhitsvalid_token_
edm::EDGetTokenT< std::vector< unsigned int > > trackhitsvalid_token_
Definition: SiStripGainFromCalibTree.cc:236
HLT_25ns10e33_v2_cff.InputTag
tuple InputTag
Definition: HLT_25ns10e33_v2_cff.py:58467
x
T x() const
Cartesian x coordinate.
Definition: Basic3DVectorLD.h:127
SiStripGainFromCalibTree::APVsCollOrdered
std::vector< stAPVGain * > APVsCollOrdered
Definition: SiStripGainFromCalibTree.cc:269
SiStripGainFromCalibTree::APVsColl
__gnu_cxx::hash_map< unsigned int, stAPVGain *, __gnu_cxx::hash< unsigned int >, isEqual > APVsColl
Definition: SiStripGainFromCalibTree.cc:270
SiStripGainFromCalibTree::path_token_
edm::EDGetTokenT< std::vector< double > > path_token_
Definition: SiStripGainFromCalibTree.cc:251
iEvent
int iEvent
Definition: GenABIO.cc:230
SiStripGainFromCalibTree::trackalgo_token_
edm::EDGetTokenT< std::vector< int > > trackalgo_token_
Definition: SiStripGainFromCalibTree.cc:237
GlobalPosition_Frontier_DevDB_cff.tag
tuple tag
Definition: GlobalPosition_Frontier_DevDB_cff.py:11
SiStripGainFromCalibTree::trackp_token_
edm::EDGetTokenT< std::vector< float > > trackp_token_
Definition: SiStripGainFromCalibTree.cc:232
SiStripGainFromCalibTree::IsGoodLandauFit
bool IsGoodLandauFit(double *FitResults)
Definition: SiStripGainFromCalibTree.cc:822
HcalObjRepresent::Fill
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
Definition: HcalObjRepresent.h:611
patZpeak.handle
tuple handle
Definition: patZpeak.py:22
SiStripGainFromCalibTree::algoBeginJob
virtual void algoBeginJob(const edm::EventSetup &iSetup) override
Definition: SiStripGainFromCalibTree.cc:451
SiStripApvGain::Range
std::pair< ContainerIterator, ContainerIterator > Range
Definition: SiStripApvGain.h:27
j
int j
Definition: DBlmapReader.cc:9
SiStripGainFromCalibTree::TrigTech_token_
edm::EDGetTokenT< std::vector< bool > > TrigTech_token_
Definition: SiStripGainFromCalibTree.cc:228
SiStripGainFromCalibTree::chargeoverpath_token_
edm::EDGetTokenT< std::vector< double > > chargeoverpath_token_
Definition: SiStripGainFromCalibTree.cc:252
Exhume::I
const std::complex< double > I
Definition: I.h:8
SiStripGainFromCalibTree::gainused
const std::vector< double > * gainused
Definition: SiStripGainFromCalibTree.cc:254
SiStripGainFromCalibTree::Charge_Vs_PathlengthTECP1
std::vector< MonitorElement * > Charge_Vs_PathlengthTECP1
Definition: SiStripGainFromCalibTree.cc:206
SiStripGainFromCalibTree::localdirz_token_
edm::EDGetTokenT< std::vector< double > > localdirz_token_
Definition: SiStripGainFromCalibTree.cc:244
SiStripGainFromCalibTree::firststrip
const std::vector< unsigned short > * firststrip
Definition: SiStripGainFromCalibTree.cc:245
SiStripGainFromCalibTree::Charge_Vs_PathlengthTIDM
std::vector< MonitorElement * > Charge_Vs_PathlengthTIDM
Definition: SiStripGainFromCalibTree.cc:205
SiStripGainFromCalibTree::localdiry
const std::vector< double > * localdiry
Definition: SiStripGainFromCalibTree.cc:243
SiStripGainFromCalibTree::charge_token_
edm::EDGetTokenT< std::vector< unsigned int > > charge_token_
Definition: SiStripGainFromCalibTree.cc:250
DetId::subdetId
int subdetId() const
get the contents of the subdetector field (not cast into any detector&#39;s numbering enum) ...
Definition: DetId.h:37
funct::integral
Integral< F, X >::type integral(const F &f)
Definition: Integral.h:69
SiStripGainFromCalibTree::algoEndJob
virtual void algoEndJob() override
Definition: SiStripGainFromCalibTree.cc:755
SiStripGainFromCalibTree::nstrips_token_
edm::EDGetTokenT< std::vector< unsigned short > > nstrips_token_
Definition: SiStripGainFromCalibTree.cc:246
SiStripGainFromCalibTree::Charge_Vs_PathlengthTECP2
std::vector< MonitorElement * > Charge_Vs_PathlengthTECP2
Definition: SiStripGainFromCalibTree.cc:207
SiStripGainFromCalibTree::trackindex
const std::vector< int > * trackindex
Definition: SiStripGainFromCalibTree.cc:240
DetId
Definition: DetId.h:18
SiStripGainFromCalibTree::statCollectionFromMode
int statCollectionFromMode(const char *tag) const
Definition: SiStripGainFromCalibTree.cc:274
SiStripGainFromCalibTree::Charge_Vs_PathlengthTIB
std::vector< MonitorElement * > Charge_Vs_PathlengthTIB
Definition: SiStripGainFromCalibTree.cc:202
SiStripGainFromCalibTree::firststrip_token_
edm::EDGetTokenT< std::vector< unsigned short > > firststrip_token_
Definition: SiStripGainFromCalibTree.cc:245
SiStripGainFromCalibTree::algoAnalyze
virtual void algoAnalyze(const edm::Event &, const edm::EventSetup &) override
Definition: SiStripGainFromCalibTree.cc:1274
SiStripGainFromCalibTree::trackpt_token_
edm::EDGetTokenT< std::vector< float > > trackpt_token_
Definition: SiStripGainFromCalibTree.cc:233
edm::Handle::product
T const * product() const
Definition: Handle.h:81
SiStripGainFromCalibTree::overlapping
const std::vector< bool > * overlapping
Definition: SiStripGainFromCalibTree.cc:248
SiStripGainFromCalibTree::saturation
const std::vector< bool > * saturation
Definition: SiStripGainFromCalibTree.cc:247
fed_dqm_sourceclient-live_cfg.path
string path
Definition: fed_dqm_sourceclient-live_cfg.py:33
edm::merge
bool merge(LuminosityBlockRange &lh, LuminosityBlockRange &rh)
Definition: LuminosityBlockRange.cc:78
SiStripGainFromCalibTree::getPeakOfLandau
void getPeakOfLandau(TH1 *InputHisto, double *FitResults, double LowRange=50, double HighRange=5400)
Definition: SiStripGainFromCalibTree.cc:795
edm::EventSetup::get
const T & get() const
Definition: EventSetup.h:56
edm::ESHandle::product
T const * product() const
Definition: ESHandle.h:86
SiStripGainFromCalibTree::farfromedge
const std::vector< bool > * farfromedge
Definition: SiStripGainFromCalibTree.cc:249
SiStripGainFromCalibTree::TrigTech
const std::vector< bool > * TrigTech
Definition: SiStripGainFromCalibTree.cc:228
ConditionDBWriter< SiStripApvGain >::setDoStore
void setDoStore(const bool doStore)
When set to false the payload will not be written to the db.
Definition: ConditionDBWriter.h:371
edm::EventBase::id
edm::EventID id() const
Definition: EventBase.h:59
a
double a
Definition: hdecay.h:121
SiStripGainFromCalibTree::localdirx_token_
edm::EDGetTokenT< std::vector< double > > localdirx_token_
Definition: SiStripGainFromCalibTree.cc:242
SiStripGainFromCalibTree::saturation_token_
edm::EDGetTokenT< std::vector< bool > > saturation_token_
Definition: SiStripGainFromCalibTree.cc:247
SiStripGainFromCalibTree::trackchi2ndof_token_
edm::EDGetTokenT< std::vector< double > > trackchi2ndof_token_
Definition: SiStripGainFromCalibTree.cc:231
SiStripGainFromCalibTree::Charge_Vs_PathlengthTIDP
std::vector< MonitorElement * > Charge_Vs_PathlengthTIDP
Definition: SiStripGainFromCalibTree.cc:204
SiStripGainFromCalibTree::localdiry_token_
edm::EDGetTokenT< std::vector< double > > localdiry_token_
Definition: SiStripGainFromCalibTree.cc:243
edm::false
volatile std::atomic< bool > shutdown_flag false
Definition: UnixSignalHandlers.cc:22
SiStripGainFromCalibTree::charge
const std::vector< unsigned int > * charge
Definition: SiStripGainFromCalibTree.cc:250
edm::ESHandleBase::isValid
bool isValid() const
Definition: ESHandle.h:47
SiStripGainFromCalibTree::path
const std::vector< double > * path
Definition: SiStripGainFromCalibTree.cc:251
T
long double T
Definition: Basic3DVectorLD.h:57
SiStripGainFromCalibTree::trackphi
const std::vector< double > * trackphi
Definition: SiStripGainFromCalibTree.cc:235
SiStripGainFromCalibTree::nstrips
const std::vector< unsigned short > * nstrips
Definition: SiStripGainFromCalibTree.cc:246
SiStripGainFromCalibTree::gainused_token_
edm::EDGetTokenT< std::vector< double > > gainused_token_
Definition: SiStripGainFromCalibTree.cc:254
SiStripGainFromCalibTree::bookDQMHistos
void bookDQMHistos(const char *dqm_dir, const char *tag)
Definition: SiStripGainFromCalibTree.cc:411
SiStripGainFromCalibTree::trackphi_token_
edm::EDGetTokenT< std::vector< double > > trackphi_token_
Definition: SiStripGainFromCalibTree.cc:235
SiStripGainFromCalibTree::tracketa
const std::vector< double > * tracketa
Definition: SiStripGainFromCalibTree.cc:234
edm::Run
Definition: Run.h:43
Generated for CMSSW Reference Manual by doxygen 1.8.5