Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-01-24 01:11:12

0001 // Original Author:  Loic QUERTENMONT
0002 //         Created:  Wed Feb  6 08:55:18 CET 2008
0003 
0004 #include <memory>
0005 
0006 #include "FWCore/Framework/interface/Frameworkfwd.h"
0007 #include "FWCore/Framework/interface/Event.h"
0008 #include "FWCore/Framework/interface/MakerMacros.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include "FWCore/Utilities/interface/Exception.h"
0012 
0013 #include "Geometry/CommonDetUnit/interface/GeomDet.h"
0014 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
0015 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0016 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
0017 #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
0018 
0019 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
0020 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
0021 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
0022 #include "Geometry/CommonDetUnit/interface/PixelGeomDetUnit.h"
0023 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
0024 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
0025 
0026 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
0027 
0028 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
0029 #include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
0030 
0031 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
0032 #include "DataFormats/TrackReco/interface/Track.h"
0033 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0034 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
0035 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
0036 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
0037 #include "DataFormats/DetId/interface/DetId.h"
0038 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
0039 #include "DataFormats/TrackerCommon/interface/TrackerTopology.h"
0040 #include "Geometry/Records/interface/TrackerTopologyRcd.h"
0041 #include "DataFormats/TrackReco/interface/DeDxHit.h"
0042 #include "DataFormats/TrackReco/interface/TrackDeDxHits.h"
0043 
0044 #include "CommonTools/ConditionDBWriter/interface/ConditionDBWriter.h"
0045 #include "CondFormats/SiStripObjects/interface/SiStripApvGain.h"
0046 
0047 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0048 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
0049 
0050 #include "DQMServices/Core/interface/DQMStore.h"
0051 
0052 #include "CalibFormats/SiStripObjects/interface/SiStripGain.h"
0053 #include "CalibTracker/Records/interface/SiStripGainRcd.h"
0054 
0055 #include "TFile.h"
0056 #include "TObjString.h"
0057 #include "TString.h"
0058 #include "TH1F.h"
0059 #include "TH2F.h"
0060 #include "TProfile.h"
0061 #include "TF1.h"
0062 #include "TROOT.h"
0063 
0064 #include <unordered_map>
0065 #include <memory>
0066 
0067 using namespace edm;
0068 using namespace reco;
0069 using namespace std;
0070 
0071 struct stAPVGain {
0072   unsigned int Index;
0073   int DetId;
0074   int APVId;
0075   int SubDet;
0076   float Eta;
0077   float R;
0078   float Phi;
0079   float Thickness;
0080   double MPV;
0081   double Gain;
0082   double PreviousGain;
0083   char Side;
0084 };
0085 
0086 class SiStripGainFromData : public ConditionDBWriter<SiStripApvGain> {
0087 public:
0088   typedef dqm::legacy::MonitorElement MonitorElement;
0089   typedef dqm::legacy::DQMStore DQMStore;
0090   explicit SiStripGainFromData(const edm::ParameterSet&);
0091   ~SiStripGainFromData() override;
0092 
0093 private:
0094   void algoBeginJob(const edm::EventSetup&) override;
0095   void algoEndJob() override;
0096   void algoBeginRun(const edm::Run&, const edm::EventSetup&) override;
0097   //      virtual void algoBeginRun(const edm::Event& iEvent, const edm::EventSetup& iSetup);
0098   void algoAnalyze(const edm::Event&, const edm::EventSetup&) override;
0099 
0100   std::unique_ptr<SiStripApvGain> getNewObject() override;
0101   DQMStore* dqmStore_;
0102 
0103   double ComputeChargeOverPath(const SiStripCluster* Cluster,
0104                                TrajectoryStateOnSurface trajState,
0105                                const edm::EventSetup* iSetup,
0106                                const Track* track,
0107                                double trajChi2OverN);
0108   bool IsFarFromBorder(TrajectoryStateOnSurface trajState, const uint32_t detid, const edm::EventSetup* iSetup);
0109 
0110   void getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange = 0, double HighRange = 5400);
0111 
0112   const edm::EventSetup* iSetup_;
0113   const edm::Event* iEvent_;
0114 
0115   bool CheckLocalAngle;
0116   unsigned int MinNrEntries;
0117   double MaxMPVError;
0118   double MaxChi2OverNDF;
0119   double MinTrackMomentum;
0120   double MaxTrackMomentum;
0121   double MinTrackEta;
0122   double MaxTrackEta;
0123   unsigned int MaxNrStrips;
0124   unsigned int MinTrackHits;
0125   double MaxTrackChiOverNdf;
0126   bool AllowSaturation;
0127   bool FirstSetOfConstants;
0128   bool Validation;
0129   int CalibrationLevel;
0130   bool CheckIfFileExist;
0131 
0132   std::string AlgoMode;
0133   std::string OutputGains;
0134   std::string OutputHistos;
0135   std::string TrajToTrackProducer;
0136   std::string TrajToTrackLabel;
0137 
0138   vector<string> VInputFiles;
0139 
0140   MonitorElement* tmp;
0141 
0142   TH2F* Tracks_P_Vs_Eta;
0143   TH2F* Tracks_Pt_Vs_Eta;
0144 
0145   TH1F* NumberOfEntriesByAPV;
0146   TH1F* HChi2OverNDF;
0147   TH1F* HTrackChi2OverNDF;
0148   TH1F* HTrackHits;
0149 
0150   TH2F* MPV_Vs_EtaTIB;
0151   TH2F* MPV_Vs_EtaTID;
0152   TH2F* MPV_Vs_EtaTOB;
0153   TH2F* MPV_Vs_EtaTEC;
0154   TH2F* MPV_Vs_EtaTEC1;
0155   TH2F* MPV_Vs_EtaTEC2;
0156 
0157   TH2F* MPV_Vs_PhiTIB;
0158   TH2F* MPV_Vs_PhiTID;
0159   TH2F* MPV_Vs_PhiTOB;
0160   TH2F* MPV_Vs_PhiTEC;
0161   TH2F* MPV_Vs_PhiTEC1;
0162   TH2F* MPV_Vs_PhiTEC2;
0163 
0164   TH2F* Charge_Vs_PathTIB;
0165   TH2F* Charge_Vs_PathTID;
0166   TH2F* Charge_Vs_PathTOB;
0167   TH2F* Charge_Vs_PathTEC;
0168   TH2F* Charge_Vs_PathTEC1;
0169   TH2F* Charge_Vs_PathTEC2;
0170 
0171   TH1F* MPV_Vs_PathTIB;
0172   TH1F* MPV_Vs_PathTID;
0173   TH1F* MPV_Vs_PathTOB;
0174   TH1F* MPV_Vs_PathTEC;
0175   TH1F* MPV_Vs_PathTEC1;
0176   TH1F* MPV_Vs_PathTEC2;
0177 
0178   TH1F* Charge_TIB;
0179   TH1F* Charge_TID;
0180   TH1F* Charge_TIDP;
0181   TH1F* Charge_TIDM;
0182   TH1F* Charge_TOB;
0183   TH1F* Charge_TEC;
0184   TH1F* Charge_TEC1;
0185   TH1F* Charge_TEC2;
0186   TH1F* Charge_TECP;
0187   TH1F* Charge_TECM;
0188 
0189   TH2F* MPV_Vs_Phi;
0190   TH2F* MPV_Vs_Eta;
0191   TH2F* MPV_Vs_R;
0192 
0193   //      TH2F*        PD_Vs_Eta;
0194   //      TH2F*        PD_Vs_R;
0195 
0196   TH1F* APV_DetId;
0197   TH1F* APV_Id;
0198   TH1F* APV_Eta;
0199   TH1F* APV_R;
0200   TH1F* APV_SubDet;
0201   TH2F* APV_Momentum;
0202   TH2F* APV_Charge;
0203   TH2F* APV_PathLength;
0204   TH1F* APV_PathLengthM;
0205   TH1F* APV_MPV;
0206   TH1F* APV_Gain;
0207   TH1F* APV_CumulGain;
0208   TH1F* APV_PrevGain;
0209   TH1F* APV_Thickness;
0210 
0211   TH1F* MPVs;
0212   TH1F* MPVs320;
0213   TH1F* MPVs500;
0214 
0215   //      TH2F*        MPV_vs_10RplusEta;
0216 
0217   TH1F* NSatStripInCluster;
0218   TH1F* NHighStripInCluster;
0219   //      TH2F*        Charge_Vs_PathLength_CS1;
0220   //      TH2F*        Charge_Vs_PathLength_CS2;
0221   //      TH2F*        Charge_Vs_PathLength_CS3;
0222   //      TH2F*        Charge_Vs_PathLength_CS4;
0223   //      TH2F*        Charge_Vs_PathLength_CS5;
0224 
0225   //      TH1F*        MPV_Vs_PathLength_CS1;
0226   //      TH1F*        MPV_Vs_PathLength_CS2;
0227   //      TH1F*        MPV_Vs_PathLength_CS3;
0228   //      TH1F*        MPV_Vs_PathLength_CS4;
0229   //      TH1F*        MPV_Vs_PathLength_CS5;
0230 
0231   //      TH1F*        FWHM_Vs_PathLength_CS1;
0232   //      TH1F*        FWHM_Vs_PathLength_CS2;
0233   //      TH1F*        FWHM_Vs_PathLength_CS3;
0234   //      TH1F*        FWHM_Vs_PathLength_CS4;
0235   //      TH1F*        FWHM_Vs_PathLength_CS5;
0236 
0237   TH2F* Charge_Vs_PathLength;
0238   TH2F* Charge_Vs_PathLength320;
0239   TH2F* Charge_Vs_PathLength500;
0240 
0241   TH1F* MPV_Vs_PathLength;
0242   TH1F* MPV_Vs_PathLength320;
0243   TH1F* MPV_Vs_PathLength500;
0244 
0245   TH1F* FWHM_Vs_PathLength;
0246   TH1F* FWHM_Vs_PathLength320;
0247   TH1F* FWHM_Vs_PathLength500;
0248 
0249   TH2F* Charge_Vs_TransversAngle;
0250   TH1F* MPV_Vs_TransversAngle;
0251   TH2F* NStrips_Vs_TransversAngle;
0252 
0253   TH2F* Charge_Vs_Alpha;
0254   TH1F* MPV_Vs_Alpha;
0255   TH2F* NStrips_Vs_Alpha;
0256 
0257   TH2F* Charge_Vs_Beta;
0258   TH1F* MPV_Vs_Beta;
0259   TH2F* NStrips_Vs_Beta;
0260 
0261   TH2F* Error_Vs_MPV;
0262   TH2F* Error_Vs_Entries;
0263   TH2F* Error_Vs_Eta;
0264   TH2F* Error_Vs_Phi;
0265 
0266   TH2F* NoMPV_Vs_EtaPhi;
0267 
0268   TH2F* HitLocalPosition;
0269   TH2F* HitLocalPositionBefCut;
0270 
0271   TH1F* JobInfo;
0272 
0273   TH1F* HFirstStrip;
0274 
0275   unsigned int NEvent;
0276   unsigned int SRun;
0277   unsigned int SEvent;
0278   TimeValue_t STimestamp;
0279   unsigned int ERun;
0280   unsigned int EEvent;
0281   TimeValue_t ETimestamp;
0282 
0283 private:
0284   class isEqual {
0285   public:
0286     template <class T>
0287     bool operator()(const T& PseudoDetId1, const T& PseudoDetId2) {
0288       return PseudoDetId1 == PseudoDetId2;
0289     }
0290   };
0291 
0292   std::vector<stAPVGain*> APVsCollOrdered;
0293   std::unordered_map<unsigned int, stAPVGain*> APVsColl;
0294 
0295   edm::ESGetToken<TrackerTopology, TrackerTopologyRcd> tTopoToken_;
0296   edm::ESGetToken<TrackerGeometry, TrackerDigiGeometryRecord> tkGeomToken_;
0297   edm::ESGetToken<SiStripGain, SiStripGainRcd> gainToken_;
0298 };
0299 
0300 SiStripGainFromData::SiStripGainFromData(const edm::ParameterSet& iConfig)
0301     : ConditionDBWriter<SiStripApvGain>(iConfig) {
0302   AlgoMode = iConfig.getParameter<std::string>("AlgoMode");
0303 
0304   OutputGains = iConfig.getParameter<std::string>("OutputGains");
0305   OutputHistos = iConfig.getParameter<std::string>("OutputHistos");
0306 
0307   TrajToTrackProducer = iConfig.getParameter<std::string>("TrajToTrackProducer");
0308   TrajToTrackLabel = iConfig.getParameter<std::string>("TrajToTrackLabel");
0309 
0310   CheckLocalAngle = iConfig.getUntrackedParameter<bool>("checkLocalAngle", false);
0311   MinNrEntries = iConfig.getUntrackedParameter<unsigned>("minNrEntries", 20);
0312   MaxMPVError = iConfig.getUntrackedParameter<double>("maxMPVError", 500.0);
0313   MaxChi2OverNDF = iConfig.getUntrackedParameter<double>("maxChi2OverNDF", 5.0);
0314   MinTrackMomentum = iConfig.getUntrackedParameter<double>("minTrackMomentum", 3.0);
0315   MaxTrackMomentum = iConfig.getUntrackedParameter<double>("maxTrackMomentum", 99999.0);
0316   MinTrackEta = iConfig.getUntrackedParameter<double>("minTrackEta", -5.0);
0317   MaxTrackEta = iConfig.getUntrackedParameter<double>("maxTrackEta", 5.0);
0318   MaxNrStrips = iConfig.getUntrackedParameter<unsigned>("maxNrStrips", 2);
0319   MinTrackHits = iConfig.getUntrackedParameter<unsigned>("MinTrackHits", 8);
0320   MaxTrackChiOverNdf = iConfig.getUntrackedParameter<double>("MaxTrackChiOverNdf", 3);
0321   AllowSaturation = iConfig.getUntrackedParameter<bool>("AllowSaturation", false);
0322   FirstSetOfConstants = iConfig.getUntrackedParameter<bool>("FirstSetOfConstants", true);
0323   Validation = iConfig.getUntrackedParameter<bool>("Validation", false);
0324   CheckIfFileExist = iConfig.getUntrackedParameter<bool>("CheckIfFileExist", false);
0325 
0326   CalibrationLevel = iConfig.getUntrackedParameter<int>("CalibrationLevel", 0);
0327 
0328   if (strcmp(AlgoMode.c_str(), "WriteOnDB") == 0)
0329     VInputFiles = iConfig.getParameter<vector<string> >("VInputFiles");
0330 
0331   dqmStore_ = edm::Service<DQMStore>().operator->();
0332 
0333   tTopoToken_ = esConsumes<edm::Transition::BeginRun>();
0334   tkGeomToken_ = esConsumes<edm::Transition::BeginRun>();
0335   gainToken_ = esConsumes<edm::Transition::BeginRun>();
0336 
0337   //if( OutputHistos!="" )
0338   //  dqmStore_->open(OutputHistos.c_str(), true);
0339 }
0340 
0341 SiStripGainFromData::~SiStripGainFromData() {}
0342 
0343 void SiStripGainFromData::algoBeginJob(const edm::EventSetup& iSetup) {
0344   const TrackerTopology* const tTopo = &iSetup.getData(tTopoToken_);
0345 
0346   iSetup_ = &iSetup;
0347 
0348   //   TH1::AddDirectory(kTRUE);
0349 
0350   tmp = dqmStore_->book1D("JobInfo", "JobInfo", 20, 0, 20);
0351   JobInfo = tmp->getTH1F();
0352 
0353   tmp = dqmStore_->book1D("APV_DetId", "APV_DetId", 72785, 0, 72784);
0354   APV_DetId = tmp->getTH1F();
0355   tmp = dqmStore_->book1D("APV_Id", "APV_Id", 72785, 0, 72784);
0356   APV_Id = tmp->getTH1F();
0357   tmp = dqmStore_->book1D("APV_Eta", "APV_Eta", 72785, 0, 72784);
0358   APV_Eta = tmp->getTH1F();
0359   tmp = dqmStore_->book1D("APV_R", "APV_R", 72785, 0, 72784);
0360   APV_R = tmp->getTH1F();
0361   tmp = dqmStore_->book1D("APV_SubDet", "APV_SubDet", 72785, 0, 72784);
0362   APV_SubDet = tmp->getTH1F();
0363   tmp = dqmStore_->book2D("APV_Momentum", "APV_Momentum", 72785, 0, 72784, 50, 0, 200);
0364   APV_Momentum = tmp->getTH2F();
0365   tmp = dqmStore_->book2D("APV_Charge", "APV_Charge", 72785, 0, 72784, 1000, 0, 2000);
0366   APV_Charge = tmp->getTH2F();
0367   tmp = dqmStore_->book2D("APV_PathLength", "APV_PathLength", 72785, 0, 72784, 100, 0.2, 1.4);
0368   APV_PathLength = tmp->getTH2F();
0369   tmp = dqmStore_->book1D("APV_PathLengthM", "APV_PathLengthM", 72785, 0, 72784);
0370   APV_PathLengthM = tmp->getTH1F();
0371   tmp = dqmStore_->book1D("APV_MPV", "APV_MPV", 72785, 0, 72784);
0372   APV_MPV = tmp->getTH1F();
0373   tmp = dqmStore_->book1D("APV_Gain", "APV_Gain", 72785, 0, 72784);
0374   APV_Gain = tmp->getTH1F();
0375   tmp = dqmStore_->book1D("APV_PrevGain", "APV_PrevGain", 72785, 0, 72784);
0376   APV_PrevGain = tmp->getTH1F();
0377   tmp = dqmStore_->book1D("APV_CumulGain", "APV_CumulGain", 72785, 0, 72784);
0378   APV_CumulGain = tmp->getTH1F();
0379   tmp = dqmStore_->book1D("APV_Thickness", "APV_Thicknes", 72785, 0, 72784);
0380   APV_Thickness = tmp->getTH1F();
0381 
0382   tmp = dqmStore_->book2D("Tracks_P_Vs_Eta", "Tracks_P_Vs_Eta", 30, 0, 3, 50, 0, 200);
0383   Tracks_P_Vs_Eta = tmp->getTH2F();
0384   tmp = dqmStore_->book2D("Tracks_Pt_Vs_Eta", "Tracks_Pt_Vs_Eta", 30, 0, 3, 50, 0, 200);
0385   Tracks_Pt_Vs_Eta = tmp->getTH2F();
0386 
0387   tmp = dqmStore_->book2D("Charge_Vs_PathTIB", "Charge_Vs_PathTIB", 100, 0.2, 1.4, 500, 0, 2000);
0388   Charge_Vs_PathTIB = tmp->getTH2F();
0389   tmp = dqmStore_->book2D("Charge_Vs_PathTID", "Charge_Vs_PathTID", 100, 0.2, 1.4, 500, 0, 2000);
0390   Charge_Vs_PathTID = tmp->getTH2F();
0391   tmp = dqmStore_->book2D("Charge_Vs_PathTOB", "Charge_Vs_PathTOB", 100, 0.2, 1.4, 500, 0, 2000);
0392   Charge_Vs_PathTOB = tmp->getTH2F();
0393   tmp = dqmStore_->book2D("Charge_Vs_PathTEC", "Charge_Vs_PathTEC", 100, 0.2, 1.4, 500, 0, 2000);
0394   Charge_Vs_PathTEC = tmp->getTH2F();
0395   tmp = dqmStore_->book2D("Charge_Vs_PathTEC1", "Charge_Vs_PathTEC1", 100, 0.2, 1.4, 500, 0, 2000);
0396   Charge_Vs_PathTEC1 = tmp->getTH2F();
0397   tmp = dqmStore_->book2D("Charge_Vs_PathTEC2", "Charge_Vs_PathTEC2", 100, 0.2, 1.4, 500, 0, 2000);
0398   Charge_Vs_PathTEC2 = tmp->getTH2F();
0399 
0400   tmp = dqmStore_->book1D("Charge_TIB", "Charge_TIB", 1000, 0, 2000);
0401   Charge_TIB = tmp->getTH1F();
0402   tmp = dqmStore_->book1D("Charge_TID", "Charge_TID", 1000, 0, 2000);
0403   Charge_TID = tmp->getTH1F();
0404   tmp = dqmStore_->book1D("Charge_TID+", "Charge_TID+", 1000, 0, 2000);
0405   Charge_TIDP = tmp->getTH1F();
0406   tmp = dqmStore_->book1D("Charge_TID-", "Charge_TID-", 1000, 0, 2000);
0407   Charge_TIDM = tmp->getTH1F();
0408   tmp = dqmStore_->book1D("Charge_TOB", "Charge_TOB", 1000, 0, 2000);
0409   Charge_TOB = tmp->getTH1F();
0410   tmp = dqmStore_->book1D("Charge_TEC", "Charge_TEC", 1000, 0, 2000);
0411   Charge_TEC = tmp->getTH1F();
0412   tmp = dqmStore_->book1D("Charge_TEC1", "Charge_TEC1", 1000, 0, 2000);
0413   Charge_TEC1 = tmp->getTH1F();
0414   tmp = dqmStore_->book1D("Charge_TEC2", "Charge_TEC2", 1000, 0, 2000);
0415   Charge_TEC2 = tmp->getTH1F();
0416   tmp = dqmStore_->book1D("Charge_TEC+", "Charge_TEC+", 1000, 0, 2000);
0417   Charge_TECP = tmp->getTH1F();
0418   tmp = dqmStore_->book1D("Charge_TEC-", "Charge_TEC-", 1000, 0, 2000);
0419   Charge_TECM = tmp->getTH1F();
0420 
0421   /*
0422    tmp  = dqmStore_->book2D ("Charge_Vs_PathLength_CS1", "Charge_Vs_PathLength_CS1"  , 250,0.2,1.4, 500,0,2000); Charge_Vs_PathLength_CS1 = tmp->getTH2F();
0423    tmp  = dqmStore_->book2D ("Charge_Vs_PathLength_CS2", "Charge_Vs_PathLength_CS2"  , 250,0.2,1.4, 500,0,2000); Charge_Vs_PathLength_CS2 = tmp->getTH2F();
0424    tmp  = dqmStore_->book2D ("Charge_Vs_PathLength_CS3", "Charge_Vs_PathLength_CS3"  , 250,0.2,1.4, 500,0,2000); Charge_Vs_PathLength_CS3 = tmp->getTH2F();
0425    tmp  = dqmStore_->book2D ("Charge_Vs_PathLength_CS4", "Charge_Vs_PathLength_CS4"  , 250,0.2,1.4, 500,0,2000); Charge_Vs_PathLength_CS4 = tmp->getTH2F();
0426    tmp  = dqmStore_->book2D ("Charge_Vs_PathLength_CS5", "Charge_Vs_PathLength_CS5"  , 250,0.2,1.4, 500,0,2000); Charge_Vs_PathLength_CS5 = tmp->getTH2F();
0427 */
0428   tmp = dqmStore_->book2D("Charge_Vs_PathLength", "Charge_Vs_PathLength", 100, 0.2, 1.4, 500, 0, 2000);
0429   Charge_Vs_PathLength = tmp->getTH2F();
0430   tmp = dqmStore_->book2D("Charge_Vs_PathLength320", "Charge_Vs_PathLength", 100, 0.2, 1.4, 500, 0, 2000);
0431   Charge_Vs_PathLength320 = tmp->getTH2F();
0432   tmp = dqmStore_->book2D("Charge_Vs_PathLength500", "Charge_Vs_PathLength", 100, 0.2, 1.4, 500, 0, 2000);
0433   Charge_Vs_PathLength500 = tmp->getTH2F();
0434 
0435   tmp = dqmStore_->book2D("Charge_Vs_TransversAngle", "Charge_Vs_TransversAngle", 220, -20, 200, 500, 0, 2000);
0436   Charge_Vs_TransversAngle = tmp->getTH2F();
0437   tmp = dqmStore_->book2D("Charge_Vs_Alpha", "Charge_Vs_Alpha", 220, -20, 200, 500, 0, 2000);
0438   Charge_Vs_Alpha = tmp->getTH2F();
0439   tmp = dqmStore_->book2D("Charge_Vs_Beta", "Charge_Vs_Beta", 220, -20, 200, 500, 0, 2000);
0440   Charge_Vs_Beta = tmp->getTH2F();
0441 
0442   tmp = dqmStore_->book2D("NStrips_Vs_TransversAngle", "NStrips_Vs_TransversAngle", 220, -20, 200, 10, 0, 10);
0443   NStrips_Vs_TransversAngle = tmp->getTH2F();
0444   tmp = dqmStore_->book2D("NStrips_Vs_Alpha", "NStrips_Vs_Alpha", 220, -20, 200, 10, 0, 10);
0445   NStrips_Vs_Alpha = tmp->getTH2F();
0446   tmp = dqmStore_->book2D("NStrips_Vs_Beta", "NStrips_Vs_Beta", 220, -20, 200, 10, 0, 10);
0447   NStrips_Vs_Beta = tmp->getTH2F();
0448   tmp = dqmStore_->book1D("NHighStripInCluster", "NHighStripInCluster", 15, 0, 14);
0449   NHighStripInCluster = tmp->getTH1F();
0450   tmp = dqmStore_->book1D("NSatStripInCluster", "NSatStripInCluster", 50, 0, 50);
0451   NSatStripInCluster = tmp->getTH1F();
0452 
0453   tmp = dqmStore_->book1D("TrackChi2OverNDF", "TrackChi2OverNDF", 500, 0, 10);
0454   HTrackChi2OverNDF = tmp->getTH1F();
0455   tmp = dqmStore_->book1D("TrackHits", "TrackHits", 40, 0, 40);
0456   HTrackHits = tmp->getTH1F();
0457 
0458   tmp = dqmStore_->book1D("FirstStrip", "FirstStrip", 800, 0, 800);
0459   HFirstStrip = tmp->getTH1F();
0460 
0461   if (strcmp(AlgoMode.c_str(), "MultiJob") != 0) {
0462     tmp = dqmStore_->book2D("MPV_Vs_EtaTIB", "MPV_Vs_EtaTIB", 50, -3.0, 3.0, 300, 0, 600);
0463     MPV_Vs_EtaTIB = tmp->getTH2F();
0464     tmp = dqmStore_->book2D("MPV_Vs_EtaTID", "MPV_Vs_EtaTID", 50, -3.0, 3.0, 300, 0, 600);
0465     MPV_Vs_EtaTID = tmp->getTH2F();
0466     tmp = dqmStore_->book2D("MPV_Vs_EtaTOB", "MPV_Vs_EtaTOB", 50, -3.0, 3.0, 300, 0, 600);
0467     MPV_Vs_EtaTOB = tmp->getTH2F();
0468     tmp = dqmStore_->book2D("MPV_Vs_EtaTEC", "MPV_Vs_EtaTEC", 50, -3.0, 3.0, 300, 0, 600);
0469     MPV_Vs_EtaTEC = tmp->getTH2F();
0470     tmp = dqmStore_->book2D("MPV_Vs_EtaTEC1", "MPV_Vs_EtaTEC1", 50, -3.0, 3.0, 300, 0, 600);
0471     MPV_Vs_EtaTEC1 = tmp->getTH2F();
0472     tmp = dqmStore_->book2D("MPV_Vs_EtaTEC2", "MPV_Vs_EtaTEC2", 50, -3.0, 3.0, 300, 0, 600);
0473     MPV_Vs_EtaTEC2 = tmp->getTH2F();
0474 
0475     tmp = dqmStore_->book2D("MPV_Vs_PhiTIB", "MPV_Vs_PhiTIB", 50, -3.2, 3.2, 300, 0, 600);
0476     MPV_Vs_PhiTIB = tmp->getTH2F();
0477     tmp = dqmStore_->book2D("MPV_Vs_PhiTID", "MPV_Vs_PhiTID", 50, -3.2, 3.2, 300, 0, 600);
0478     MPV_Vs_PhiTID = tmp->getTH2F();
0479     tmp = dqmStore_->book2D("MPV_Vs_PhiTOB", "MPV_Vs_PhiTOB", 50, -3.2, 3.2, 300, 0, 600);
0480     MPV_Vs_PhiTOB = tmp->getTH2F();
0481     tmp = dqmStore_->book2D("MPV_Vs_PhiTEC", "MPV_Vs_PhiTEC", 50, -3.2, 3.2, 300, 0, 600);
0482     MPV_Vs_PhiTEC = tmp->getTH2F();
0483     tmp = dqmStore_->book2D("MPV_Vs_PhiTEC1", "MPV_Vs_PhiTEC1", 50, -3.2, 3.2, 300, 0, 600);
0484     MPV_Vs_PhiTEC1 = tmp->getTH2F();
0485     tmp = dqmStore_->book2D("MPV_Vs_PhiTEC2", "MPV_Vs_PhiTEC2", 50, -3.2, 3.2, 300, 0, 600);
0486     MPV_Vs_PhiTEC2 = tmp->getTH2F();
0487 
0488     tmp = dqmStore_->book1D("MPV_Vs_PathTIB", "MPV_Vs_PathTIB", 100, 0.2, 1.4);
0489     MPV_Vs_PathTIB = tmp->getTH1F();
0490     tmp = dqmStore_->book1D("MPV_Vs_PathTID", "MPV_Vs_PathTID", 100, 0.2, 1.4);
0491     MPV_Vs_PathTID = tmp->getTH1F();
0492     tmp = dqmStore_->book1D("MPV_Vs_PathTOB", "MPV_Vs_PathTOB", 100, 0.2, 1.4);
0493     MPV_Vs_PathTOB = tmp->getTH1F();
0494     tmp = dqmStore_->book1D("MPV_Vs_PathTEC", "MPV_Vs_PathTEC", 100, 0.2, 1.4);
0495     MPV_Vs_PathTEC = tmp->getTH1F();
0496     tmp = dqmStore_->book1D("MPV_Vs_PathTEC1", "MPV_Vs_PathTEC1", 100, 0.2, 1.4);
0497     MPV_Vs_PathTEC1 = tmp->getTH1F();
0498     tmp = dqmStore_->book1D("MPV_Vs_PathTEC2", "MPV_Vs_PathTEC2", 100, 0.2, 1.4);
0499     MPV_Vs_PathTEC2 = tmp->getTH1F();
0500 
0501     tmp = dqmStore_->book2D("MPV_Vs_Phi", "MPV_Vs_Phi", 50, -3.2, 3.2, 300, 0, 600);
0502     MPV_Vs_Phi = tmp->getTH2F();
0503     tmp = dqmStore_->book2D("MPV_Vs_Eta", "MPV_Vs_Eta", 50, -3.0, 3.0, 300, 0, 600);
0504     MPV_Vs_Eta = tmp->getTH2F();
0505     tmp = dqmStore_->book2D("MPV_Vs_R", "MPV_Vs_R", 150, 0.0, 150.0, 300, 0, 600);
0506     MPV_Vs_R = tmp->getTH2F();
0507     /*   
0508       tmp  = dqmStore_->book1D ("MPV_Vs_PathLength_CS1"   , "MPV_Vs_PathLength_CS1" , 250, 0.2, 1.4); MPV_Vs_PathLength_CS1 = tmp->getTH1F();
0509       tmp  = dqmStore_->book1D ("MPV_Vs_PathLength_CS2"   , "MPV_Vs_PathLength_CS2" , 250, 0.2, 1.4); MPV_Vs_PathLength_CS2 = tmp->getTH1F();
0510       tmp  = dqmStore_->book1D ("MPV_Vs_PathLength_CS3"   , "MPV_Vs_PathLength_CS3" , 250, 0.2, 1.4); MPV_Vs_PathLength_CS3 = tmp->getTH1F();
0511       tmp  = dqmStore_->book1D ("MPV_Vs_PathLength_CS4"   , "MPV_Vs_PathLength_CS4" , 250, 0.2, 1.4); MPV_Vs_PathLength_CS4 = tmp->getTH1F();
0512       tmp  = dqmStore_->book1D ("MPV_Vs_PathLength_CS5"   , "MPV_Vs_PathLength_CS5" , 250, 0.2, 1.4); MPV_Vs_PathLength_CS5 = tmp->getTH1F();
0513 
0514       tmp  = dqmStore_->book1D ("FWHM_Vs_PathLength_CS1"  , "FWHM_Vs_PathLength_CS1", 250, 0.2, 1.4); FWHM_Vs_PathLength_CS1 = tmp->getTH1F();
0515       tmp  = dqmStore_->book1D ("FWHM_Vs_PathLength_CS2"  , "FWHM_Vs_PathLength_CS2", 250, 0.2, 1.4); FWHM_Vs_PathLength_CS2 = tmp->getTH1F();
0516       tmp  = dqmStore_->book1D ("FWHM_Vs_PathLength_CS3"  , "FWHM_Vs_PathLength_CS3", 250, 0.2, 1.4); FWHM_Vs_PathLength_CS3 = tmp->getTH1F();
0517       tmp  = dqmStore_->book1D ("FWHM_Vs_PathLength_CS4"  , "FWHM_Vs_PathLength_CS4", 250, 0.2, 1.4); FWHM_Vs_PathLength_CS4 = tmp->getTH1F();
0518       tmp  = dqmStore_->book1D ("FWHM_Vs_PathLength_CS5"  , "FWHM_Vs_PathLength_CS5", 250, 0.2, 1.4); FWHM_Vs_PathLength_CS5 = tmp->getTH1F();
0519 */
0520     tmp = dqmStore_->book1D("MPV_Vs_PathLength", "MPV_Vs_PathLength", 100, 0.2, 1.4);
0521     MPV_Vs_PathLength = tmp->getTH1F();
0522     tmp = dqmStore_->book1D("MPV_Vs_PathLength320", "MPV_Vs_PathLength", 100, 0.2, 1.4);
0523     MPV_Vs_PathLength320 = tmp->getTH1F();
0524     tmp = dqmStore_->book1D("MPV_Vs_PathLength500", "MPV_Vs_PathLength", 100, 0.2, 1.4);
0525     MPV_Vs_PathLength500 = tmp->getTH1F();
0526 
0527     tmp = dqmStore_->book1D("FWHM_Vs_PathLength", "FWHM_Vs_PathLength", 100, 0.2, 1.4);
0528     FWHM_Vs_PathLength = tmp->getTH1F();
0529     tmp = dqmStore_->book1D("FWHM_Vs_PathLength320", "FWHM_Vs_PathLength", 100, 0.2, 1.4);
0530     FWHM_Vs_PathLength320 = tmp->getTH1F();
0531     tmp = dqmStore_->book1D("FWHM_Vs_PathLength500", "FWHM_Vs_PathLength", 100, 0.2, 1.4);
0532     FWHM_Vs_PathLength500 = tmp->getTH1F();
0533 
0534     tmp = dqmStore_->book1D("MPV_Vs_TransversAngle", "MPV_Vs_TransversAngle", 220, -20, 200);
0535     MPV_Vs_TransversAngle = tmp->getTH1F();
0536     tmp = dqmStore_->book1D("MPV_Vs_Alpha", "MPV_Vs_Alpha", 220, -20, 200);
0537     MPV_Vs_Alpha = tmp->getTH1F();
0538     tmp = dqmStore_->book1D("MPV_Vs_Beta", "MPV_Vs_Beta", 220, -20, 200);
0539     MPV_Vs_Beta = tmp->getTH1F();
0540 
0541     tmp = dqmStore_->book2D("Error_Vs_MPV", "Error_Vs_MPV", 600, 0, 600, 100, 0, 50);
0542     Error_Vs_MPV = tmp->getTH2F();
0543     tmp = dqmStore_->book2D("Error_Vs_Entries", "Error_Vs_Entries", 500, 0, 10000, 100, 0, 50);
0544     Error_Vs_Entries = tmp->getTH2F();
0545     tmp = dqmStore_->book2D("Error_Vs_Eta", "Error_Vs_Eta", 50, -3.0, 3.0, 100, 0, 50);
0546     Error_Vs_Eta = tmp->getTH2F();
0547     tmp = dqmStore_->book2D("Error_Vs_Phi", "Error_Vs_Phi", 50, -3.2, 3.2, 100, 0, 50);
0548     Error_Vs_Phi = tmp->getTH2F();
0549 
0550     tmp = dqmStore_->book2D("NoMPV_Vs_EtaPhi", "NoMPV_Vs_EtaPhi", 50, -3.0, 3.0, 50, -3.2, 3.2);
0551     NoMPV_Vs_EtaPhi = tmp->getTH2F();
0552 
0553     tmp = dqmStore_->book1D("NumberOfEntriesByAPV", "NumberOfEntriesByAPV", 1000, 0, 10000);
0554     NumberOfEntriesByAPV = tmp->getTH1F();
0555     tmp = dqmStore_->book1D("Chi2OverNDF", "Chi2OverNDF", 500, 0, 25);
0556     HChi2OverNDF = tmp->getTH1F();
0557 
0558     tmp = dqmStore_->book1D("MPVs", "MPVs", 600, 0, 600);
0559     MPVs = tmp->getTH1F();
0560     tmp = dqmStore_->book1D("MPVs320", "MPVs320", 600, 0, 600);
0561     MPVs320 = tmp->getTH1F();
0562     tmp = dqmStore_->book1D("MPVs500", "MPVs500", 600, 0, 600);
0563     MPVs500 = tmp->getTH1F();
0564 
0565     //      MPV_vs_10RplusEta          tmp  = dqmStore_->book2D ("MPV_vs_10RplusEta","MPV_vs_10RplusEta", 48000,0,2400, 800,100,500);
0566   }
0567 
0568   gROOT->cd();
0569 
0570   auto const& Det = iSetup.getData(tkGeomToken_).dets();
0571 
0572   //   if(strcmp(AlgoMode.c_str(),"MultiJob")!=0 && !FirstSetOfConstants){
0573   if (!iSetup.getHandle(gainToken_)) {
0574     printf("\n#####################\n\nERROR --> gainHandle is not valid\n\n#####################\n\n");
0575     exit(0);
0576   }
0577   //   }
0578 
0579   unsigned int Id = 0;
0580   for (unsigned int i = 0; i < Det.size(); i++) {
0581     DetId Detid = Det[i]->geographicalId();
0582     int SubDet = Detid.subdetId();
0583 
0584     if (SubDet == StripSubdetector::TIB || SubDet == StripSubdetector::TID || SubDet == StripSubdetector::TOB ||
0585         SubDet == StripSubdetector::TEC) {
0586       auto DetUnit = dynamic_cast<StripGeomDetUnit const*>(Det[i]);
0587       if (!DetUnit)
0588         continue;
0589 
0590       const StripTopology& Topo = DetUnit->specificTopology();
0591       unsigned int NAPV = Topo.nstrips() / 128;
0592 
0593       double Phi = DetUnit->position().basicVector().phi();
0594       double Eta = DetUnit->position().basicVector().eta();
0595       double R = DetUnit->position().basicVector().transverse();
0596       double Thick = DetUnit->surface().bounds().thickness();
0597 
0598       for (unsigned int j = 0; j < NAPV; j++) {
0599         stAPVGain* APV = new stAPVGain;
0600         APV->Index = Id;
0601         APV->DetId = Detid.rawId();
0602         APV->APVId = j;
0603         APV->SubDet = SubDet;
0604         APV->MPV = -1;
0605         APV->Gain = -1;
0606         APV->PreviousGain = 1;
0607         APV->Eta = Eta;
0608         APV->Phi = Phi;
0609         APV->R = R;
0610         APV->Thickness = Thick;
0611         APV->Side = 0;
0612 
0613         if (SubDet == StripSubdetector::TID) {
0614           APV->Side = tTopo->tecSide(Detid);
0615         } else if (SubDet == StripSubdetector::TEC) {
0616           APV->Side = tTopo->tecSide(Detid);
0617         }
0618 
0619         APVsCollOrdered.push_back(APV);
0620         APVsColl[(APV->DetId << 3) | APV->APVId] = APV;
0621         Id++;
0622 
0623         APV_DetId->Fill(Id, APV->DetId);
0624         APV_Id->Fill(Id, APV->APVId);
0625         APV_Eta->Fill(Id, APV->Eta);
0626         APV_R->Fill(Id, APV->R);
0627         APV_SubDet->Fill(Id, APV->SubDet);
0628         APV_Thickness->Fill(Id, APV->Thickness);
0629       }
0630     }
0631   }
0632 
0633   NEvent = 0;
0634   SRun = 0;
0635   SEvent = 0;
0636   STimestamp = 0;
0637   ERun = 0;
0638   EEvent = 0;
0639   ETimestamp = 0;
0640 }
0641 
0642 void SiStripGainFromData::algoEndJob() {
0643   unsigned int I = 0;
0644 
0645   if (strcmp(AlgoMode.c_str(), "WriteOnDB") == 0 || strcmp(AlgoMode.c_str(), "Merge") == 0) {
0646     TH1::AddDirectory(kTRUE);
0647 
0648     TFile* file = nullptr;
0649     for (unsigned int f = 0; f < VInputFiles.size(); f++) {
0650       printf("Loading New Input File : %s\n", VInputFiles[f].c_str());
0651       if (CheckIfFileExist) {
0652         FILE* doesFileExist = fopen(VInputFiles[f].c_str(), "r");
0653         if (!doesFileExist) {
0654           printf("File %s doesn't exist\n", VInputFiles[f].c_str());
0655           continue;
0656         } else {
0657           fclose(doesFileExist);
0658         }
0659       }
0660       file = new TFile(VInputFiles[f].c_str());
0661       if (!file || file->IsZombie()) {
0662         printf("### Bug With File %s\n### File will be skipped \n", VInputFiles[f].c_str());
0663         continue;
0664       }
0665       APV_Charge->Add((TH1*)file->FindObjectAny("APV_Charge"), 1);
0666       APV_Momentum->Add((TH1*)file->FindObjectAny("APV_Momentum"), 1);
0667       APV_PathLength->Add((TH1*)file->FindObjectAny("APV_PathLength"), 1);
0668 
0669       Tracks_P_Vs_Eta->Add((TH1*)file->FindObjectAny("Tracks_P_Vs_Eta"), 1);
0670       Tracks_Pt_Vs_Eta->Add((TH1*)file->FindObjectAny("Tracks_Pt_Vs_Eta"), 1);
0671 
0672       Charge_Vs_PathTIB->Add((TH1*)file->FindObjectAny("Charge_Vs_PathTIB"), 1);
0673       Charge_Vs_PathTID->Add((TH1*)file->FindObjectAny("Charge_Vs_PathTID"), 1);
0674       Charge_Vs_PathTOB->Add((TH1*)file->FindObjectAny("Charge_Vs_PathTOB"), 1);
0675       Charge_Vs_PathTEC->Add((TH1*)file->FindObjectAny("Charge_Vs_PathTEC"), 1);
0676       Charge_Vs_PathTEC1->Add((TH1*)file->FindObjectAny("Charge_Vs_PathTEC1"), 1);
0677       Charge_Vs_PathTEC2->Add((TH1*)file->FindObjectAny("Charge_Vs_PathTEC2"), 1);
0678 
0679       HTrackChi2OverNDF->Add((TH1*)file->FindObjectAny("TrackChi2OverNDF"), 1);
0680       HTrackHits->Add((TH1*)file->FindObjectAny("TrackHits"), 1);
0681 
0682       NHighStripInCluster->Add((TH1*)file->FindObjectAny("NHighStripInCluster"), 1);
0683       NSatStripInCluster->Add((TH1*)file->FindObjectAny("NSatStripInCluster"), 1);
0684       Charge_Vs_PathLength->Add((TH1*)file->FindObjectAny("Charge_Vs_PathLength"), 1);
0685       Charge_Vs_PathLength320->Add((TH1*)file->FindObjectAny("Charge_Vs_PathLength320"), 1);
0686       Charge_Vs_PathLength500->Add((TH1*)file->FindObjectAny("Charge_Vs_PathLength500"), 1);
0687       Charge_Vs_TransversAngle->Add((TH1*)file->FindObjectAny("Charge_Vs_TransversAngle"), 1);
0688       NStrips_Vs_TransversAngle->Add((TH1*)file->FindObjectAny("NStrips_Vs_TransversAngle"), 1);
0689       Charge_Vs_Alpha->Add((TH1*)file->FindObjectAny("Charge_Vs_Alpha"), 1);
0690       NStrips_Vs_Alpha->Add((TH1*)file->FindObjectAny("NStrips_Vs_Alpha"), 1);
0691       HFirstStrip->Add((TH1*)file->FindObjectAny("FirstStrip"), 1);
0692 
0693       TH1F* JobInfo_tmp = (TH1F*)file->FindObjectAny("JobInfo");
0694       NEvent += (unsigned int)JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(1));
0695       unsigned int tmp_SRun = (unsigned int)JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(3));
0696       unsigned int tmp_SEvent = (unsigned int)JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(4));
0697       unsigned int tmp_ERun = (unsigned int)JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(6));
0698       unsigned int tmp_EEvent = (unsigned int)JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(7));
0699 
0700       if (SRun == 0)
0701         SRun = tmp_SRun;
0702 
0703       if (tmp_SRun < SRun) {
0704         SRun = tmp_SRun;
0705         SEvent = tmp_SEvent;
0706       } else if (tmp_SRun == SRun && tmp_SEvent < SEvent) {
0707         SEvent = tmp_SEvent;
0708       }
0709 
0710       if (tmp_ERun > ERun) {
0711         ERun = tmp_ERun;
0712         EEvent = tmp_EEvent;
0713       } else if (tmp_ERun == ERun && tmp_EEvent > EEvent) {
0714         EEvent = tmp_EEvent;
0715       }
0716 
0717       printf("Deleting Current Input File\n");
0718       file->Close();
0719       delete file;
0720     }
0721   }
0722 
0723   JobInfo->Fill(1, NEvent);
0724   JobInfo->Fill(3, SRun);
0725   JobInfo->Fill(4, SEvent);
0726   JobInfo->Fill(6, ERun);
0727   JobInfo->Fill(7, EEvent);
0728 
0729   if (strcmp(AlgoMode.c_str(), "MultiJob") != 0) {
0730     TH1D* Proj = nullptr;
0731     double FitResults[5];
0732     I = 0;
0733     for (auto it = APVsColl.begin(); it != APVsColl.end(); it++) {
0734       if (I % 3650 == 0)
0735         printf("Fitting Histograms \t %6.2f%%\n", (100.0 * I) / APVsColl.size());
0736       I++;
0737       stAPVGain* APV = it->second;
0738 
0739       int bin = APV_Charge->GetXaxis()->FindBin(APV->Index);
0740       Proj = APV_Charge->ProjectionY(" ", bin, bin, "e");
0741       Proj = (TH1D*)Proj->Clone();
0742       if (Proj == nullptr)
0743         continue;
0744 
0745       // ADD PROJECTTIONS COMMING FROM THE SECOND APV IN THE PAIR
0746       if (CalibrationLevel == 1) {
0747         int SecondAPVId = APV->APVId;
0748         if (SecondAPVId % 2 == 0) {
0749           SecondAPVId = SecondAPVId + 1;
0750         } else {
0751           SecondAPVId = SecondAPVId - 1;
0752         }
0753         stAPVGain* APV2 = APVsColl[(APV->DetId << 3) | SecondAPVId];
0754 
0755         int bin2 = APV_Charge->GetXaxis()->FindBin(APV2->Index);
0756         TH1D* Proj2 = APV_Charge->ProjectionY(" ", bin2, bin2, "e");
0757         if (Proj2 != nullptr) {
0758           Proj->Add(Proj2, 1);
0759         }
0760       } else if (CalibrationLevel > 1) {
0761         //       printf("%8i %i--> %4.0f + %4.0f\n",APV->DetId, APV->APVId, 0.0, Proj->GetEntries());
0762         for (auto it2 = APVsColl.begin(); it2 != APVsColl.end(); it2++) {
0763           stAPVGain* APV2 = it2->second;
0764 
0765           if (APV2->DetId != APV->DetId)
0766             continue;
0767           if (APV2->APVId == APV->APVId)
0768             continue;
0769 
0770           int bin2 = APV_Charge->GetXaxis()->FindBin(APV2->Index);
0771           TH1D* Proj2 = APV_Charge->ProjectionY(" ", bin2, bin2, "e");
0772           if (Proj2 != nullptr) {
0773             //                   printf("%8i %i--> %4.0f + %4.0f\n",APV2->DetId, APV2->APVId, Proj->GetEntries(), Proj2->GetEntries());
0774             Proj->Add(Proj2, 1);
0775           }
0776         }
0777         //             printf("%8i %i--> %4.0f Full\n",APV->DetId, APV->APVId, Proj->GetEntries());
0778       }
0779 
0780       //std::cout << "Proj->GetEntries(): " << Proj->GetEntries() << ", Proj->GetMean(): " << Proj->GetMean() << std::endl;
0781 
0782       getPeakOfLandau(Proj, FitResults);
0783       APV->MPV = FitResults[0];
0784       //         printf("MPV = %f - %f\n",FitResults[0], FitResults[1]);
0785       if (FitResults[0] != -0.5 && FitResults[1] < MaxMPVError) {
0786         APV_MPV->Fill(APV->Index, APV->MPV);
0787         MPVs->Fill(APV->MPV);
0788         if (APV->Thickness < 0.04)
0789           MPVs320->Fill(APV->MPV);
0790         if (APV->Thickness > 0.04)
0791           MPVs500->Fill(APV->MPV);
0792 
0793         MPV_Vs_R->Fill(APV->R, APV->MPV);
0794         MPV_Vs_Eta->Fill(APV->Eta, APV->MPV);
0795         if (APV->SubDet == StripSubdetector::TIB)
0796           MPV_Vs_EtaTIB->Fill(APV->Eta, APV->MPV);
0797         if (APV->SubDet == StripSubdetector::TID)
0798           MPV_Vs_EtaTID->Fill(APV->Eta, APV->MPV);
0799         if (APV->SubDet == StripSubdetector::TOB)
0800           MPV_Vs_EtaTOB->Fill(APV->Eta, APV->MPV);
0801         if (APV->SubDet == StripSubdetector::TEC) {
0802           MPV_Vs_EtaTEC->Fill(APV->Eta, APV->MPV);
0803           if (APV->Thickness < 0.04)
0804             MPV_Vs_EtaTEC1->Fill(APV->Eta, APV->MPV);
0805           if (APV->Thickness > 0.04)
0806             MPV_Vs_EtaTEC2->Fill(APV->Eta, APV->MPV);
0807         }
0808         MPV_Vs_Phi->Fill(APV->Phi, APV->MPV);
0809         if (APV->SubDet == StripSubdetector::TIB)
0810           MPV_Vs_PhiTIB->Fill(APV->Phi, APV->MPV);
0811         if (APV->SubDet == StripSubdetector::TID)
0812           MPV_Vs_PhiTID->Fill(APV->Phi, APV->MPV);
0813         if (APV->SubDet == StripSubdetector::TOB)
0814           MPV_Vs_PhiTOB->Fill(APV->Phi, APV->MPV);
0815         if (APV->SubDet == StripSubdetector::TEC) {
0816           MPV_Vs_PhiTEC->Fill(APV->Phi, APV->MPV);
0817           if (APV->Thickness < 0.04)
0818             MPV_Vs_PhiTEC1->Fill(APV->Phi, APV->MPV);
0819           if (APV->Thickness > 0.04)
0820             MPV_Vs_PhiTEC2->Fill(APV->Phi, APV->MPV);
0821         }
0822 
0823         if (APV->SubDet == StripSubdetector::TIB)
0824           Charge_TIB->Add(Proj, 1);
0825         if (APV->SubDet == StripSubdetector::TID) {
0826           Charge_TID->Add(Proj, 1);
0827           if (APV->Side == 1)
0828             Charge_TIDM->Add(Proj, 1);
0829           if (APV->Side == 2)
0830             Charge_TIDP->Add(Proj, 1);
0831         }
0832         if (APV->SubDet == StripSubdetector::TOB)
0833           Charge_TOB->Add(Proj, 1);
0834         if (APV->SubDet == StripSubdetector::TEC) {
0835           Charge_TEC->Add(Proj, 1);
0836           if (APV->Thickness < 0.04)
0837             Charge_TEC1->Add(Proj, 1);
0838           if (APV->Thickness > 0.04)
0839             Charge_TEC2->Add(Proj, 1);
0840           if (APV->Side == 1)
0841             Charge_TECM->Add(Proj, 1);
0842           if (APV->Side == 2)
0843             Charge_TECP->Add(Proj, 1);
0844         }
0845       }
0846 
0847       if (APV->SubDet == StripSubdetector::TIB)
0848         Charge_TIB->Add(Proj, 1);
0849       if (APV->SubDet == StripSubdetector::TID) {
0850         Charge_TID->Add(Proj, 1);
0851         if (APV->Side == 1)
0852           Charge_TIDM->Add(Proj, 1);
0853         if (APV->Side == 2)
0854           Charge_TIDP->Add(Proj, 1);
0855       }
0856       if (APV->SubDet == StripSubdetector::TOB)
0857         Charge_TOB->Add(Proj, 1);
0858       if (APV->SubDet == StripSubdetector::TEC) {
0859         Charge_TEC->Add(Proj, 1);
0860         if (APV->Thickness < 0.04)
0861           Charge_TEC1->Add(Proj, 1);
0862         if (APV->Thickness > 0.04)
0863           Charge_TEC2->Add(Proj, 1);
0864         if (APV->Side == 1)
0865           Charge_TECM->Add(Proj, 1);
0866         if (APV->Side == 2)
0867           Charge_TECP->Add(Proj, 1);
0868       }
0869 
0870       if (FitResults[0] != -0.5) {
0871         HChi2OverNDF->Fill(FitResults[4]);
0872         Error_Vs_MPV->Fill(FitResults[0], FitResults[1]);
0873         Error_Vs_Entries->Fill(Proj->GetEntries(), FitResults[1]);
0874         Error_Vs_Eta->Fill(APV->Eta, FitResults[1]);
0875         Error_Vs_Phi->Fill(APV->Phi, FitResults[1]);
0876       }
0877       NumberOfEntriesByAPV->Fill(Proj->GetEntries());
0878       delete Proj;
0879 
0880       Proj = APV_PathLength->ProjectionY(" ", bin, bin, "e");
0881       if (Proj == nullptr)
0882         continue;
0883 
0884       APV_PathLengthM->SetBinContent(APV->Index, Proj->GetMean(1));
0885       APV_PathLengthM->SetBinError(APV->Index, Proj->GetMeanError(1));
0886       //         delete Proj;
0887     }
0888 
0889     unsigned int GOOD = 0;
0890     unsigned int BAD = 0;
0891     double MPVmean = 300.;  //MPVs->GetMean();
0892     for (auto it = APVsColl.begin(); it != APVsColl.end(); it++) {
0893       stAPVGain* APV = it->second;
0894       if (APV->MPV > 0) {
0895         APV->Gain = APV->MPV / MPVmean;  // APV->MPV;
0896         GOOD++;
0897       } else {
0898         NoMPV_Vs_EtaPhi->Fill(APV->Eta, APV->Phi);
0899         APV->Gain = 1;
0900         BAD++;
0901       }
0902       if (APV->Gain <= 0)
0903         APV->Gain = 1;
0904       APV_Gain->Fill(APV->Index, APV->Gain);
0905 
0906       if (!FirstSetOfConstants)
0907         APV->Gain *= APV->PreviousGain;
0908       APV_CumulGain->Fill(APV->Index, APV->Gain);
0909     }
0910 
0911     for (int j = 0; j < Charge_Vs_PathLength->GetXaxis()->GetNbins(); j++) {
0912       Proj = Charge_Vs_PathLength->ProjectionY(" ", j, j, "e");
0913       getPeakOfLandau(Proj, FitResults);
0914       if (FitResults[0] == -0.5)
0915         continue;
0916       MPV_Vs_PathLength->SetBinContent(j, FitResults[0] / Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j));
0917       MPV_Vs_PathLength->SetBinError(j, FitResults[1] / Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j));
0918       FWHM_Vs_PathLength->SetBinContent(
0919           j, FitResults[2] / (FitResults[0] / Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j)));
0920       FWHM_Vs_PathLength->SetBinError(
0921           j, FitResults[3] / (FitResults[0] / Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j)));
0922       delete Proj;
0923     }
0924 
0925     for (int j = 0; j < Charge_Vs_PathLength320->GetXaxis()->GetNbins(); j++) {
0926       Proj = Charge_Vs_PathLength320->ProjectionY(" ", j, j, "e");
0927       getPeakOfLandau(Proj, FitResults);
0928       if (FitResults[0] == -0.5)
0929         continue;
0930       MPV_Vs_PathLength320->SetBinContent(j, FitResults[0] / Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j));
0931       MPV_Vs_PathLength320->SetBinError(j, FitResults[1] / Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j));
0932       FWHM_Vs_PathLength320->SetBinContent(
0933           j, FitResults[2] / (FitResults[0] / Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j)));
0934       FWHM_Vs_PathLength320->SetBinError(
0935           j, FitResults[3] / (FitResults[0] / Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j)));
0936       delete Proj;
0937     }
0938 
0939     for (int j = 0; j < Charge_Vs_PathLength500->GetXaxis()->GetNbins(); j++) {
0940       Proj = Charge_Vs_PathLength500->ProjectionY(" ", j, j, "e");
0941       getPeakOfLandau(Proj, FitResults);
0942       if (FitResults[0] == -0.5)
0943         continue;
0944       MPV_Vs_PathLength500->SetBinContent(j, FitResults[0] / Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j));
0945       MPV_Vs_PathLength500->SetBinError(j, FitResults[1] / Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j));
0946       FWHM_Vs_PathLength500->SetBinContent(
0947           j, FitResults[2] / (FitResults[0] / Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j)));
0948       FWHM_Vs_PathLength500->SetBinError(
0949           j, FitResults[3] / (FitResults[0] / Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j)));
0950       delete Proj;
0951     }
0952     /*
0953       for(int j=0;j<Charge_Vs_PathLength_CS1->GetXaxis()->GetNbins();j++){
0954          Proj      = Charge_Vs_PathLength_CS1->ProjectionY(" ",j,j,"e");
0955          getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
0956          MPV_Vs_PathLength_CS1->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength_CS1->GetXaxis()->GetBinCenter(j));
0957          MPV_Vs_PathLength_CS1->SetBinError   (j, FitResults[1]/Charge_Vs_PathLength_CS1->GetXaxis()->GetBinCenter(j));
0958          FWHM_Vs_PathLength_CS1->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength_CS1->GetXaxis()->GetBinCenter(j) ));
0959          FWHM_Vs_PathLength_CS1->SetBinError  (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength_CS1->GetXaxis()->GetBinCenter(j) ));
0960          delete Proj;
0961       }
0962 
0963       for(int j=0;j<Charge_Vs_PathLength_CS2->GetXaxis()->GetNbins();j++){
0964          Proj      = Charge_Vs_PathLength_CS2->ProjectionY(" ",j,j,"e");
0965          getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
0966          MPV_Vs_PathLength_CS2->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength_CS2->GetXaxis()->GetBinCenter(j));
0967          MPV_Vs_PathLength_CS2->SetBinError   (j, FitResults[1]/Charge_Vs_PathLength_CS2->GetXaxis()->GetBinCenter(j));
0968          FWHM_Vs_PathLength_CS2->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength_CS2->GetXaxis()->GetBinCenter(j) ));
0969          FWHM_Vs_PathLength_CS2->SetBinError  (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength_CS2->GetXaxis()->GetBinCenter(j) ));
0970          delete Proj;
0971       }
0972 
0973       for(int j=0;j<Charge_Vs_PathLength_CS3->GetXaxis()->GetNbins();j++){
0974          Proj      = Charge_Vs_PathLength_CS3->ProjectionY(" ",j,j,"e");
0975          getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
0976          MPV_Vs_PathLength_CS3->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength_CS3->GetXaxis()->GetBinCenter(j));
0977          MPV_Vs_PathLength_CS3->SetBinError   (j, FitResults[1]/Charge_Vs_PathLength_CS3->GetXaxis()->GetBinCenter(j));
0978          FWHM_Vs_PathLength_CS3->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength_CS3->GetXaxis()->GetBinCenter(j) ));
0979          FWHM_Vs_PathLength_CS3->SetBinError  (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength_CS3->GetXaxis()->GetBinCenter(j) ));
0980          delete Proj;
0981       }
0982 
0983       for(int j=0;j<Charge_Vs_PathLength_CS4->GetXaxis()->GetNbins();j++){
0984          Proj      = Charge_Vs_PathLength_CS4->ProjectionY(" ",j,j,"e");
0985          getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
0986          MPV_Vs_PathLength_CS4->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength_CS4->GetXaxis()->GetBinCenter(j));
0987          MPV_Vs_PathLength_CS4->SetBinError   (j, FitResults[1]/Charge_Vs_PathLength_CS4->GetXaxis()->GetBinCenter(j));
0988          FWHM_Vs_PathLength_CS4->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength_CS4->GetXaxis()->GetBinCenter(j) ));
0989          FWHM_Vs_PathLength_CS4->SetBinError  (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength_CS4->GetXaxis()->GetBinCenter(j) ));
0990          delete Proj;
0991       }
0992 
0993       for(int j=0;j<Charge_Vs_PathLength_CS5->GetXaxis()->GetNbins();j++){
0994          Proj      = Charge_Vs_PathLength_CS5->ProjectionY(" ",j,j,"e");
0995          getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
0996          MPV_Vs_PathLength_CS5->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength_CS5->GetXaxis()->GetBinCenter(j));
0997          MPV_Vs_PathLength_CS5->SetBinError   (j, FitResults[1]/Charge_Vs_PathLength_CS5->GetXaxis()->GetBinCenter(j));
0998          FWHM_Vs_PathLength_CS5->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength_CS5->GetXaxis()->GetBinCenter(j) ));
0999          FWHM_Vs_PathLength_CS5->SetBinError  (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength_CS5->GetXaxis()->GetBinCenter(j) ));
1000          delete Proj;
1001       }
1002 */
1003 
1004     for (int j = 0; j < Charge_Vs_PathTIB->GetXaxis()->GetNbins(); j++) {
1005       Proj = Charge_Vs_PathTIB->ProjectionY(" ", j, j, "e");
1006       getPeakOfLandau(Proj, FitResults);
1007       if (FitResults[0] == -0.5)
1008         continue;
1009       MPV_Vs_PathTIB->SetBinContent(j, FitResults[0] / Charge_Vs_PathTIB->GetXaxis()->GetBinCenter(j));
1010       MPV_Vs_PathTIB->SetBinError(j, FitResults[1] / Charge_Vs_PathTIB->GetXaxis()->GetBinCenter(j));
1011       delete Proj;
1012     }
1013 
1014     for (int j = 0; j < Charge_Vs_PathTID->GetXaxis()->GetNbins(); j++) {
1015       Proj = Charge_Vs_PathTID->ProjectionY(" ", j, j, "e");
1016       getPeakOfLandau(Proj, FitResults);
1017       if (FitResults[0] == -0.5)
1018         continue;
1019       MPV_Vs_PathTID->SetBinContent(j, FitResults[0] / Charge_Vs_PathTID->GetXaxis()->GetBinCenter(j));
1020       MPV_Vs_PathTID->SetBinError(j, FitResults[1] / Charge_Vs_PathTID->GetXaxis()->GetBinCenter(j));
1021       delete Proj;
1022     }
1023 
1024     for (int j = 0; j < Charge_Vs_PathTOB->GetXaxis()->GetNbins(); j++) {
1025       Proj = Charge_Vs_PathTOB->ProjectionY(" ", j, j, "e");
1026       getPeakOfLandau(Proj, FitResults);
1027       if (FitResults[0] == -0.5)
1028         continue;
1029       MPV_Vs_PathTOB->SetBinContent(j, FitResults[0] / Charge_Vs_PathTOB->GetXaxis()->GetBinCenter(j));
1030       MPV_Vs_PathTOB->SetBinError(j, FitResults[1] / Charge_Vs_PathTOB->GetXaxis()->GetBinCenter(j));
1031       delete Proj;
1032     }
1033 
1034     for (int j = 0; j < Charge_Vs_PathTEC->GetXaxis()->GetNbins(); j++) {
1035       Proj = Charge_Vs_PathTEC->ProjectionY(" ", j, j, "e");
1036       getPeakOfLandau(Proj, FitResults);
1037       if (FitResults[0] == -0.5)
1038         continue;
1039       MPV_Vs_PathTEC->SetBinContent(j, FitResults[0] / Charge_Vs_PathTEC->GetXaxis()->GetBinCenter(j));
1040       MPV_Vs_PathTEC->SetBinError(j, FitResults[1] / Charge_Vs_PathTEC->GetXaxis()->GetBinCenter(j));
1041       delete Proj;
1042     }
1043 
1044     for (int j = 0; j < Charge_Vs_PathTEC1->GetXaxis()->GetNbins(); j++) {
1045       Proj = Charge_Vs_PathTEC1->ProjectionY(" ", j, j, "e");
1046       getPeakOfLandau(Proj, FitResults);
1047       if (FitResults[0] == -0.5)
1048         continue;
1049       MPV_Vs_PathTEC1->SetBinContent(j, FitResults[0] / Charge_Vs_PathTEC1->GetXaxis()->GetBinCenter(j));
1050       MPV_Vs_PathTEC1->SetBinError(j, FitResults[1] / Charge_Vs_PathTEC1->GetXaxis()->GetBinCenter(j));
1051       delete Proj;
1052     }
1053 
1054     for (int j = 0; j < Charge_Vs_PathTEC2->GetXaxis()->GetNbins(); j++) {
1055       Proj = Charge_Vs_PathTEC2->ProjectionY(" ", j, j, "e");
1056       getPeakOfLandau(Proj, FitResults);
1057       if (FitResults[0] == -0.5)
1058         continue;
1059       MPV_Vs_PathTEC2->SetBinContent(j, FitResults[0] / Charge_Vs_PathTEC2->GetXaxis()->GetBinCenter(j));
1060       MPV_Vs_PathTEC2->SetBinError(j, FitResults[1] / Charge_Vs_PathTEC2->GetXaxis()->GetBinCenter(j));
1061       delete Proj;
1062     }
1063 
1064     for (int j = 1; j < Charge_Vs_TransversAngle->GetXaxis()->GetNbins(); j++) {
1065       Proj = Charge_Vs_TransversAngle->ProjectionY(" ", j, j, "e");
1066       getPeakOfLandau(Proj, FitResults);
1067       if (FitResults[0] == -0.5)
1068         continue;
1069       MPV_Vs_TransversAngle->SetBinContent(j, FitResults[0]);
1070       MPV_Vs_TransversAngle->SetBinError(j, FitResults[1]);
1071       delete Proj;
1072     }
1073 
1074     for (int j = 1; j < Charge_Vs_Alpha->GetXaxis()->GetNbins(); j++) {
1075       Proj = Charge_Vs_Alpha->ProjectionY(" ", j, j, "e");
1076       getPeakOfLandau(Proj, FitResults);
1077       if (FitResults[0] == -0.5)
1078         continue;
1079       MPV_Vs_Alpha->SetBinContent(j, FitResults[0]);
1080       MPV_Vs_Alpha->SetBinError(j, FitResults[1]);
1081       delete Proj;
1082     }
1083 
1084     for (int j = 1; j < Charge_Vs_Beta->GetXaxis()->GetNbins(); j++) {
1085       Proj = Charge_Vs_Beta->ProjectionY(" ", j, j, "e");
1086       getPeakOfLandau(Proj, FitResults);
1087       if (FitResults[0] == -0.5)
1088         continue;
1089       MPV_Vs_Beta->SetBinContent(j, FitResults[0]);
1090       MPV_Vs_Beta->SetBinError(j, FitResults[1]);
1091       delete Proj;
1092     }
1093 
1094     FILE* Gains = fopen(OutputGains.c_str(), "w");
1095     fprintf(Gains, "NEvents = %i\n", NEvent);
1096     fprintf(Gains, "Number of APVs = %lu\n", static_cast<unsigned long>(APVsColl.size()));
1097     fprintf(Gains, "GoodFits = %i BadFits = %i ratio = %f\n", GOOD, BAD, (100.0 * GOOD) / (GOOD + BAD));
1098     for (std::vector<stAPVGain*>::iterator it = APVsCollOrdered.begin(); it != APVsCollOrdered.end(); it++) {
1099       stAPVGain* APV = *it;
1100       fprintf(Gains,
1101               "%i | %i | PreviousGain = %7.5f NewGain = %7.5f\n",
1102               APV->DetId,
1103               APV->APVId,
1104               APV->PreviousGain,
1105               APV->Gain);
1106     }
1107 
1108     std::vector<int> DetIdOfBuggedAPV;
1109     fprintf(Gains, "----------------------------------------------------------------------\n");
1110     for (std::vector<stAPVGain*>::iterator it = APVsCollOrdered.begin(); it != APVsCollOrdered.end(); it++) {
1111       stAPVGain* APV = *it;
1112       if (APV->MPV > 0 && APV->MPV < 200) {
1113         bool tmpBug = false;
1114         for (unsigned int b = 0; b < DetIdOfBuggedAPV.size() && !tmpBug; b++) {
1115           if (DetIdOfBuggedAPV[b] == APV->DetId)
1116             tmpBug = true;
1117         }
1118         if (!tmpBug) {
1119           fprintf(Gains, "%i,\n", APV->DetId);
1120           DetIdOfBuggedAPV.push_back(APV->DetId);
1121         }
1122       }
1123     }
1124 
1125     fclose(Gains);
1126 
1127     //      delete [] FitResults;
1128     //      delete Proj;
1129   }
1130 
1131   dqmStore_->cd();
1132   dqmStore_->save(OutputHistos);
1133 }
1134 
1135 void SiStripGainFromData::algoBeginRun(const edm::Run&, const edm::EventSetup& iSetup) {
1136   if ((strcmp(AlgoMode.c_str(), "MultiJob") != 0 && !FirstSetOfConstants) || Validation) {
1137     const auto gainHandle = iSetup.getHandle(gainToken_);
1138     if (!gainHandle.isValid()) {
1139       printf("\n#####################\n\nERROR --> gainHandle is not valid\n\n#####################\n\n");
1140       exit(0);
1141     }
1142 
1143     for (std::vector<stAPVGain*>::iterator it = APVsCollOrdered.begin(); it != APVsCollOrdered.end(); it++) {
1144       stAPVGain* APV = *it;
1145 
1146       if (gainHandle.isValid()) {
1147         SiStripApvGain::Range detGainRange = gainHandle->getRange(APV->DetId);
1148         APV->PreviousGain = *(detGainRange.first + APV->APVId);
1149         //             APV_PrevGain->Fill(APV->Index,APV->PreviousGain);
1150         APV_PrevGain->SetBinContent(APV_PrevGain->GetXaxis()->FindBin(APV->Index), APV->PreviousGain);
1151         if (APV->PreviousGain < 0)
1152           APV->PreviousGain = 1;
1153       } else {
1154         printf("GAIN HANDLE IS NOT VALID\n");
1155       }
1156     }
1157   }
1158 }
1159 
1160 void SiStripGainFromData::algoAnalyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
1161   if (strcmp(AlgoMode.c_str(), "WriteOnDB") == 0)
1162     return;
1163 
1164   if (NEvent == 0) {
1165     SRun = iEvent.id().run();
1166     SEvent = iEvent.id().event();
1167     STimestamp = iEvent.time().value();
1168   }
1169   ERun = iEvent.id().run();
1170   EEvent = iEvent.id().event();
1171   ETimestamp = iEvent.time().value();
1172   NEvent++;
1173 
1174   iEvent_ = &iEvent;
1175 
1176   Handle<TrajTrackAssociationCollection> trajTrackAssociationHandle;
1177   iEvent.getByLabel(TrajToTrackProducer, TrajToTrackLabel, trajTrackAssociationHandle);
1178   const TrajTrackAssociationCollection TrajToTrackMap = *trajTrackAssociationHandle.product();
1179 
1180   for (TrajTrackAssociationCollection::const_iterator it = TrajToTrackMap.begin(); it != TrajToTrackMap.end(); ++it) {
1181     const Track track = *it->val;
1182     const Trajectory traj = *it->key;
1183 
1184     if (track.p() < MinTrackMomentum || track.p() > MaxTrackMomentum || track.eta() < MinTrackEta ||
1185         track.eta() > MaxTrackEta)
1186       continue;
1187 
1188     Tracks_Pt_Vs_Eta->Fill(fabs(track.eta()), track.pt());
1189     Tracks_P_Vs_Eta->Fill(fabs(track.eta()), track.p());
1190 
1191     //BEGIN TO COMPUTE NDOF FOR TRACKS NO IMPLEMENTED BEFORE 200pre3
1192     int ndof = 0;
1193     const Trajectory::RecHitContainer transRecHits = traj.recHits();
1194 
1195     for (Trajectory::RecHitContainer::const_iterator rechit = transRecHits.begin(); rechit != transRecHits.end();
1196          ++rechit)
1197       if ((*rechit)->isValid())
1198         ndof += (*rechit)->dimension();
1199     ndof -= 5;
1200     //END TO COMPUTE NDOF FOR TRACKS NO IMPLEMENTED BEFORE 200pre3
1201 
1202     HTrackChi2OverNDF->Fill(traj.chiSquared() / ndof);
1203     if (traj.chiSquared() / ndof > MaxTrackChiOverNdf)
1204       continue;
1205 
1206     vector<TrajectoryMeasurement> measurements = traj.measurements();
1207     HTrackHits->Fill(traj.foundHits());
1208     if (traj.foundHits() < (int)MinTrackHits)
1209       continue;
1210     /*
1211       //BEGIN TO COMPUTE #MATCHEDRECHIT IN THE TRACK
1212       int NMatchedHit = 0;
1213       for(vector<TrajectoryMeasurement>::const_iterator measurement_it = measurements.begin(); measurement_it!=measurements.end(); measurement_it++){
1214          TrajectoryStateOnSurface trajState = measurement_it->updatedState();
1215          if( !trajState.isValid() ) continue;
1216          const TrackingRecHit*         hit               = (*measurement_it->recHit()).hit();
1217          const SiStripMatchedRecHit2D* sistripmatchedhit = dynamic_cast<const SiStripMatchedRecHit2D*>(hit);
1218      if(sistripmatchedhit)NMatchedHit++;
1219 //       NMatchedHit++;
1220 
1221       }
1222       //END TO COMPUTE #MATCHEDRECHIT IN THE TRACK
1223 
1224       if(NMatchedHit<2){
1225      printf("NOT ENOUGH MATCHED RECHITS : %i\n",NMatchedHit);
1226      continue;
1227       }
1228 */
1229 
1230     for (vector<TrajectoryMeasurement>::const_iterator measurement_it = measurements.begin();
1231          measurement_it != measurements.end();
1232          measurement_it++) {
1233       TrajectoryStateOnSurface trajState = measurement_it->updatedState();
1234       if (!trajState.isValid())
1235         continue;
1236 
1237       const TrackingRecHit* hit = (*measurement_it->recHit()).hit();
1238       const SiStripRecHit1D* sistripsimple1dhit = dynamic_cast<const SiStripRecHit1D*>(hit);
1239       const SiStripRecHit2D* sistripsimplehit = dynamic_cast<const SiStripRecHit2D*>(hit);
1240       const SiStripMatchedRecHit2D* sistripmatchedhit = dynamic_cast<const SiStripMatchedRecHit2D*>(hit);
1241 
1242       if (sistripsimplehit) {
1243         ComputeChargeOverPath(
1244             (sistripsimplehit->cluster()).get(), trajState, &iSetup, &track, traj.chiSquared() / ndof);
1245       } else if (sistripmatchedhit) {
1246         ComputeChargeOverPath(&sistripmatchedhit->monoCluster(), trajState, &iSetup, &track, traj.chiSquared() / ndof);
1247         ComputeChargeOverPath(
1248             &sistripmatchedhit->stereoCluster(), trajState, &iSetup, &track, traj.chiSquared() / ndof);
1249       } else if (sistripsimple1dhit) {
1250         ComputeChargeOverPath(
1251             (sistripsimple1dhit->cluster()).get(), trajState, &iSetup, &track, traj.chiSquared() / ndof);
1252       } else {
1253       }
1254     }
1255   }
1256 }
1257 
1258 double SiStripGainFromData::ComputeChargeOverPath(const SiStripCluster* Cluster,
1259                                                   TrajectoryStateOnSurface trajState,
1260                                                   const edm::EventSetup* iSetup,
1261                                                   const Track* track,
1262                                                   double trajChi2OverN) {
1263   LocalVector trackDirection = trajState.localDirection();
1264   double cosine = trackDirection.z() / trackDirection.mag();
1265   auto const& Ampls = Cluster->amplitudes();
1266   uint32_t DetId = 0;  // is 0 since long time Cluster->geographicalId();
1267   int FirstStrip = Cluster->firstStrip();
1268   int APVId = FirstStrip / 128;
1269   stAPVGain* APV = APVsColl[(DetId << 3) | APVId];
1270   int Saturation = 0;
1271   bool Overlaping = false;
1272   int Charge = 0;
1273   unsigned int NHighStrip = 0;
1274 
1275   if (!IsFarFromBorder(trajState, DetId, iSetup))
1276     return -1;
1277 
1278   if (FirstStrip == 0)
1279     Overlaping = true;
1280   if (FirstStrip == 128)
1281     Overlaping = true;
1282   if (FirstStrip == 256)
1283     Overlaping = true;
1284   if (FirstStrip == 384)
1285     Overlaping = true;
1286   if (FirstStrip == 512)
1287     Overlaping = true;
1288   if (FirstStrip == 640)
1289     Overlaping = true;
1290 
1291   if (FirstStrip <= 127 && FirstStrip + Ampls.size() > 127)
1292     Overlaping = true;
1293   if (FirstStrip <= 255 && FirstStrip + Ampls.size() > 255)
1294     Overlaping = true;
1295   if (FirstStrip <= 383 && FirstStrip + Ampls.size() > 383)
1296     Overlaping = true;
1297   if (FirstStrip <= 511 && FirstStrip + Ampls.size() > 511)
1298     Overlaping = true;
1299   if (FirstStrip <= 639 && FirstStrip + Ampls.size() > 639)
1300     Overlaping = true;
1301 
1302   if (FirstStrip + Ampls.size() == 127)
1303     Overlaping = true;
1304   if (FirstStrip + Ampls.size() == 255)
1305     Overlaping = true;
1306   if (FirstStrip + Ampls.size() == 383)
1307     Overlaping = true;
1308   if (FirstStrip + Ampls.size() == 511)
1309     Overlaping = true;
1310   if (FirstStrip + Ampls.size() == 639)
1311     Overlaping = true;
1312   if (FirstStrip + Ampls.size() == 767)
1313     Overlaping = true;
1314   if (Overlaping)
1315     return -1;
1316 
1317   /*
1318    if(FirstStrip==0                                 )Overlaping=true;
1319    if(FirstStrip<=255 && FirstStrip+Ampls.size()>255)Overlaping=true;
1320    if(FirstStrip<=511 && FirstStrip+Ampls.size()>511)Overlaping=true;
1321    if(FirstStrip+Ampls.size()==511                  )Overlaping=true;
1322    if(FirstStrip+Ampls.size()==767                  )Overlaping=true;
1323    if(Overlaping)return -1;
1324 */
1325 
1326   for (unsigned int a = 0; a < Ampls.size(); a++) {
1327     Charge += Ampls[a];
1328     if (Ampls[a] >= 254)
1329       Saturation++;
1330     if (Ampls[a] >= 20)
1331       NHighStrip++;
1332   }
1333   double path = (10.0 * APV->Thickness) / fabs(cosine);
1334   double ClusterChargeOverPath = (double)Charge / path;
1335 
1336   NSatStripInCluster->Fill(Saturation);
1337 
1338   if (Ampls.size() > MaxNrStrips)
1339     return -1;
1340   if (Saturation > 0 && !AllowSaturation)
1341     return -1;
1342   Charge_Vs_PathLength->Fill(path, Charge);
1343   if (APV->Thickness < 0.04)
1344     Charge_Vs_PathLength320->Fill(path, Charge);
1345   if (APV->Thickness > 0.04)
1346     Charge_Vs_PathLength500->Fill(path, Charge);
1347   if (APV->SubDet == StripSubdetector::TIB)
1348     Charge_Vs_PathTIB->Fill(path, Charge);
1349   if (APV->SubDet == StripSubdetector::TID)
1350     Charge_Vs_PathTID->Fill(path, Charge);
1351   if (APV->SubDet == StripSubdetector::TOB)
1352     Charge_Vs_PathTOB->Fill(path, Charge);
1353   if (APV->SubDet == StripSubdetector::TEC) {
1354     Charge_Vs_PathTEC->Fill(path, Charge);
1355     if (APV->Thickness < 0.04)
1356       Charge_Vs_PathTEC1->Fill(path, Charge);
1357     if (APV->Thickness > 0.04)
1358       Charge_Vs_PathTEC2->Fill(path, Charge);
1359   }
1360 
1361   double trans = atan2(trackDirection.y(), trackDirection.x()) * (180 / 3.14159265);
1362   double alpha =
1363       acos(trackDirection.x() / sqrt(pow(trackDirection.x(), 2) + pow(trackDirection.z(), 2))) * (180 / 3.14159265);
1364   double beta =
1365       acos(trackDirection.y() / sqrt(pow(trackDirection.x(), 2) + pow(trackDirection.z(), 2))) * (180 / 3.14159265);
1366 
1367   if (path > 0.4 && path < 0.45) {
1368     Charge_Vs_TransversAngle->Fill(trans, Charge / path);
1369     Charge_Vs_Alpha->Fill(alpha, Charge / path);
1370     Charge_Vs_Beta->Fill(beta, Charge / path);
1371   }
1372 
1373   NStrips_Vs_TransversAngle->Fill(trans, Ampls.size());
1374   NStrips_Vs_Alpha->Fill(alpha, Ampls.size());
1375   NStrips_Vs_Beta->Fill(beta, Ampls.size());
1376 
1377   NHighStripInCluster->Fill(NHighStrip);
1378   //   if(NHighStrip==1)   Charge_Vs_PathLength_CS1->Fill(path, Charge );
1379   //   if(NHighStrip==2)   Charge_Vs_PathLength_CS2->Fill(path, Charge );
1380   //   if(NHighStrip==3)   Charge_Vs_PathLength_CS3->Fill(path, Charge );
1381   //   if(NHighStrip==4)   Charge_Vs_PathLength_CS4->Fill(path, Charge );
1382   //   if(NHighStrip==5)   Charge_Vs_PathLength_CS5->Fill(path, Charge );
1383 
1384   HFirstStrip->Fill(FirstStrip);
1385 
1386   if (Validation) {
1387     ClusterChargeOverPath = ClusterChargeOverPath / APV->PreviousGain;
1388   }
1389 
1390   APV_Charge->Fill(APV->Index, ClusterChargeOverPath);
1391   APV_Momentum->Fill(APV->Index, trajState.globalMomentum().mag());
1392   APV_PathLength->Fill(APV->Index, path);
1393 
1394   return ClusterChargeOverPath;
1395 }
1396 
1397 bool SiStripGainFromData::IsFarFromBorder(TrajectoryStateOnSurface trajState,
1398                                           const uint32_t detid,
1399                                           const edm::EventSetup* iSetup) {
1400   const auto& tkGeom = iSetup->getData(tkGeomToken_);
1401 
1402   LocalPoint HitLocalPos = trajState.localPosition();
1403   LocalError HitLocalError = trajState.localError().positionError();
1404 
1405   const GeomDetUnit* it = tkGeom.idToDetUnit(DetId(detid));
1406   if (dynamic_cast<const StripGeomDetUnit*>(it) == nullptr && dynamic_cast<const PixelGeomDetUnit*>(it) == nullptr) {
1407     std::cout << "this detID doesn't seem to belong to the Tracker" << std::endl;
1408     return false;
1409   }
1410 
1411   const BoundPlane plane = it->surface();
1412   const TrapezoidalPlaneBounds* trapezoidalBounds(dynamic_cast<const TrapezoidalPlaneBounds*>(&(plane.bounds())));
1413   const RectangularPlaneBounds* rectangularBounds(dynamic_cast<const RectangularPlaneBounds*>(&(plane.bounds())));
1414 
1415   double DistFromBorder = 1.0;
1416   //double HalfWidth      = it->surface().bounds().width()  /2.0;
1417   double HalfLength;
1418 
1419   if (trapezoidalBounds) {
1420     std::array<const float, 4> const& parameters = (*trapezoidalBounds).parameters();
1421     HalfLength = parameters[3];
1422     //double t       = (HalfLength + HitLocalPos.y()) / (2*HalfLength) ;
1423     //HalfWidth      = parameters[0] + (parameters[1]-parameters[0]) * t;
1424   } else if (rectangularBounds) {
1425     //HalfWidth      = it->surface().bounds().width()  /2.0;
1426     HalfLength = it->surface().bounds().length() / 2.0;
1427   } else {
1428     return false;
1429   }
1430 
1431   //  if (fabs(HitLocalPos.x())+HitLocalError.xx() >= (HalfWidth  - DistFromBorder) ) return false;//Don't think is really necessary
1432   if (fabs(HitLocalPos.y()) + HitLocalError.yy() >= (HalfLength - DistFromBorder))
1433     return false;
1434 
1435   return true;
1436 }
1437 
1438 void SiStripGainFromData::getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange, double HighRange) {
1439   double adcs = -0.5;
1440   double adcs_err = 0.;
1441   double width = -0.5;
1442   double width_err = 0;
1443   double chi2overndf = -0.5;
1444 
1445   double nr_of_entries = InputHisto->GetEntries();
1446 
1447   if ((unsigned int)nr_of_entries < MinNrEntries) {
1448     FitResults[0] = adcs;
1449     FitResults[1] = adcs_err;
1450     FitResults[2] = width;
1451     FitResults[3] = width_err;
1452     FitResults[4] = chi2overndf;
1453     return;
1454   }
1455 
1456   // perform fit with standard landau
1457   TF1* MyLandau = new TF1("MyLandau", "landau", LowRange, HighRange);
1458   MyLandau->SetParameter("MPV", 300);
1459 
1460   InputHisto->Fit(MyLandau, "QR WW");
1461 
1462   // MPV is parameter 1 (0=constant, 1=MPV, 2=Sigma)
1463   adcs = MyLandau->GetParameter("MPV");
1464   adcs_err = MyLandau->GetParError(1);
1465   width = MyLandau->GetParameter(2);
1466   width_err = MyLandau->GetParError(2);
1467   chi2overndf = MyLandau->GetChisquare() / MyLandau->GetNDF();
1468 
1469   // if still wrong, give up
1470   if (adcs < 2. || chi2overndf > MaxChi2OverNDF) {
1471     adcs = -0.5;
1472     adcs_err = 0.;
1473     width = -0.5;
1474     width_err = 0;
1475     chi2overndf = -0.5;
1476   }
1477 
1478   FitResults[0] = adcs;
1479   FitResults[1] = adcs_err;
1480   FitResults[2] = width;
1481   FitResults[3] = width_err;
1482   FitResults[4] = chi2overndf;
1483 
1484   delete MyLandau;
1485 }
1486 
1487 std::unique_ptr<SiStripApvGain> SiStripGainFromData::getNewObject() {
1488   cout << "START getNewObject\n";
1489 
1490   //  if( !(strcmp(AlgoMode.c_str(),"WriteOnDB")==0 || strcmp(AlgoMode.c_str(),"SingleJob")==0) )return NULL;
1491   if (!(strcmp(AlgoMode.c_str(), "WriteOnDB") == 0 || strcmp(AlgoMode.c_str(), "SingleJob") == 0))
1492     return std::make_unique<SiStripApvGain>();
1493 
1494   auto obj = std::make_unique<SiStripApvGain>();
1495   std::vector<float>* theSiStripVector = nullptr;
1496   int PreviousDetId = -1;
1497   for (unsigned int a = 0; a < APVsCollOrdered.size(); a++) {
1498     stAPVGain* APV = APVsCollOrdered[a];
1499     if (APV == nullptr) {
1500       printf("Bug\n");
1501       continue;
1502     }
1503     if (APV->DetId != PreviousDetId) {
1504       if (theSiStripVector != nullptr) {
1505         SiStripApvGain::Range range(theSiStripVector->begin(), theSiStripVector->end());
1506         if (!obj->put(PreviousDetId, range))
1507           printf("Bug to put detId = %i\n", PreviousDetId);
1508       }
1509 
1510       theSiStripVector = new std::vector<float>;
1511       PreviousDetId = APV->DetId;
1512     }
1513     printf("%i | %i | PreviousGain = %7.5f NewGain = %7.5f\n", APV->DetId, APV->APVId, APV->PreviousGain, APV->Gain);
1514     if (theSiStripVector != nullptr) {
1515       theSiStripVector->push_back(APV->Gain);
1516     }
1517     //      theSiStripVector->push_back(APV->Gain);
1518   }
1519 
1520   if (theSiStripVector != nullptr) {
1521     SiStripApvGain::Range range(theSiStripVector->begin(), theSiStripVector->end());
1522     if (!obj->put(PreviousDetId, range))
1523       printf("Bug to put detId = %i\n", PreviousDetId);
1524   }
1525 
1526   cout << "END getNewObject\n";
1527   return obj;
1528 }
1529 
1530 DEFINE_FWK_MODULE(SiStripGainFromData);