00001
00002
00003
00004 #include <memory>
00005
00006 #include "FWCore/Framework/interface/Frameworkfwd.h"
00007 #include "FWCore/Framework/interface/EDAnalyzer.h"
00008 #include "FWCore/Framework/interface/Event.h"
00009 #include "FWCore/Framework/interface/MakerMacros.h"
00010 #include "FWCore/Framework/interface/ESHandle.h"
00011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00013 #include "FWCore/Utilities/interface/Exception.h"
00014
00015 #include "Geometry/CommonDetUnit/interface/GeomDetUnit.h"
00016 #include "Geometry/CommonDetUnit/interface/GeomDetType.h"
00017 #include "Geometry/CommonTopologies/interface/StripTopology.h"
00018 #include "DataFormats/GeometrySurface/interface/TrapezoidalPlaneBounds.h"
00019 #include "DataFormats/GeometrySurface/interface/RectangularPlaneBounds.h"
00020
00021 #include "Geometry/TrackerGeometryBuilder/interface/TrackerGeometry.h"
00022 #include "Geometry/Records/interface/TrackerDigiGeometryRecord.h"
00023 #include "Geometry/TrackerGeometryBuilder/interface/StripGeomDetUnit.h"
00024 #include "Geometry/TrackerGeometryBuilder/interface/PixelGeomDetUnit.h"
00025 #include "Geometry/TrackerNumberingBuilder/interface/GeometricDet.h"
00026 #include "Geometry/CommonDetUnit/interface/TrackingGeometry.h"
00027
00028 #include "DataFormats/SiStripCluster/interface/SiStripClusterCollection.h"
00029
00030 #include "CalibFormats/SiStripObjects/interface/SiStripDetCabling.h"
00031 #include "CalibTracker/Records/interface/SiStripDetCablingRcd.h"
00032
00033 #include "DataFormats/FEDRawData/interface/FEDNumbering.h"
00034 #include "DataFormats/TrackReco/interface/Track.h"
00035 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00036 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit1D.h"
00037 #include "DataFormats/TrackerRecHit2D/interface/SiStripRecHit2D.h"
00038 #include "DataFormats/TrackerRecHit2D/interface/SiStripMatchedRecHit2D.h"
00039 #include "DataFormats/SiStripDetId/interface/SiStripSubStructure.h"
00040 #include "DataFormats/DetId/interface/DetId.h"
00041 #include "DataFormats/SiStripDetId/interface/StripSubdetector.h"
00042 #include "DataFormats/SiStripDetId/interface/TECDetId.h"
00043 #include "DataFormats/SiStripDetId/interface/TIBDetId.h"
00044 #include "DataFormats/SiStripDetId/interface/TIDDetId.h"
00045 #include "DataFormats/SiStripDetId/interface/TOBDetId.h"
00046 #include "DataFormats/TrackReco/interface/DeDxHit.h"
00047 #include "DataFormats/TrackReco/interface/TrackDeDxHits.h"
00048
00049 #include "CommonTools/ConditionDBWriter/interface/ConditionDBWriter.h"
00050 #include "CondFormats/SiStripObjects/interface/SiStripApvGain.h"
00051
00052 #include "TrackingTools/PatternTools/interface/Trajectory.h"
00053 #include "TrackingTools/PatternTools/interface/TrajTrackAssociation.h"
00054
00055 #include "DQMServices/Core/interface/DQMStore.h"
00056 #include "DQMServices/Core/interface/MonitorElement.h"
00057
00058 #include "CalibFormats/SiStripObjects/interface/SiStripGain.h"
00059 #include "CalibTracker/Records/interface/SiStripGainRcd.h"
00060
00061
00062 #include "TFile.h"
00063 #include "TObjString.h"
00064 #include "TString.h"
00065 #include "TH1F.h"
00066 #include "TH2F.h"
00067 #include "TProfile.h"
00068 #include "TF1.h"
00069 #include "TROOT.h"
00070
00071 #include <ext/hash_map>
00072
00073
00074
00075 using namespace edm;
00076 using namespace reco;
00077 using namespace std;
00078 using namespace __gnu_cxx;
00079 using __gnu_cxx::hash_map;
00080 using __gnu_cxx::hash;
00081
00082
00083 struct stAPVGain{unsigned int Index; int DetId; int APVId; int SubDet; float Eta; float R; float Phi; float Thickness; double MPV; double Gain; double PreviousGain; char Side;};
00084
00085 class SiStripGainFromData : public ConditionDBWriter<SiStripApvGain> {
00086 public:
00087 explicit SiStripGainFromData(const edm::ParameterSet&);
00088 ~SiStripGainFromData();
00089
00090
00091 private:
00092 virtual void algoBeginJob(const edm::EventSetup&) ;
00093 virtual void algoEndJob() ;
00094 virtual void algoBeginRun(const edm::Run &, const edm::EventSetup &);
00095
00096 virtual void algoAnalyze(const edm::Event &, const edm::EventSetup &);
00097
00098 SiStripApvGain* getNewObject();
00099 DQMStore* dqmStore_;
00100 DQMStore* dqmStore_infile;
00101
00102 double ComputeChargeOverPath(const SiStripCluster* Cluster,TrajectoryStateOnSurface trajState, const edm::EventSetup* iSetup, const Track* track, double trajChi2OverN);
00103 bool IsFarFromBorder(TrajectoryStateOnSurface trajState, const uint32_t detid, const edm::EventSetup* iSetup);
00104
00105 void getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange=0, double HighRange=5400);
00106
00107 const edm::EventSetup* iSetup_;
00108 const edm::Event* iEvent_;
00109
00110 bool CheckLocalAngle;
00111 unsigned int MinNrEntries;
00112 double MaxMPVError;
00113 double MaxChi2OverNDF;
00114 double MinTrackMomentum;
00115 double MaxTrackMomentum;
00116 double MinTrackEta;
00117 double MaxTrackEta;
00118 unsigned int MaxNrStrips;
00119 unsigned int MinTrackHits;
00120 double MaxTrackChiOverNdf;
00121 bool AllowSaturation;
00122 bool FirstSetOfConstants;
00123 bool Validation;
00124 int CalibrationLevel;
00125 bool CheckIfFileExist;
00126
00127 std::string AlgoMode;
00128 std::string OutputGains;
00129 std::string OutputHistos;
00130 std::string TrajToTrackProducer;
00131 std::string TrajToTrackLabel;
00132
00133 vector<string> VInputFiles;
00134
00135 MonitorElement* tmp;
00136
00137 TH2F* Tracks_P_Vs_Eta;
00138 TH2F* Tracks_Pt_Vs_Eta;
00139
00140 TH1F* NumberOfEntriesByAPV;
00141 TH1F* HChi2OverNDF;
00142 TH1F* HTrackChi2OverNDF;
00143 TH1F* HTrackHits;
00144
00145 TH2F* MPV_Vs_EtaTIB;
00146 TH2F* MPV_Vs_EtaTID;
00147 TH2F* MPV_Vs_EtaTOB;
00148 TH2F* MPV_Vs_EtaTEC;
00149 TH2F* MPV_Vs_EtaTEC1;
00150 TH2F* MPV_Vs_EtaTEC2;
00151
00152 TH2F* MPV_Vs_PhiTIB;
00153 TH2F* MPV_Vs_PhiTID;
00154 TH2F* MPV_Vs_PhiTOB;
00155 TH2F* MPV_Vs_PhiTEC;
00156 TH2F* MPV_Vs_PhiTEC1;
00157 TH2F* MPV_Vs_PhiTEC2;
00158
00159 TH2F* Charge_Vs_PathTIB;
00160 TH2F* Charge_Vs_PathTID;
00161 TH2F* Charge_Vs_PathTOB;
00162 TH2F* Charge_Vs_PathTEC;
00163 TH2F* Charge_Vs_PathTEC1;
00164 TH2F* Charge_Vs_PathTEC2;
00165
00166 TH1F* MPV_Vs_PathTIB;
00167 TH1F* MPV_Vs_PathTID;
00168 TH1F* MPV_Vs_PathTOB;
00169 TH1F* MPV_Vs_PathTEC;
00170 TH1F* MPV_Vs_PathTEC1;
00171 TH1F* MPV_Vs_PathTEC2;
00172
00173
00174 TH1F* Charge_TIB;
00175 TH1F* Charge_TID;
00176 TH1F* Charge_TIDP;
00177 TH1F* Charge_TIDM;
00178 TH1F* Charge_TOB;
00179 TH1F* Charge_TEC;
00180 TH1F* Charge_TEC1;
00181 TH1F* Charge_TEC2;
00182 TH1F* Charge_TECP;
00183 TH1F* Charge_TECM;
00184
00185 TH2F* MPV_Vs_Phi;
00186 TH2F* MPV_Vs_Eta;
00187 TH2F* MPV_Vs_R;
00188
00189
00190
00191
00192 TH1F* APV_DetId;
00193 TH1F* APV_Id;
00194 TH1F* APV_Eta;
00195 TH1F* APV_R;
00196 TH1F* APV_SubDet;
00197 TH2F* APV_Momentum;
00198 TH2F* APV_Charge;
00199 TH2F* APV_PathLength;
00200 TH1F* APV_PathLengthM;
00201 TH1F* APV_MPV;
00202 TH1F* APV_Gain;
00203 TH1F* APV_CumulGain;
00204 TH1F* APV_PrevGain;
00205 TH1F* APV_Thickness;
00206
00207 TH1F* MPVs;
00208 TH1F* MPVs320;
00209 TH1F* MPVs500;
00210
00211
00212
00213
00214 TH1F* NSatStripInCluster;
00215 TH1F* NHighStripInCluster;
00216
00217
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232
00233
00234
00235 TH2F* Charge_Vs_PathLength;
00236 TH2F* Charge_Vs_PathLength320;
00237 TH2F* Charge_Vs_PathLength500;
00238
00239 TH1F* MPV_Vs_PathLength;
00240 TH1F* MPV_Vs_PathLength320;
00241 TH1F* MPV_Vs_PathLength500;
00242
00243 TH1F* FWHM_Vs_PathLength;
00244 TH1F* FWHM_Vs_PathLength320;
00245 TH1F* FWHM_Vs_PathLength500;
00246
00247
00248 TH2F* Charge_Vs_TransversAngle;
00249 TH1F* MPV_Vs_TransversAngle;
00250 TH2F* NStrips_Vs_TransversAngle;
00251
00252 TH2F* Charge_Vs_Alpha;
00253 TH1F* MPV_Vs_Alpha;
00254 TH2F* NStrips_Vs_Alpha;
00255
00256 TH2F* Charge_Vs_Beta;
00257 TH1F* MPV_Vs_Beta;
00258 TH2F* NStrips_Vs_Beta;
00259
00260 TH2F* Error_Vs_MPV;
00261 TH2F* Error_Vs_Entries;
00262 TH2F* Error_Vs_Eta;
00263 TH2F* Error_Vs_Phi;
00264
00265 TH2F* NoMPV_Vs_EtaPhi;
00266
00267 TH2F* HitLocalPosition;
00268 TH2F* HitLocalPositionBefCut;
00269
00270 TH1F* JobInfo;
00271
00272 TH1F* HFirstStrip;
00273
00274 unsigned int NEvent;
00275 unsigned int SRun;
00276 unsigned int SEvent;
00277 TimeValue_t STimestamp;
00278 unsigned int ERun;
00279 unsigned int EEvent;
00280 TimeValue_t ETimestamp;
00281
00282 private :
00283 class isEqual{
00284 public:
00285 template <class T> bool operator () (const T& PseudoDetId1, const T& PseudoDetId2) { return PseudoDetId1==PseudoDetId2; }
00286 };
00287
00288 std::vector<stAPVGain*> APVsCollOrdered;
00289 __gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual > APVsColl;
00290 };
00291
00292 SiStripGainFromData::SiStripGainFromData(const edm::ParameterSet& iConfig) : ConditionDBWriter<SiStripApvGain>(iConfig)
00293 {
00294 AlgoMode = iConfig.getParameter<std::string>("AlgoMode");
00295
00296 OutputGains = iConfig.getParameter<std::string>("OutputGains");
00297 OutputHistos = iConfig.getParameter<std::string>("OutputHistos");
00298
00299 TrajToTrackProducer = iConfig.getParameter<std::string>("TrajToTrackProducer");
00300 TrajToTrackLabel = iConfig.getParameter<std::string>("TrajToTrackLabel");
00301
00302 CheckLocalAngle = iConfig.getUntrackedParameter<bool> ("checkLocalAngle" , false);
00303 MinNrEntries = iConfig.getUntrackedParameter<unsigned>("minNrEntries" , 20);
00304 MaxMPVError = iConfig.getUntrackedParameter<double> ("maxMPVError" , 500.0);
00305 MaxChi2OverNDF = iConfig.getUntrackedParameter<double> ("maxChi2OverNDF" , 5.0);
00306 MinTrackMomentum = iConfig.getUntrackedParameter<double> ("minTrackMomentum" , 3.0);
00307 MaxTrackMomentum = iConfig.getUntrackedParameter<double> ("maxTrackMomentum" , 99999.0);
00308 MinTrackEta = iConfig.getUntrackedParameter<double> ("minTrackEta" , -5.0);
00309 MaxTrackEta = iConfig.getUntrackedParameter<double> ("maxTrackEta" , 5.0);
00310 MaxNrStrips = iConfig.getUntrackedParameter<unsigned>("maxNrStrips" , 2);
00311 MinTrackHits = iConfig.getUntrackedParameter<unsigned>("MinTrackHits" , 8);
00312 MaxTrackChiOverNdf = iConfig.getUntrackedParameter<double> ("MaxTrackChiOverNdf" , 3);
00313 AllowSaturation = iConfig.getUntrackedParameter<bool> ("AllowSaturation" , false);
00314 FirstSetOfConstants = iConfig.getUntrackedParameter<bool> ("FirstSetOfConstants", true);
00315 Validation = iConfig.getUntrackedParameter<bool> ("Validation" , false);
00316 CheckIfFileExist = iConfig.getUntrackedParameter<bool> ("CheckIfFileExist" , false);
00317
00318 CalibrationLevel = iConfig.getUntrackedParameter<int> ("CalibrationLevel" , 0);
00319
00320
00321 if( strcmp(AlgoMode.c_str(),"WriteOnDB")==0 )
00322 VInputFiles = iConfig.getParameter<vector<string> >("VInputFiles");
00323
00324 dqmStore_ = edm::Service<DQMStore>().operator->();
00325 dqmStore_infile = edm::Service<DQMStore>().operator->();
00326
00327
00328
00329 }
00330
00331
00332 SiStripGainFromData::~SiStripGainFromData()
00333 {
00334 }
00335
00336
00337 void
00338 SiStripGainFromData::algoBeginJob(const edm::EventSetup& iSetup)
00339 {
00340 iSetup_ = &iSetup;
00341
00342
00343
00344 tmp = dqmStore_->book1D ("JobInfo" , "JobInfo", 20,0,20); JobInfo = tmp->getTH1F();
00345
00346 tmp = dqmStore_->book1D ("APV_DetId" , "APV_DetId" , 72785,0,72784); APV_DetId = tmp->getTH1F();
00347 tmp = dqmStore_->book1D ("APV_Id" , "APV_Id" , 72785,0,72784); APV_Id = tmp->getTH1F();
00348 tmp = dqmStore_->book1D ("APV_Eta" , "APV_Eta" , 72785,0,72784); APV_Eta = tmp->getTH1F();
00349 tmp = dqmStore_->book1D ("APV_R" , "APV_R" , 72785,0,72784); APV_R = tmp->getTH1F();
00350 tmp = dqmStore_->book1D ("APV_SubDet" , "APV_SubDet" , 72785,0,72784); APV_SubDet = tmp->getTH1F();
00351 tmp = dqmStore_->book2D ("APV_Momentum" , "APV_Momentum" , 72785,0,72784, 50,0,200); APV_Momentum = tmp->getTH2F();
00352 tmp = dqmStore_->book2D ("APV_Charge" , "APV_Charge" , 72785,0,72784, 1000,0,2000); APV_Charge = tmp->getTH2F();
00353 tmp = dqmStore_->book2D ("APV_PathLength" , "APV_PathLength" , 72785,0,72784, 100,0.2,1.4); APV_PathLength = tmp->getTH2F();
00354 tmp = dqmStore_->book1D ("APV_PathLengthM", "APV_PathLengthM", 72785,0,72784); APV_PathLengthM = tmp->getTH1F();
00355 tmp = dqmStore_->book1D ("APV_MPV" , "APV_MPV" , 72785,0,72784); APV_MPV = tmp->getTH1F();
00356 tmp = dqmStore_->book1D ("APV_Gain" , "APV_Gain" , 72785,0,72784); APV_Gain = tmp->getTH1F();
00357 tmp = dqmStore_->book1D ("APV_PrevGain" , "APV_PrevGain" , 72785,0,72784); APV_PrevGain = tmp->getTH1F();
00358 tmp = dqmStore_->book1D ("APV_CumulGain" , "APV_CumulGain" , 72785,0,72784); APV_CumulGain = tmp->getTH1F();
00359 tmp = dqmStore_->book1D ("APV_Thickness" , "APV_Thicknes" , 72785,0,72784); APV_Thickness = tmp->getTH1F();
00360
00361
00362 tmp = dqmStore_->book2D ("Tracks_P_Vs_Eta" , "Tracks_P_Vs_Eta" , 30, 0,3,50,0,200); Tracks_P_Vs_Eta = tmp->getTH2F();
00363 tmp = dqmStore_->book2D ("Tracks_Pt_Vs_Eta" , "Tracks_Pt_Vs_Eta", 30, 0,3,50,0,200); Tracks_Pt_Vs_Eta = tmp->getTH2F();
00364
00365 tmp = dqmStore_->book2D ("Charge_Vs_PathTIB" , "Charge_Vs_PathTIB" ,100,0.2,1.4, 500,0,2000); Charge_Vs_PathTIB = tmp->getTH2F();
00366 tmp = dqmStore_->book2D ("Charge_Vs_PathTID" , "Charge_Vs_PathTID" ,100,0.2,1.4, 500,0,2000); Charge_Vs_PathTID = tmp->getTH2F();
00367 tmp = dqmStore_->book2D ("Charge_Vs_PathTOB" , "Charge_Vs_PathTOB" ,100,0.2,1.4, 500,0,2000); Charge_Vs_PathTOB = tmp->getTH2F();
00368 tmp = dqmStore_->book2D ("Charge_Vs_PathTEC" , "Charge_Vs_PathTEC" ,100,0.2,1.4, 500,0,2000); Charge_Vs_PathTEC = tmp->getTH2F();
00369 tmp = dqmStore_->book2D ("Charge_Vs_PathTEC1", "Charge_Vs_PathTEC1",100,0.2,1.4, 500,0,2000); Charge_Vs_PathTEC1 = tmp->getTH2F();
00370 tmp = dqmStore_->book2D ("Charge_Vs_PathTEC2", "Charge_Vs_PathTEC2",100,0.2,1.4, 500,0,2000); Charge_Vs_PathTEC2 = tmp->getTH2F();
00371
00372
00373 tmp = dqmStore_->book1D ("Charge_TIB" , "Charge_TIB" ,1000,0,2000); Charge_TIB = tmp->getTH1F();
00374 tmp = dqmStore_->book1D ("Charge_TID" , "Charge_TID" ,1000,0,2000); Charge_TID = tmp->getTH1F();
00375 tmp = dqmStore_->book1D ("Charge_TID+", "Charge_TID+",1000,0,2000); Charge_TIDP = tmp->getTH1F();
00376 tmp = dqmStore_->book1D ("Charge_TID-", "Charge_TID-",1000,0,2000); Charge_TIDM = tmp->getTH1F();
00377 tmp = dqmStore_->book1D ("Charge_TOB" , "Charge_TOB" ,1000,0,2000); Charge_TOB = tmp->getTH1F();
00378 tmp = dqmStore_->book1D ("Charge_TEC" , "Charge_TEC" ,1000,0,2000); Charge_TEC = tmp->getTH1F();
00379 tmp = dqmStore_->book1D ("Charge_TEC1", "Charge_TEC1",1000,0,2000); Charge_TEC1 = tmp->getTH1F();
00380 tmp = dqmStore_->book1D ("Charge_TEC2", "Charge_TEC2",1000,0,2000); Charge_TEC2 = tmp->getTH1F();
00381 tmp = dqmStore_->book1D ("Charge_TEC+", "Charge_TEC+",1000,0,2000); Charge_TECP = tmp->getTH1F();
00382 tmp = dqmStore_->book1D ("Charge_TEC-", "Charge_TEC-",1000,0,2000); Charge_TECM = tmp->getTH1F();
00383
00384
00385
00386
00387
00388
00389
00390
00391
00392 tmp = dqmStore_->book2D ("Charge_Vs_PathLength" , "Charge_Vs_PathLength" , 100,0.2,1.4, 500,0,2000); Charge_Vs_PathLength = tmp->getTH2F();
00393 tmp = dqmStore_->book2D ("Charge_Vs_PathLength320" , "Charge_Vs_PathLength" , 100,0.2,1.4, 500,0,2000); Charge_Vs_PathLength320 = tmp->getTH2F();
00394 tmp = dqmStore_->book2D ("Charge_Vs_PathLength500" , "Charge_Vs_PathLength" , 100,0.2,1.4, 500,0,2000); Charge_Vs_PathLength500 = tmp->getTH2F();
00395
00396 tmp = dqmStore_->book2D ("Charge_Vs_TransversAngle" , "Charge_Vs_TransversAngle" , 220,-20,200, 500,0,2000); Charge_Vs_TransversAngle = tmp->getTH2F();
00397 tmp = dqmStore_->book2D ("Charge_Vs_Alpha" , "Charge_Vs_Alpha" , 220,-20,200, 500,0,2000); Charge_Vs_Alpha = tmp->getTH2F();
00398 tmp = dqmStore_->book2D ("Charge_Vs_Beta" , "Charge_Vs_Beta" , 220,-20,200, 500,0,2000); Charge_Vs_Beta = tmp->getTH2F();
00399
00400 tmp = dqmStore_->book2D ("NStrips_Vs_TransversAngle", "NStrips_Vs_TransversAngle", 220,-20,200, 10,0,10); NStrips_Vs_TransversAngle = tmp->getTH2F();
00401 tmp = dqmStore_->book2D ("NStrips_Vs_Alpha" , "NStrips_Vs_Alpha" , 220,-20,200, 10,0,10); NStrips_Vs_Alpha = tmp->getTH2F();
00402 tmp = dqmStore_->book2D ("NStrips_Vs_Beta" , "NStrips_Vs_Beta" , 220,-20,200, 10,0,10); NStrips_Vs_Beta = tmp->getTH2F();
00403 tmp = dqmStore_->book1D ("NHighStripInCluster" , "NHighStripInCluster" , 15,0,14); NHighStripInCluster = tmp->getTH1F();
00404 tmp = dqmStore_->book1D ("NSatStripInCluster" , "NSatStripInCluster" , 50,0,50); NSatStripInCluster = tmp->getTH1F();
00405
00406 tmp = dqmStore_->book1D ("TrackChi2OverNDF","TrackChi2OverNDF", 500, 0,10); HTrackChi2OverNDF = tmp->getTH1F();
00407 tmp = dqmStore_->book1D ("TrackHits","TrackHits", 40, 0,40); HTrackHits = tmp->getTH1F();
00408
00409 tmp = dqmStore_->book1D ("FirstStrip","FirstStrip", 800, 0,800); HFirstStrip = tmp->getTH1F();
00410
00411 if( strcmp(AlgoMode.c_str(),"MultiJob")!=0 ){
00412
00413 tmp = dqmStore_->book2D ("MPV_Vs_EtaTIB" , "MPV_Vs_EtaTIB" , 50, -3.0, 3.0, 300, 0, 600); MPV_Vs_EtaTIB = tmp->getTH2F();
00414 tmp = dqmStore_->book2D ("MPV_Vs_EtaTID" , "MPV_Vs_EtaTID" , 50, -3.0, 3.0, 300, 0, 600); MPV_Vs_EtaTID = tmp->getTH2F();
00415 tmp = dqmStore_->book2D ("MPV_Vs_EtaTOB" , "MPV_Vs_EtaTOB" , 50, -3.0, 3.0, 300, 0, 600); MPV_Vs_EtaTOB = tmp->getTH2F();
00416 tmp = dqmStore_->book2D ("MPV_Vs_EtaTEC" , "MPV_Vs_EtaTEC" , 50, -3.0, 3.0, 300, 0, 600); MPV_Vs_EtaTEC = tmp->getTH2F();
00417 tmp = dqmStore_->book2D ("MPV_Vs_EtaTEC1" , "MPV_Vs_EtaTEC1", 50, -3.0, 3.0, 300, 0, 600); MPV_Vs_EtaTEC1 = tmp->getTH2F();
00418 tmp = dqmStore_->book2D ("MPV_Vs_EtaTEC2" , "MPV_Vs_EtaTEC2", 50, -3.0, 3.0, 300, 0, 600); MPV_Vs_EtaTEC2 = tmp->getTH2F();
00419
00420 tmp = dqmStore_->book2D ("MPV_Vs_PhiTIB" , "MPV_Vs_PhiTIB" , 50, -3.2, 3.2, 300, 0, 600); MPV_Vs_PhiTIB = tmp->getTH2F();
00421 tmp = dqmStore_->book2D ("MPV_Vs_PhiTID" , "MPV_Vs_PhiTID" , 50, -3.2, 3.2, 300, 0, 600); MPV_Vs_PhiTID = tmp->getTH2F();
00422 tmp = dqmStore_->book2D ("MPV_Vs_PhiTOB" , "MPV_Vs_PhiTOB" , 50, -3.2, 3.2, 300, 0, 600); MPV_Vs_PhiTOB = tmp->getTH2F();
00423 tmp = dqmStore_->book2D ("MPV_Vs_PhiTEC" , "MPV_Vs_PhiTEC" , 50, -3.2, 3.2, 300, 0, 600); MPV_Vs_PhiTEC = tmp->getTH2F();
00424 tmp = dqmStore_->book2D ("MPV_Vs_PhiTEC1" , "MPV_Vs_PhiTEC1", 50, -3.2, 3.2, 300, 0, 600); MPV_Vs_PhiTEC1 = tmp->getTH2F();
00425 tmp = dqmStore_->book2D ("MPV_Vs_PhiTEC2" , "MPV_Vs_PhiTEC2", 50, -3.2, 3.2, 300, 0, 600); MPV_Vs_PhiTEC2 = tmp->getTH2F();
00426
00427
00428 tmp = dqmStore_->book1D ("MPV_Vs_PathTIB" , "MPV_Vs_PathTIB" ,100,0.2,1.4); MPV_Vs_PathTIB = tmp->getTH1F();
00429 tmp = dqmStore_->book1D ("MPV_Vs_PathTID" , "MPV_Vs_PathTID" ,100,0.2,1.4); MPV_Vs_PathTID = tmp->getTH1F();
00430 tmp = dqmStore_->book1D ("MPV_Vs_PathTOB" , "MPV_Vs_PathTOB" ,100,0.2,1.4); MPV_Vs_PathTOB = tmp->getTH1F();
00431 tmp = dqmStore_->book1D ("MPV_Vs_PathTEC" , "MPV_Vs_PathTEC" ,100,0.2,1.4); MPV_Vs_PathTEC = tmp->getTH1F();
00432 tmp = dqmStore_->book1D ("MPV_Vs_PathTEC1" , "MPV_Vs_PathTEC1" ,100,0.2,1.4); MPV_Vs_PathTEC1 = tmp->getTH1F();
00433 tmp = dqmStore_->book1D ("MPV_Vs_PathTEC2" , "MPV_Vs_PathTEC2" ,100,0.2,1.4); MPV_Vs_PathTEC2 = tmp->getTH1F();
00434
00435 tmp = dqmStore_->book2D ("MPV_Vs_Phi", "MPV_Vs_Phi", 50, -3.2, 3.2 , 300, 0, 600); MPV_Vs_Phi = tmp->getTH2F();
00436 tmp = dqmStore_->book2D ("MPV_Vs_Eta", "MPV_Vs_Eta", 50, -3.0, 3.0 , 300, 0, 600); MPV_Vs_Eta = tmp->getTH2F();
00437 tmp = dqmStore_->book2D ("MPV_Vs_R" , "MPV_Vs_R" , 150, 0.0, 150.0, 300, 0, 600); MPV_Vs_R = tmp->getTH2F();
00438
00439
00440
00441
00442
00443
00444
00445
00446
00447
00448
00449
00450
00451 tmp = dqmStore_->book1D ("MPV_Vs_PathLength" , "MPV_Vs_PathLength" , 100, 0.2, 1.4); MPV_Vs_PathLength = tmp->getTH1F();
00452 tmp = dqmStore_->book1D ("MPV_Vs_PathLength320" , "MPV_Vs_PathLength" , 100, 0.2, 1.4); MPV_Vs_PathLength320 = tmp->getTH1F();
00453 tmp = dqmStore_->book1D ("MPV_Vs_PathLength500" , "MPV_Vs_PathLength" , 100, 0.2, 1.4); MPV_Vs_PathLength500 = tmp->getTH1F();
00454
00455 tmp = dqmStore_->book1D ("FWHM_Vs_PathLength" , "FWHM_Vs_PathLength" , 100, 0.2, 1.4); FWHM_Vs_PathLength = tmp->getTH1F();
00456 tmp = dqmStore_->book1D ("FWHM_Vs_PathLength320" , "FWHM_Vs_PathLength" , 100, 0.2, 1.4); FWHM_Vs_PathLength320 = tmp->getTH1F();
00457 tmp = dqmStore_->book1D ("FWHM_Vs_PathLength500" , "FWHM_Vs_PathLength" , 100, 0.2, 1.4); FWHM_Vs_PathLength500 = tmp->getTH1F();
00458
00459 tmp = dqmStore_->book1D ("MPV_Vs_TransversAngle" , "MPV_Vs_TransversAngle" , 220, -20, 200); MPV_Vs_TransversAngle = tmp->getTH1F();
00460 tmp = dqmStore_->book1D ("MPV_Vs_Alpha" , "MPV_Vs_Alpha" , 220, -20, 200); MPV_Vs_Alpha = tmp->getTH1F();
00461 tmp = dqmStore_->book1D ("MPV_Vs_Beta" , "MPV_Vs_Beta" , 220, -20, 200); MPV_Vs_Beta = tmp->getTH1F();
00462
00463 tmp = dqmStore_->book2D ("Error_Vs_MPV" , "Error_Vs_MPV" ,600,0,600 ,100 ,0 ,50); Error_Vs_MPV = tmp->getTH2F();
00464 tmp = dqmStore_->book2D ("Error_Vs_Entries","Error_Vs_Entries",500,0,10000 ,100 ,0 ,50); Error_Vs_Entries = tmp->getTH2F();
00465 tmp = dqmStore_->book2D ("Error_Vs_Eta" , "Error_Vs_Eta" ,50 ,-3.0,3.0 ,100 ,0 ,50 ); Error_Vs_Eta = tmp->getTH2F();
00466 tmp = dqmStore_->book2D ("Error_Vs_Phi" , "Error_Vs_Phi" ,50 ,-3.2,3.2 ,100 ,0 ,50); Error_Vs_Phi = tmp->getTH2F();
00467
00468
00469 tmp = dqmStore_->book2D ("NoMPV_Vs_EtaPhi" , "NoMPV_Vs_EtaPhi" ,50,-3.0,3.0 ,50 ,-3.2,3.2); NoMPV_Vs_EtaPhi = tmp->getTH2F();
00470
00471
00472
00473 tmp = dqmStore_->book1D ("NumberOfEntriesByAPV" , "NumberOfEntriesByAPV" , 1000, 0,10000); NumberOfEntriesByAPV = tmp->getTH1F();
00474 tmp = dqmStore_->book1D ("Chi2OverNDF","Chi2OverNDF", 500, 0,25); HChi2OverNDF = tmp->getTH1F();
00475
00476 tmp = dqmStore_->book1D ("MPVs", "MPVs", 600,0,600); MPVs = tmp->getTH1F();
00477 tmp = dqmStore_->book1D ("MPVs320", "MPVs320", 600,0,600); MPVs320 = tmp->getTH1F();
00478 tmp = dqmStore_->book1D ("MPVs500", "MPVs500", 600,0,600); MPVs500 = tmp->getTH1F();
00479
00480
00481 }
00482
00483 gROOT->cd();
00484
00485 edm::ESHandle<TrackerGeometry> tkGeom;
00486 iSetup.get<TrackerDigiGeometryRecord>().get( tkGeom );
00487 vector<GeomDet*> Det = tkGeom->dets();
00488
00489
00490 edm::ESHandle<SiStripGain> gainHandle;
00491
00492 iSetup.get<SiStripGainRcd>().get(gainHandle);
00493 if(!gainHandle.isValid()){printf("\n#####################\n\nERROR --> gainHandle is not valid\n\n#####################\n\n");exit(0);}
00494
00495
00496 unsigned int Id=0;
00497 for(unsigned int i=0;i<Det.size();i++){
00498 DetId Detid = Det[i]->geographicalId();
00499 int SubDet = Detid.subdetId();
00500
00501 if( SubDet == StripSubdetector::TIB || SubDet == StripSubdetector::TID ||
00502 SubDet == StripSubdetector::TOB || SubDet == StripSubdetector::TEC ){
00503
00504 StripGeomDetUnit* DetUnit = dynamic_cast<StripGeomDetUnit*> (Det[i]);
00505 if(!DetUnit)continue;
00506
00507 const StripTopology& Topo = DetUnit->specificTopology();
00508 unsigned int NAPV = Topo.nstrips()/128;
00509
00510 double Phi = DetUnit->position().basicVector().phi();
00511 double Eta = DetUnit->position().basicVector().eta();
00512 double R = DetUnit->position().basicVector().transverse();
00513 double Thick = DetUnit->surface().bounds().thickness();
00514
00515 for(unsigned int j=0;j<NAPV;j++){
00516 stAPVGain* APV = new stAPVGain;
00517 APV->Index = Id;
00518 APV->DetId = Detid.rawId();
00519 APV->APVId = j;
00520 APV->SubDet = SubDet;
00521 APV->MPV = -1;
00522 APV->Gain = -1;
00523 APV->PreviousGain = 1;
00524 APV->Eta = Eta;
00525 APV->Phi = Phi;
00526 APV->R = R;
00527 APV->Thickness = Thick;
00528 APV->Side = 0;
00529
00530 if(SubDet==StripSubdetector::TID){
00531 TIDDetId detid = TIDDetId(Detid);
00532 APV->Side = detid.side();
00533 }else if(SubDet==StripSubdetector::TEC){
00534 TECDetId detid = TECDetId(Detid);
00535 APV->Side = detid.side();
00536 }
00537
00538 APVsCollOrdered.push_back(APV);
00539 APVsColl[(APV->DetId<<3) | APV->APVId] = APV;
00540 Id++;
00541
00542 APV_DetId ->Fill(Id,APV->DetId);
00543 APV_Id ->Fill(Id,APV->APVId);
00544 APV_Eta ->Fill(Id,APV->Eta);
00545 APV_R ->Fill(Id,APV->R);
00546 APV_SubDet ->Fill(Id,APV->SubDet);
00547 APV_Thickness->Fill(Id,APV->Thickness);
00548 }
00549 }
00550 }
00551
00552 NEvent = 0;
00553 SRun = 0;
00554 SEvent = 0;
00555 STimestamp = 0;
00556 ERun = 0;
00557 EEvent = 0;
00558 ETimestamp = 0;
00559 }
00560
00561 void
00562 SiStripGainFromData::algoEndJob() {
00563
00564
00565 unsigned int I=0;
00566
00567 if( strcmp(AlgoMode.c_str(),"WriteOnDB")==0 || strcmp(AlgoMode.c_str(),"Merge")==0){
00568 TH1::AddDirectory(kTRUE);
00569
00570 TFile* file = NULL;
00571 for(unsigned int f=0;f<VInputFiles.size();f++){
00572 printf("Loading New Input File : %s\n", VInputFiles[f].c_str());
00573 if(CheckIfFileExist){
00574 FILE* doesFileExist = fopen( VInputFiles[f].c_str(), "r" );
00575 if(!doesFileExist){
00576 printf("File %s doesn't exist\n",VInputFiles[f].c_str());
00577 continue;
00578 }else{
00579 fclose(doesFileExist);
00580 }
00581 }
00582 file = new TFile( VInputFiles[f].c_str() ); if(!file || file->IsZombie() ){printf("### Bug With File %s\n### File will be skipped \n",VInputFiles[f].c_str()); continue;}
00583 APV_Charge ->Add( (TH1*) file->FindObjectAny("APV_Charge") , 1);
00584 APV_Momentum ->Add( (TH1*) file->FindObjectAny("APV_Momentum") , 1);
00585 APV_PathLength ->Add( (TH1*) file->FindObjectAny("APV_PathLength") , 1);
00586
00587 Tracks_P_Vs_Eta ->Add( (TH1*) file->FindObjectAny("Tracks_P_Vs_Eta") , 1);
00588 Tracks_Pt_Vs_Eta ->Add( (TH1*) file->FindObjectAny("Tracks_Pt_Vs_Eta") , 1);
00589
00590 Charge_Vs_PathTIB ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathTIB") , 1);
00591 Charge_Vs_PathTID ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathTID") , 1);
00592 Charge_Vs_PathTOB ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathTOB") , 1);
00593 Charge_Vs_PathTEC ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathTEC") , 1);
00594 Charge_Vs_PathTEC1 ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathTEC1") , 1);
00595 Charge_Vs_PathTEC2 ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathTEC2") , 1);
00596
00597 HTrackChi2OverNDF ->Add( (TH1*) file->FindObjectAny("TrackChi2OverNDF") , 1);
00598 HTrackHits ->Add( (TH1*) file->FindObjectAny("TrackHits") , 1);
00599
00600 NHighStripInCluster ->Add( (TH1*) file->FindObjectAny("NHighStripInCluster") , 1);
00601 NSatStripInCluster ->Add( (TH1*) file->FindObjectAny("NSatStripInCluster") , 1);
00602 Charge_Vs_PathLength ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathLength") , 1);
00603 Charge_Vs_PathLength320 ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathLength320") , 1);
00604 Charge_Vs_PathLength500 ->Add( (TH1*) file->FindObjectAny("Charge_Vs_PathLength500") , 1);
00605 Charge_Vs_TransversAngle ->Add( (TH1*) file->FindObjectAny("Charge_Vs_TransversAngle") , 1);
00606 NStrips_Vs_TransversAngle->Add( (TH1*) file->FindObjectAny("NStrips_Vs_TransversAngle"), 1);
00607 Charge_Vs_Alpha ->Add( (TH1*) file->FindObjectAny("Charge_Vs_Alpha") , 1);
00608 NStrips_Vs_Alpha ->Add( (TH1*) file->FindObjectAny("NStrips_Vs_Alpha") , 1);
00609 HFirstStrip ->Add( (TH1*) file->FindObjectAny("FirstStrip") , 1);
00610
00611 TH1F* JobInfo_tmp = (TH1F*) file->FindObjectAny("JobInfo");
00612 NEvent += (unsigned int) JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(1));
00613 unsigned int tmp_SRun = (unsigned int) JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(3));
00614 unsigned int tmp_SEvent = (unsigned int) JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(4));
00615 unsigned int tmp_ERun = (unsigned int) JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(6));
00616 unsigned int tmp_EEvent = (unsigned int) JobInfo_tmp->GetBinContent(JobInfo_tmp->GetXaxis()->FindBin(7));
00617
00618 if(SRun==0)SRun = tmp_SRun;
00619
00620 if(tmp_SRun< SRun){SRun=tmp_SRun; SEvent=tmp_SEvent;}
00621 else if(tmp_SRun==SRun && tmp_SEvent<SEvent){SEvent=tmp_SEvent;}
00622
00623 if(tmp_ERun> ERun){ERun=tmp_ERun; EEvent=tmp_EEvent;}
00624 else if(tmp_ERun==ERun && tmp_EEvent>EEvent){EEvent=tmp_EEvent;}
00625
00626 printf("Deleting Current Input File\n");
00627 file->Close();
00628 delete file;
00629 }
00630 }
00631
00632 JobInfo->Fill(1,NEvent);
00633 JobInfo->Fill(3,SRun);
00634 JobInfo->Fill(4,SEvent);
00635 JobInfo->Fill(6,ERun);
00636 JobInfo->Fill(7,EEvent);
00637
00638
00639 if( strcmp(AlgoMode.c_str(),"MultiJob")!=0 ){
00640 TH1D* Proj = NULL;
00641 double* FitResults = new double[5];
00642 I=0;
00643 for(__gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual >::iterator it = APVsColl.begin();it!=APVsColl.end();it++){
00644 if( I%3650==0 ) printf("Fitting Histograms \t %6.2f%%\n",(100.0*I)/APVsColl.size());I++;
00645 stAPVGain* APV = it->second;
00646
00647 int bin = APV_Charge->GetXaxis()->FindBin(APV->Index);
00648 Proj = APV_Charge->ProjectionY(" ",bin,bin,"e");
00649 Proj = (TH1D*)Proj->Clone();
00650 if(Proj==NULL)continue;
00651
00652
00653 if(CalibrationLevel==1){
00654 int SecondAPVId = APV->APVId;
00655 if(SecondAPVId%2==0){
00656 SecondAPVId = SecondAPVId+1;
00657 }else{
00658 SecondAPVId = SecondAPVId-1;
00659 }
00660 stAPVGain* APV2 = APVsColl[(APV->DetId<<3) | SecondAPVId];
00661
00662 int bin2 = APV_Charge->GetXaxis()->FindBin(APV2->Index);
00663 TH1D* Proj2 = APV_Charge->ProjectionY(" ",bin2,bin2,"e");
00664 if(Proj2!=NULL){
00665 Proj->Add(Proj2,1);
00666 }
00667 }else if(CalibrationLevel>1){
00668
00669 for(__gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual >::iterator it2 = APVsColl.begin();it2!=APVsColl.end();it2++){
00670 stAPVGain* APV2 = it2->second;
00671
00672 if(APV2->DetId != APV->DetId)continue;
00673 if(APV2->APVId == APV->APVId)continue;
00674
00675 int bin2 = APV_Charge->GetXaxis()->FindBin(APV2->Index);
00676 TH1D* Proj2 = APV_Charge->ProjectionY(" ",bin2,bin2,"e");
00677 if(Proj2!=NULL){
00678
00679 Proj->Add(Proj2,1);
00680 }
00681 }
00682
00683 }
00684
00685
00686
00687
00688 getPeakOfLandau(Proj,FitResults);
00689 APV->MPV = FitResults[0];
00690
00691 if(FitResults[0]!=-0.5 && FitResults[1]<MaxMPVError){
00692 APV_MPV->Fill(APV->Index,APV->MPV);
00693 MPVs ->Fill(APV->MPV);
00694 if(APV->Thickness<0.04) MPVs320->Fill(APV->MPV);
00695 if(APV->Thickness>0.04) MPVs500->Fill(APV->MPV);
00696
00697 MPV_Vs_R->Fill(APV->R,APV->MPV);
00698 MPV_Vs_Eta->Fill(APV->Eta,APV->MPV);
00699 if(APV->SubDet==StripSubdetector::TIB) MPV_Vs_EtaTIB ->Fill(APV->Eta,APV->MPV);
00700 if(APV->SubDet==StripSubdetector::TID) MPV_Vs_EtaTID ->Fill(APV->Eta,APV->MPV);
00701 if(APV->SubDet==StripSubdetector::TOB) MPV_Vs_EtaTOB ->Fill(APV->Eta,APV->MPV);
00702 if(APV->SubDet==StripSubdetector::TEC){ MPV_Vs_EtaTEC ->Fill(APV->Eta,APV->MPV);
00703 if(APV->Thickness<0.04) MPV_Vs_EtaTEC1->Fill(APV->Eta,APV->MPV);
00704 if(APV->Thickness>0.04) MPV_Vs_EtaTEC2->Fill(APV->Eta,APV->MPV);
00705 }
00706 MPV_Vs_Phi->Fill(APV->Phi,APV->MPV);
00707 if(APV->SubDet==StripSubdetector::TIB) MPV_Vs_PhiTIB ->Fill(APV->Phi,APV->MPV);
00708 if(APV->SubDet==StripSubdetector::TID) MPV_Vs_PhiTID ->Fill(APV->Phi,APV->MPV);
00709 if(APV->SubDet==StripSubdetector::TOB) MPV_Vs_PhiTOB ->Fill(APV->Phi,APV->MPV);
00710 if(APV->SubDet==StripSubdetector::TEC){ MPV_Vs_PhiTEC ->Fill(APV->Phi,APV->MPV);
00711 if(APV->Thickness<0.04) MPV_Vs_PhiTEC1->Fill(APV->Phi,APV->MPV);
00712 if(APV->Thickness>0.04) MPV_Vs_PhiTEC2->Fill(APV->Phi,APV->MPV);
00713 }
00714
00715 if(APV->SubDet==StripSubdetector::TIB) Charge_TIB ->Add(Proj,1);
00716 if(APV->SubDet==StripSubdetector::TID){ Charge_TID ->Add(Proj,1);
00717 if(APV->Side==1) Charge_TIDM->Add(Proj,1);
00718 if(APV->Side==2) Charge_TIDP->Add(Proj,1);
00719 }
00720 if(APV->SubDet==StripSubdetector::TOB) Charge_TOB ->Add(Proj,1);
00721 if(APV->SubDet==StripSubdetector::TEC){ Charge_TEC ->Add(Proj,1);
00722 if(APV->Thickness<0.04) Charge_TEC1->Add(Proj,1);
00723 if(APV->Thickness>0.04) Charge_TEC2->Add(Proj,1);
00724 if(APV->Side==1) Charge_TECM->Add(Proj,1);
00725 if(APV->Side==2) Charge_TECP->Add(Proj,1);
00726 }
00727 }
00728
00729 if(APV->SubDet==StripSubdetector::TIB) Charge_TIB ->Add(Proj,1);
00730 if(APV->SubDet==StripSubdetector::TID){ Charge_TID ->Add(Proj,1);
00731 if(APV->Side==1) Charge_TIDM->Add(Proj,1);
00732 if(APV->Side==2) Charge_TIDP->Add(Proj,1);
00733 }
00734 if(APV->SubDet==StripSubdetector::TOB) Charge_TOB ->Add(Proj,1);
00735 if(APV->SubDet==StripSubdetector::TEC){ Charge_TEC ->Add(Proj,1);
00736 if(APV->Thickness<0.04) Charge_TEC1->Add(Proj,1);
00737 if(APV->Thickness>0.04) Charge_TEC2->Add(Proj,1);
00738 if(APV->Side==1) Charge_TECM->Add(Proj,1);
00739 if(APV->Side==2) Charge_TECP->Add(Proj,1);
00740 }
00741
00742 if(FitResults[0]!=-0.5){
00743 HChi2OverNDF->Fill(FitResults[4]);
00744 Error_Vs_MPV->Fill(FitResults[0],FitResults[1]);
00745 Error_Vs_Entries->Fill(Proj->GetEntries(),FitResults[1]);
00746 Error_Vs_Eta->Fill(APV->Eta,FitResults[1]);
00747 Error_Vs_Phi->Fill(APV->Phi,FitResults[1]);
00748 }
00749 NumberOfEntriesByAPV->Fill(Proj->GetEntries());
00750 delete Proj;
00751
00752
00753 Proj = APV_PathLength->ProjectionY(" ",bin,bin,"e");
00754 if(Proj==NULL)continue;
00755
00756 APV_PathLengthM->SetBinContent(APV->Index, Proj->GetMean(1) );
00757 APV_PathLengthM->SetBinError (APV->Index, Proj->GetMeanError(1) );
00758
00759 }
00760
00761 unsigned int GOOD = 0;
00762 unsigned int BAD = 0;
00763 double MPVmean = MPVs->GetMean();
00764 MPVmean = 300;
00765 for(__gnu_cxx::hash_map<unsigned int, stAPVGain*, __gnu_cxx::hash<unsigned int>, isEqual >::iterator it = APVsColl.begin();it!=APVsColl.end();it++){
00766
00767 stAPVGain* APV = it->second;
00768 if(APV->MPV>0){
00769 APV->Gain = APV->MPV / MPVmean;
00770 GOOD++;
00771 }else{
00772 NoMPV_Vs_EtaPhi->Fill(APV->Eta, APV->Phi);
00773 APV->Gain = 1;
00774 BAD++;
00775 }
00776 if(APV->Gain<=0) APV->Gain = 1;
00777 APV_Gain->Fill(APV->Index,APV->Gain);
00778
00779 if(!FirstSetOfConstants) APV->Gain *= APV->PreviousGain;
00780 APV_CumulGain->Fill(APV->Index,APV->Gain);
00781 }
00782
00783 for(int j=0;j<Charge_Vs_PathLength->GetXaxis()->GetNbins();j++){
00784 Proj = Charge_Vs_PathLength->ProjectionY(" ",j,j,"e");
00785 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00786 MPV_Vs_PathLength->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j));
00787 MPV_Vs_PathLength->SetBinError (j, FitResults[1]/Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j));
00788 FWHM_Vs_PathLength->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j)) );
00789 FWHM_Vs_PathLength->SetBinError (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength->GetXaxis()->GetBinCenter(j)) );
00790 delete Proj;
00791 }
00792
00793 for(int j=0;j<Charge_Vs_PathLength320->GetXaxis()->GetNbins();j++){
00794 Proj = Charge_Vs_PathLength320->ProjectionY(" ",j,j,"e");
00795 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00796 MPV_Vs_PathLength320->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j));
00797 MPV_Vs_PathLength320->SetBinError (j, FitResults[1]/Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j));
00798 FWHM_Vs_PathLength320->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j)) );
00799 FWHM_Vs_PathLength320->SetBinError (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength320->GetXaxis()->GetBinCenter(j)) );
00800 delete Proj;
00801 }
00802
00803 for(int j=0;j<Charge_Vs_PathLength500->GetXaxis()->GetNbins();j++){
00804 Proj = Charge_Vs_PathLength500->ProjectionY(" ",j,j,"e");
00805 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00806 MPV_Vs_PathLength500->SetBinContent (j, FitResults[0]/Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j));
00807 MPV_Vs_PathLength500->SetBinError (j, FitResults[1]/Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j));
00808 FWHM_Vs_PathLength500->SetBinContent(j, FitResults[2]/(FitResults[0]/Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j) ));
00809 FWHM_Vs_PathLength500->SetBinError (j, FitResults[3]/(FitResults[0]/Charge_Vs_PathLength500->GetXaxis()->GetBinCenter(j) ));
00810 delete Proj;
00811 }
00812
00813
00814
00815
00816
00817
00818
00819
00820
00821
00822
00823
00824
00825
00826
00827
00828
00829
00830
00831
00832
00833
00834
00835
00836
00837
00838
00839
00840
00841
00842
00843
00844
00845
00846
00847
00848
00849
00850
00851
00852
00853
00854
00855
00856
00857
00858
00859
00860
00861
00862
00863
00864
00865 for(int j=0;j<Charge_Vs_PathTIB->GetXaxis()->GetNbins();j++){
00866 Proj = Charge_Vs_PathTIB->ProjectionY(" ",j,j,"e");
00867 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00868 MPV_Vs_PathTIB->SetBinContent(j, FitResults[0]/Charge_Vs_PathTIB->GetXaxis()->GetBinCenter(j));
00869 MPV_Vs_PathTIB->SetBinError (j, FitResults[1]/Charge_Vs_PathTIB->GetXaxis()->GetBinCenter(j));
00870 delete Proj;
00871 }
00872
00873 for(int j=0;j<Charge_Vs_PathTID->GetXaxis()->GetNbins();j++){
00874 Proj = Charge_Vs_PathTID->ProjectionY(" ",j,j,"e");
00875 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00876 MPV_Vs_PathTID->SetBinContent(j, FitResults[0]/Charge_Vs_PathTID->GetXaxis()->GetBinCenter(j));
00877 MPV_Vs_PathTID->SetBinError (j, FitResults[1]/Charge_Vs_PathTID->GetXaxis()->GetBinCenter(j));
00878 delete Proj;
00879 }
00880
00881 for(int j=0;j<Charge_Vs_PathTOB->GetXaxis()->GetNbins();j++){
00882 Proj = Charge_Vs_PathTOB->ProjectionY(" ",j,j,"e");
00883 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00884 MPV_Vs_PathTOB->SetBinContent(j, FitResults[0]/Charge_Vs_PathTOB->GetXaxis()->GetBinCenter(j));
00885 MPV_Vs_PathTOB->SetBinError (j, FitResults[1]/Charge_Vs_PathTOB->GetXaxis()->GetBinCenter(j));
00886 delete Proj;
00887 }
00888
00889 for(int j=0;j<Charge_Vs_PathTEC->GetXaxis()->GetNbins();j++){
00890 Proj = Charge_Vs_PathTEC->ProjectionY(" ",j,j,"e");
00891 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00892 MPV_Vs_PathTEC->SetBinContent(j, FitResults[0]/Charge_Vs_PathTEC->GetXaxis()->GetBinCenter(j));
00893 MPV_Vs_PathTEC->SetBinError (j, FitResults[1]/Charge_Vs_PathTEC->GetXaxis()->GetBinCenter(j));
00894 delete Proj;
00895 }
00896
00897 for(int j=0;j<Charge_Vs_PathTEC1->GetXaxis()->GetNbins();j++){
00898 Proj = Charge_Vs_PathTEC1->ProjectionY(" ",j,j,"e");
00899 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00900 MPV_Vs_PathTEC1->SetBinContent(j, FitResults[0]/Charge_Vs_PathTEC1->GetXaxis()->GetBinCenter(j));
00901 MPV_Vs_PathTEC1->SetBinError (j, FitResults[1]/Charge_Vs_PathTEC1->GetXaxis()->GetBinCenter(j));
00902 delete Proj;
00903 }
00904
00905 for(int j=0;j<Charge_Vs_PathTEC2->GetXaxis()->GetNbins();j++){
00906 Proj = Charge_Vs_PathTEC2->ProjectionY(" ",j,j,"e");
00907 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00908 MPV_Vs_PathTEC2->SetBinContent(j, FitResults[0]/Charge_Vs_PathTEC2->GetXaxis()->GetBinCenter(j));
00909 MPV_Vs_PathTEC2->SetBinError (j, FitResults[1]/Charge_Vs_PathTEC2->GetXaxis()->GetBinCenter(j));
00910 delete Proj;
00911 }
00912
00913
00914 for(int j=1;j<Charge_Vs_TransversAngle->GetXaxis()->GetNbins();j++){
00915 Proj = Charge_Vs_TransversAngle->ProjectionY(" ",j,j,"e");
00916 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00917 MPV_Vs_TransversAngle->SetBinContent(j, FitResults[0]);
00918 MPV_Vs_TransversAngle->SetBinError (j, FitResults[1]);
00919 delete Proj;
00920 }
00921
00922 for(int j=1;j<Charge_Vs_Alpha->GetXaxis()->GetNbins();j++){
00923 Proj = Charge_Vs_Alpha->ProjectionY(" ",j,j,"e");
00924 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00925 MPV_Vs_Alpha->SetBinContent(j, FitResults[0]);
00926 MPV_Vs_Alpha->SetBinError (j, FitResults[1]);
00927 delete Proj;
00928 }
00929
00930 for(int j=1;j<Charge_Vs_Beta->GetXaxis()->GetNbins();j++){
00931 Proj = Charge_Vs_Beta->ProjectionY(" ",j,j,"e");
00932 getPeakOfLandau(Proj,FitResults); if(FitResults[0] ==-0.5)continue;
00933 MPV_Vs_Beta->SetBinContent(j, FitResults[0]);
00934 MPV_Vs_Beta->SetBinError (j, FitResults[1]);
00935 delete Proj;
00936 }
00937
00938 FILE* Gains = fopen(OutputGains.c_str(),"w");
00939 fprintf(Gains,"NEvents = %i\n",NEvent);
00940 fprintf(Gains,"Number of APVs = %lu\n",static_cast<unsigned long>(APVsColl.size()));
00941 fprintf(Gains,"GoodFits = %i BadFits = %i ratio = %f\n",GOOD,BAD,(100.0*GOOD)/(GOOD+BAD));
00942 for(std::vector<stAPVGain*>::iterator it = APVsCollOrdered.begin();it!=APVsCollOrdered.end();it++){
00943 stAPVGain* APV = *it;
00944 fprintf(Gains,"%i | %i | PreviousGain = %7.5f NewGain = %7.5f\n", APV->DetId,APV->APVId,APV->PreviousGain,APV->Gain);
00945 }
00946
00947 std::vector<int> DetIdOfBuggedAPV;
00948 fprintf(Gains,"----------------------------------------------------------------------\n");
00949 for(std::vector<stAPVGain*>::iterator it = APVsCollOrdered.begin();it!=APVsCollOrdered.end();it++){
00950 stAPVGain* APV = *it;
00951 if(APV->MPV>0 && APV->MPV<200){
00952 bool tmpBug = false;
00953 for(unsigned int b=0;b<DetIdOfBuggedAPV.size()&&!tmpBug;b++){if(DetIdOfBuggedAPV[b]==APV->DetId)tmpBug=true;}
00954 if(!tmpBug){fprintf(Gains,"%i,\n",APV->DetId);DetIdOfBuggedAPV.push_back(APV->DetId);}
00955 }
00956 }
00957
00958
00959 fclose(Gains);
00960
00961
00962
00963 }
00964
00965 dqmStore_->cd();
00966 dqmStore_->save(OutputHistos.c_str());
00967
00968 }
00969
00970
00971 void SiStripGainFromData::algoBeginRun(const edm::Run &, const edm::EventSetup &iSetup){
00972
00973 edm::ESHandle<SiStripGain> gainHandle;
00974 if((strcmp(AlgoMode.c_str(),"MultiJob")!=0 && !FirstSetOfConstants) || Validation){
00975 iSetup.get<SiStripGainRcd>().get(gainHandle);
00976 if(!gainHandle.isValid()){printf("\n#####################\n\nERROR --> gainHandle is not valid\n\n#####################\n\n");exit(0);}
00977
00978
00979 for(std::vector<stAPVGain*>::iterator it = APVsCollOrdered.begin();it!=APVsCollOrdered.end();it++){
00980 stAPVGain* APV = *it;
00981
00982 if(gainHandle.isValid()){
00983 SiStripApvGain::Range detGainRange = gainHandle->getRange(APV->DetId);
00984 APV->PreviousGain = *(detGainRange.first + APV->APVId);
00985
00986 APV_PrevGain->SetBinContent(APV_PrevGain->GetXaxis()->FindBin(APV->Index),APV->PreviousGain);
00987 if(APV->PreviousGain<0)APV->PreviousGain = 1;
00988 }else{
00989 printf("GAIN HANDLE IS NOT VALID\n");
00990 }
00991 }
00992 }
00993
00994
00995 }
00996
00997
00998
00999 void
01000 SiStripGainFromData::algoAnalyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
01001 {
01002
01003 if( strcmp(AlgoMode.c_str(),"WriteOnDB")==0 ) return;
01004
01005 if(NEvent==0){
01006 SRun = iEvent.id().run();
01007 SEvent = iEvent.id().event();
01008 STimestamp = iEvent.time().value();
01009 }
01010 ERun = iEvent.id().run();
01011 EEvent = iEvent.id().event();
01012 ETimestamp = iEvent.time().value();
01013 NEvent++;
01014
01015 iEvent_ = &iEvent;
01016
01017 Handle<TrajTrackAssociationCollection> trajTrackAssociationHandle;
01018 iEvent.getByLabel(TrajToTrackProducer, TrajToTrackLabel, trajTrackAssociationHandle);
01019 const TrajTrackAssociationCollection TrajToTrackMap = *trajTrackAssociationHandle.product();
01020
01021
01022 for(TrajTrackAssociationCollection::const_iterator it = TrajToTrackMap.begin(); it!=TrajToTrackMap.end(); ++it) {
01023 const Track track = *it->val;
01024 const Trajectory traj = *it->key;
01025
01026
01027 if(track.p()<MinTrackMomentum || track.p()>MaxTrackMomentum || track.eta()<MinTrackEta || track.eta()>MaxTrackEta) continue;
01028
01029 Tracks_Pt_Vs_Eta->Fill(fabs(track.eta()),track.pt());
01030 Tracks_P_Vs_Eta ->Fill(fabs(track.eta()),track.p());
01031
01032
01033 int ndof =0;
01034 const Trajectory::RecHitContainer transRecHits = traj.recHits();
01035
01036 for(Trajectory::RecHitContainer::const_iterator rechit = transRecHits.begin(); rechit != transRecHits.end(); ++rechit)
01037 if ((*rechit)->isValid()) ndof += (*rechit)->dimension();
01038 ndof -= 5;
01039
01040
01041 HTrackChi2OverNDF->Fill(traj.chiSquared()/ndof);
01042 if(traj.chiSquared()/ndof>MaxTrackChiOverNdf)continue;
01043
01044 vector<TrajectoryMeasurement> measurements = traj.measurements();
01045 HTrackHits->Fill(traj.foundHits());
01046 if(traj.foundHits()<(int)MinTrackHits)continue;
01047
01048
01049
01050
01051
01052
01053
01054
01055
01056
01057
01058
01059
01060
01061
01062
01063
01064
01065
01066
01067 for(vector<TrajectoryMeasurement>::const_iterator measurement_it = measurements.begin(); measurement_it!=measurements.end(); measurement_it++){
01068
01069 TrajectoryStateOnSurface trajState = measurement_it->updatedState();
01070 if( !trajState.isValid() ) continue;
01071
01072 const TrackingRecHit* hit = (*measurement_it->recHit()).hit();
01073 const SiStripRecHit1D* sistripsimple1dhit = dynamic_cast<const SiStripRecHit1D*>(hit);
01074 const SiStripRecHit2D* sistripsimplehit = dynamic_cast<const SiStripRecHit2D*>(hit);
01075 const SiStripMatchedRecHit2D* sistripmatchedhit = dynamic_cast<const SiStripMatchedRecHit2D*>(hit);
01076
01077 if(sistripsimplehit){
01078 ComputeChargeOverPath((sistripsimplehit->cluster()).get(), trajState, &iSetup, &track, traj.chiSquared()/ndof);
01079 }else if(sistripmatchedhit){
01080 ComputeChargeOverPath((sistripmatchedhit->monoHit() ->cluster()).get(),trajState, &iSetup, &track, traj.chiSquared()/ndof);
01081 ComputeChargeOverPath((sistripmatchedhit->stereoHit()->cluster()).get(),trajState, &iSetup, &track, traj.chiSquared()/ndof);
01082 }else if(sistripsimple1dhit){
01083 ComputeChargeOverPath((sistripsimple1dhit->cluster()).get(), trajState, &iSetup, &track, traj.chiSquared()/ndof);
01084 }else{
01085 }
01086
01087 }
01088
01089 }
01090 }
01091
01092 double
01093 SiStripGainFromData::ComputeChargeOverPath(const SiStripCluster* Cluster ,TrajectoryStateOnSurface trajState, const edm::EventSetup* iSetup, const Track* track, double trajChi2OverN)
01094 {
01095 LocalVector trackDirection = trajState.localDirection();
01096 double cosine = trackDirection.z()/trackDirection.mag();
01097
01098
01099 const vector<uint8_t>& Ampls = Cluster->amplitudes();
01100 uint32_t DetId = Cluster->geographicalId();
01101 int FirstStrip = Cluster->firstStrip();
01102 int APVId = FirstStrip/128;
01103 stAPVGain* APV = APVsColl[(DetId<<3) | APVId];
01104 int Saturation = 0;
01105 bool Overlaping = false;
01106 int Charge = 0;
01107 unsigned int NHighStrip = 0;
01108
01109 if(!IsFarFromBorder(trajState,DetId, iSetup))return -1;
01110
01111
01112 if(FirstStrip==0 )Overlaping=true;
01113 if(FirstStrip==128 )Overlaping=true;
01114 if(FirstStrip==256 )Overlaping=true;
01115 if(FirstStrip==384 )Overlaping=true;
01116 if(FirstStrip==512 )Overlaping=true;
01117 if(FirstStrip==640 )Overlaping=true;
01118
01119 if(FirstStrip<=127 && FirstStrip+Ampls.size()>127)Overlaping=true;
01120 if(FirstStrip<=255 && FirstStrip+Ampls.size()>255)Overlaping=true;
01121 if(FirstStrip<=383 && FirstStrip+Ampls.size()>383)Overlaping=true;
01122 if(FirstStrip<=511 && FirstStrip+Ampls.size()>511)Overlaping=true;
01123 if(FirstStrip<=639 && FirstStrip+Ampls.size()>639)Overlaping=true;
01124
01125 if(FirstStrip+Ampls.size()==127 )Overlaping=true;
01126 if(FirstStrip+Ampls.size()==255 )Overlaping=true;
01127 if(FirstStrip+Ampls.size()==383 )Overlaping=true;
01128 if(FirstStrip+Ampls.size()==511 )Overlaping=true;
01129 if(FirstStrip+Ampls.size()==639 )Overlaping=true;
01130 if(FirstStrip+Ampls.size()==767 )Overlaping=true;
01131 if(Overlaping)return -1;
01132
01133
01134
01135
01136
01137
01138
01139
01140
01141
01142
01143 for(unsigned int a=0;a<Ampls.size();a++){Charge+=Ampls[a];if(Ampls[a]>=254)Saturation++;if(Ampls[a]>=20)NHighStrip++;}
01144 double path = (10.0*APV->Thickness)/fabs(cosine);
01145 double ClusterChargeOverPath = (double)Charge / path ;
01146
01147 NSatStripInCluster->Fill(Saturation);
01148
01149 if(Ampls.size()>MaxNrStrips) return -1;
01150 if(Saturation>0 && !AllowSaturation)return -1;
01151 Charge_Vs_PathLength ->Fill(path,Charge);
01152 if(APV->Thickness<0.04) Charge_Vs_PathLength320->Fill(path,Charge);
01153 if(APV->Thickness>0.04) Charge_Vs_PathLength500->Fill(path,Charge);
01154 if(APV->SubDet==StripSubdetector::TIB) Charge_Vs_PathTIB ->Fill(path,Charge);
01155 if(APV->SubDet==StripSubdetector::TID) Charge_Vs_PathTID ->Fill(path,Charge);
01156 if(APV->SubDet==StripSubdetector::TOB) Charge_Vs_PathTOB ->Fill(path,Charge);
01157 if(APV->SubDet==StripSubdetector::TEC){ Charge_Vs_PathTEC ->Fill(path,Charge);
01158 if(APV->Thickness<0.04) Charge_Vs_PathTEC1 ->Fill(path,Charge);
01159 if(APV->Thickness>0.04) Charge_Vs_PathTEC2 ->Fill(path,Charge); }
01160
01161 double trans = atan2(trackDirection.y(),trackDirection.x())*(180/3.14159265);
01162 double alpha = acos(trackDirection.x() / sqrt( pow(trackDirection.x(),2) + pow(trackDirection.z(),2) ) ) * (180/3.14159265);
01163 double beta = acos(trackDirection.y() / sqrt( pow(trackDirection.x(),2) + pow(trackDirection.z(),2) ) ) * (180/3.14159265);
01164
01165 if(path>0.4 && path<0.45){
01166 Charge_Vs_TransversAngle->Fill(trans,Charge/path);
01167 Charge_Vs_Alpha->Fill(alpha,Charge/path);
01168 Charge_Vs_Beta->Fill(beta,Charge/path);
01169 }
01170
01171 NStrips_Vs_TransversAngle->Fill(trans,Ampls.size());
01172 NStrips_Vs_Alpha ->Fill(alpha,Ampls.size());
01173 NStrips_Vs_Beta ->Fill(beta ,Ampls.size());
01174
01175 NHighStripInCluster->Fill(NHighStrip);
01176
01177
01178
01179
01180
01181
01182 HFirstStrip ->Fill(FirstStrip);
01183
01184
01185 if(Validation){ClusterChargeOverPath=ClusterChargeOverPath/APV->PreviousGain;}
01186
01187 APV_Charge ->Fill(APV->Index,ClusterChargeOverPath);
01188 APV_Momentum ->Fill(APV->Index,trajState.globalMomentum().mag());
01189 APV_PathLength->Fill(APV->Index,path);
01190
01191 return ClusterChargeOverPath;
01192 }
01193
01194 bool SiStripGainFromData::IsFarFromBorder(TrajectoryStateOnSurface trajState, const uint32_t detid, const edm::EventSetup* iSetup)
01195 {
01196 edm::ESHandle<TrackerGeometry> tkGeom; iSetup->get<TrackerDigiGeometryRecord>().get( tkGeom );
01197
01198 LocalPoint HitLocalPos = trajState.localPosition();
01199 LocalError HitLocalError = trajState.localError().positionError() ;
01200
01201 const GeomDetUnit* it = tkGeom->idToDetUnit(DetId(detid));
01202 if (dynamic_cast<const StripGeomDetUnit*>(it)==0 && dynamic_cast<const PixelGeomDetUnit*>(it)==0) {
01203 std::cout << "this detID doesn't seem to belong to the Tracker" << std::endl;
01204 return false;
01205 }
01206
01207 const BoundPlane plane = it->surface();
01208 const TrapezoidalPlaneBounds* trapezoidalBounds( dynamic_cast<const TrapezoidalPlaneBounds*>(&(plane.bounds())));
01209 const RectangularPlaneBounds* rectangularBounds( dynamic_cast<const RectangularPlaneBounds*>(&(plane.bounds())));
01210
01211 double DistFromBorder = 1.0;
01212 double HalfWidth = it->surface().bounds().width() /2.0;
01213 double HalfLength = it->surface().bounds().length() /2.0;
01214
01215 if(trapezoidalBounds)
01216 {
01217 std::vector<float> const & parameters = (*trapezoidalBounds).parameters();
01218 HalfLength = parameters[3];
01219 double t = (HalfLength + HitLocalPos.y()) / (2*HalfLength) ;
01220 HalfWidth = parameters[0] + (parameters[1]-parameters[0]) * t;
01221 }else if(rectangularBounds){
01222 HalfWidth = it->surface().bounds().width() /2.0;
01223 HalfLength = it->surface().bounds().length() /2.0;
01224 }else{return false;}
01225
01226
01227 if (fabs(HitLocalPos.y())+HitLocalError.yy() >= (HalfLength - DistFromBorder) ) return false;
01228
01229 return true;
01230 }
01231
01232 void SiStripGainFromData::getPeakOfLandau(TH1* InputHisto, double* FitResults, double LowRange, double HighRange)
01233 {
01234 double adcs = -0.5;
01235 double adcs_err = 0.;
01236 double width = -0.5;
01237 double width_err = 0;
01238 double chi2overndf = -0.5;
01239
01240 double nr_of_entries = InputHisto->GetEntries();
01241
01242 if( (unsigned int)nr_of_entries < MinNrEntries){
01243 FitResults[0] = adcs;
01244 FitResults[1] = adcs_err;
01245 FitResults[2] = width;
01246 FitResults[3] = width_err;
01247 FitResults[4] = chi2overndf;
01248 return;
01249 }
01250
01251
01252 TF1* MyLandau = new TF1("MyLandau","landau",LowRange, HighRange);
01253 MyLandau->SetParameter("MPV",300);
01254
01255 InputHisto->Fit("MyLandau","QR WW");
01256 TF1 * fitfunction = (TF1*) InputHisto->GetListOfFunctions()->First();
01257
01258
01259 adcs = fitfunction->GetParameter("MPV");
01260 adcs_err = fitfunction->GetParError(1);
01261 width = fitfunction->GetParameter(2);
01262 width_err = fitfunction->GetParError(2);
01263 chi2overndf = fitfunction->GetChisquare() / fitfunction->GetNDF();
01264
01265
01266 if(adcs<2. || chi2overndf>MaxChi2OverNDF){
01267 adcs = -0.5; adcs_err = 0.;
01268 width = -0.5; width_err = 0;
01269 chi2overndf = -0.5;
01270 }
01271
01272 FitResults[0] = adcs;
01273 FitResults[1] = adcs_err;
01274 FitResults[2] = width;
01275 FitResults[3] = width_err;
01276 FitResults[4] = chi2overndf;
01277
01278 delete MyLandau;
01279 }
01280
01281
01282
01283
01284 SiStripApvGain* SiStripGainFromData::getNewObject()
01285 {
01286 cout << "START getNewObject\n";
01287
01288
01289 if( !(strcmp(AlgoMode.c_str(),"WriteOnDB")==0 || strcmp(AlgoMode.c_str(),"SingleJob")==0) )return new SiStripApvGain();
01290
01291
01292 SiStripApvGain * obj = new SiStripApvGain();
01293 std::vector<float>* theSiStripVector = NULL;
01294 int PreviousDetId = -1;
01295 for(unsigned int a=0;a<APVsCollOrdered.size();a++)
01296 {
01297 stAPVGain* APV = APVsCollOrdered[a];
01298 if(APV==NULL){ printf("Bug\n"); continue; }
01299 if(APV->DetId != PreviousDetId){
01300 if(theSiStripVector!=NULL){
01301 SiStripApvGain::Range range(theSiStripVector->begin(),theSiStripVector->end());
01302 if ( !obj->put(PreviousDetId,range) ) printf("Bug to put detId = %i\n",PreviousDetId);
01303 }
01304
01305 theSiStripVector = new std::vector<float>;
01306 PreviousDetId = APV->DetId;
01307 }
01308 printf("%i | %i | PreviousGain = %7.5f NewGain = %7.5f\n", APV->DetId,APV->APVId,APV->PreviousGain,APV->Gain);
01309 theSiStripVector->push_back(APV->Gain);
01310
01311 }
01312
01313 if(theSiStripVector!=NULL){
01314 SiStripApvGain::Range range(theSiStripVector->begin(),theSiStripVector->end());
01315 if ( !obj->put(PreviousDetId,range) ) printf("Bug to put detId = %i\n",PreviousDetId);
01316 }
01317
01318 cout << "END getNewObject\n";
01319 return obj;
01320 }
01321
01322
01323
01324 DEFINE_FWK_MODULE(SiStripGainFromData);
1.7.3