CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/CalibTracker/SiStripChannelGain/plugins/SiStripGainFromData.cc

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