Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:31:55 

0001 #include "Validation/CaloTowers/interface/CaloTowersValidation.h"
0002 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
0003 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0004 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0005 
0006 CaloTowersValidation::CaloTowersValidation(edm::ParameterSet const& conf) {
0007   tok_calo_ = consumes<CaloTowerCollection>(conf.getUntrackedParameter<edm::InputTag>("CaloTowerCollectionLabel"));
0008   tok_evt_ = consumes<edm::HepMCProduct>(edm::InputTag("generatorSmeared"));
0009 
0010   // DQM ROOT output
0011   outputFile_ = conf.getUntrackedParameter<std::string>("outputFile", "myfile.root");
0012 
0013   hcalselector_ = conf.getUntrackedParameter<std::string>("hcalselector", "all");
0014 
0015   mc_ = conf.getUntrackedParameter<bool>("mc", true);
0016   useAllHistos_ = conf.getUntrackedParameter<bool>("useAllHistos", false);
0017 
0018   etaMin[0] = 0.;
0019   etaMax[0] = 1.4;
0020   etaMin[1] = 1.4;
0021   etaMax[1] = 2.9;
0022   etaMin[2] = 2.9;
0023   etaMax[2] = 5.2;
0024 
0025   isub = 0;
0026   if (hcalselector_ == "HB")
0027     isub = 1;
0028   if (hcalselector_ == "HE")
0029     isub = 2;
0030   if (hcalselector_ == "HF")
0031     isub = 3;
0032 
0033   imc = 1;
0034   if (!mc_)
0035     imc = 0;
0036 
0037   if (!outputFile_.empty()) {
0038     edm::LogInfo("OutputInfo") << " Hcal RecHit Task histograms will be saved to '" << outputFile_.c_str() << "'";
0039   } else {
0040     edm::LogInfo("OutputInfo") << " Hcal RecHit Task histograms will NOT be saved";
0041   }
0042 
0043   nevent = 0;
0044   // const char * sub = hcalselector_.c_str();
0045 }
0046 
0047 CaloTowersValidation::~CaloTowersValidation() {}
0048 
0049 void CaloTowersValidation::bookHistograms(DQMStore::IBooker& ibooker,
0050                                           edm::Run const& irun,
0051                                           edm::EventSetup const& isetup) {
0052   Char_t histo[100];
0053 
0054   ibooker.setCurrentFolder("CaloTowersV/CaloTowersTask");
0055 
0056   //These two histos are not drawn by our macros, but they are used
0057   //in the EndJob for norms and such so I am leaving them alone for now
0058   //-------------------------------------------------------------------------------------------
0059   sprintf(histo, "Ntowers_per_event_vs_ieta");
0060   Ntowers_vs_ieta = ibooker.book1D(histo, histo, 83, -41.5, 41.5);
0061 
0062   sprintf(histo, "CaloTowersTask_map_Nentries");
0063   mapEnergy_N = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0064   //-------------------------------------------------------------------------------------------
0065 
0066   //These the single pion scan histos
0067   //-------------------------------------------------------------------------------------------
0068   sprintf(histo, "emean_vs_ieta_E");
0069   emean_vs_ieta_E = ibooker.bookProfile(histo, histo, 83, -41.5, 41.5, -100., 2000., " ");
0070   sprintf(histo, "emean_vs_ieta_H");
0071   emean_vs_ieta_H = ibooker.bookProfile(histo, histo, 83, -41.5, 41.5, -100., 2000., " ");
0072   sprintf(histo, "emean_vs_ieta_EH");
0073   emean_vs_ieta_EH = ibooker.bookProfile(histo, histo, 83, -41.5, 41.5, -100., 2000., " ");
0074   //These are drawn
0075   sprintf(histo, "emean_vs_ieta_E1");
0076   emean_vs_ieta_E1 = ibooker.bookProfile(histo, histo, 83, -41.5, 41.5, -100., 2000., " ");
0077   sprintf(histo, "emean_vs_ieta_H1");
0078   emean_vs_ieta_H1 = ibooker.bookProfile(histo, histo, 83, -41.5, 41.5, -100., 2000., " ");
0079   sprintf(histo, "emean_vs_ieta_EH1");
0080   emean_vs_ieta_EH1 = ibooker.bookProfile(histo, histo, 83, -41.5, 41.5, -100., 2000., " ");
0081   //-------------------------------------------------------------------------------------------
0082 
0083   //Map energy histos are not drawn
0084   if (useAllHistos_) {
0085     sprintf(histo, "CaloTowersTask_map_energy_E");
0086     mapEnergy_E = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0087     sprintf(histo, "CaloTowersTask_map_energy_H");
0088     mapEnergy_H = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0089     sprintf(histo, "CaloTowersTask_map_energy_EH");
0090     mapEnergy_EH = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0091   }
0092 
0093   //All ECAL cell histos are used
0094   // XXX: ECAL 0-25 [0-26, 26 bins]   HCAL 0-4 [0-5, 5 bins]
0095   sprintf(histo, "number_of_bad_cells_Ecal_EB");
0096   numBadCellsEcal_EB = ibooker.book1D(histo, histo, 26, 0, 26);
0097   sprintf(histo, "number_of_bad_cells_Ecal_EE");
0098   numBadCellsEcal_EE = ibooker.book1D(histo, histo, 26, 0, 26);
0099   sprintf(histo, "number_of_recovered_cells_Ecal_EB");
0100   numRcvCellsEcal_EB = ibooker.book1D(histo, histo, 26, 0, 26);
0101   sprintf(histo, "number_of_recovered_cells_Ecal_EE");
0102   numRcvCellsEcal_EE = ibooker.book1D(histo, histo, 26, 0, 26);
0103   sprintf(histo, "number_of_problematic_cells_Ecal_EB");
0104   numPrbCellsEcal_EB = ibooker.book1D(histo, histo, 26, 0, 26);
0105   sprintf(histo, "number_of_problematic_cells_Ecal_EE");
0106   numPrbCellsEcal_EE = ibooker.book1D(histo, histo, 26, 0, 26);
0107 
0108   //Occupancy vs. ieta is drawn, occupancy map is needed to draw it
0109   sprintf(histo, "CaloTowersTask_map_occupancy");
0110   occupancy_map = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0111 
0112   sprintf(histo, "CaloTowersTask_occupancy_vs_ieta");
0113   occupancy_vs_ieta = ibooker.book1D(histo, histo, 83, -41.5, 41.5);
0114 
0115   if (isub == 1 || isub == 0) {
0116     //All cell histos are used
0117     sprintf(histo, "number_of_bad_cells_Hcal_HB");
0118     numBadCellsHcal_HB = ibooker.book1D(histo, histo, 5, 0, 5);
0119     sprintf(histo, "number_of_recovered_cells_Hcal_HB");
0120     numRcvCellsHcal_HB = ibooker.book1D(histo, histo, 5, 0, 5);
0121     sprintf(histo, "number_of_problematic_cells_Hcal_HB");
0122     numPrbCellsHcal_HB = ibooker.book1D(histo, histo, 5, 0, 5);
0123 
0124     //These are the five oldest CaloTower histos used: NTowers, E in HCAL/ECAL, MET and SET
0125     //-------------------------------------------------------------------------------------------
0126     sprintf(histo, "CaloTowersTask_energy_HCAL_HB");
0127     meEnergyHcal_HB = ibooker.book1D(histo, histo, 4100, -200, 8000);
0128 
0129     sprintf(histo, "CaloTowersTask_energy_ECAL_HB");
0130     meEnergyEcal_HB = ibooker.book1D(histo, histo, 3100, -200, 6000);
0131 
0132     sprintf(histo, "CaloTowersTask_number_of_fired_towers_HB");
0133     meNumFiredTowers_HB = ibooker.book1D(histo, histo, 1000, 0, 2000);
0134 
0135     sprintf(histo, "CaloTowersTask_MET_HB");
0136     MET_HB = ibooker.book1D(histo, histo, 3000, 0., 3000.);
0137 
0138     sprintf(histo, "CaloTowersTask_SET_HB");
0139     SET_HB = ibooker.book1D(histo, histo, 8000, 0., 8000.);
0140     //-------------------------------------------------------------------------------------------
0141 
0142     //Timing histos and profiles -- all six are necessary
0143     //-------------------------------------------------------------------------------------------
0144     sprintf(histo, "CaloTowersTask_EM_Timing_HB");
0145     emTiming_HB = ibooker.book1D(histo, histo, 110, -120., 100.);
0146 
0147     sprintf(histo, "CaloTowersTask_HAD_Timing_HB");
0148     hadTiming_HB = ibooker.book1D(histo, histo, 70, -48., 92.);
0149 
0150     //Energy-Timing histos are divided into low, medium and high to reduce memory usage
0151     //EM
0152     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_Low_HB");
0153     emEnergyTiming_Low_HB = ibooker.book2D(histo, histo, 40, 0., 40., 110, -120., 100.);
0154 
0155     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_HB");
0156     emEnergyTiming_HB = ibooker.book2D(histo, histo, 200, 0., 400., 110, -120., 100.);
0157 
0158     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_High_HB");
0159     emEnergyTiming_High_HB = ibooker.book2D(histo, histo, 200, 0., 3000., 110, -120., 100.);
0160 
0161     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_Low_HB");
0162     emEnergyTiming_profile_Low_HB = ibooker.bookProfile(histo, histo, 40, 0., 40., 110, -120., 100.);
0163 
0164     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_HB");
0165     emEnergyTiming_profile_HB = ibooker.bookProfile(histo, histo, 200, 0., 400., 110, -120., 100.);
0166 
0167     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_High_HB");
0168     emEnergyTiming_profile_High_HB = ibooker.bookProfile(histo, histo, 200, 0., 3000., 110, -120., 100.);
0169 
0170     //HAD
0171     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_Low_HB");
0172     hadEnergyTiming_Low_HB = ibooker.book2D(histo, histo, 40, 0., 40., 70, -48., 92.);
0173 
0174     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_HB");
0175     hadEnergyTiming_HB = ibooker.book2D(histo, histo, 100, 0., 200., 70, -48., 92.);
0176 
0177     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_High_HB");
0178     hadEnergyTiming_High_HB = ibooker.book2D(histo, histo, 300, 0., 3000., 70, -48., 92.);
0179 
0180     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_Low_HB");
0181     hadEnergyTiming_profile_Low_HB = ibooker.bookProfile(histo, histo, 40, 0., 40., 70, -48., 92.);
0182 
0183     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_HB");
0184     hadEnergyTiming_profile_HB = ibooker.bookProfile(histo, histo, 100, 0., 200., 70, -48., 92.);
0185 
0186     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_High_HB");
0187     hadEnergyTiming_profile_High_HB = ibooker.bookProfile(histo, histo, 300, 0., 3000., 70, -48., 92.);
0188     //-------------------------------------------------------------------------------------------
0189 
0190     //Everything else is not drawn
0191     if (useAllHistos_) {
0192       sprintf(histo, "CaloTowersTask_sum_of_energy_HCAL_vs_ECAL_HB");
0193       meEnergyHcalvsEcal_HB = ibooker.book2D(histo, histo, 500, 0., 500., 500, 0., 500.);
0194 
0195       sprintf(histo, "CaloTowersTask_energy_OUTER_HB");
0196       meEnergyHO_HB = ibooker.book1D(histo, histo, 1640, -200, 8000);
0197 
0198       sprintf(histo, "CaloTowersTask_energy_of_ECAL_component_of_tower_HB");
0199       meEnergyEcalTower_HB = ibooker.book1D(histo, histo, 440, -200, 2000);
0200 
0201       sprintf(histo, "CaloTowersTask_energy_of_HCAL_component_of_tower_HB");
0202       meEnergyHcalTower_HB = ibooker.book1D(histo, histo, 440, -200, 2000);
0203 
0204       sprintf(histo, "CaloTowersTask_energy_HcalPlusEcalPlusHO_HB");
0205       meTotEnergy_HB = ibooker.book1D(histo, histo, 400, 0., 2000.);
0206 
0207       sprintf(histo, "CaloTowersTask_map_energy_HB");
0208       mapEnergy_HB = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0209       sprintf(histo, "CaloTowersTask_map_energy_HCAL_HB");
0210       mapEnergyHcal_HB = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0211       sprintf(histo, "CaloTowersTask_map_energy_ECAL_HB");
0212       mapEnergyEcal_HB = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0213 
0214       sprintf(histo, "CaloTowersTask_phi_MET_HB");
0215       phiMET_HB = ibooker.book1D(histo, histo, 72, -3.1415926535898, 3.1415926535898);
0216     }
0217   }
0218 
0219   if (isub == 2 || isub == 0) {
0220     //All cell histos are used
0221     sprintf(histo, "number_of_bad_cells_Hcal_HE");
0222     numBadCellsHcal_HE = ibooker.book1D(histo, histo, 5, 0, 5);
0223     sprintf(histo, "number_of_recovered_cells_Hcal_HE");
0224     numRcvCellsHcal_HE = ibooker.book1D(histo, histo, 5, 0, 5);
0225     sprintf(histo, "number_of_problematic_cells_Hcal_HE");
0226     numPrbCellsHcal_HE = ibooker.book1D(histo, histo, 5, 0, 5);
0227 
0228     //These are the five oldest CaloTower histos used: NTowers, E in HCAL/ECAL, MET and SET
0229     //-------------------------------------------------------------------------------------------
0230     sprintf(histo, "CaloTowersTask_energy_HCAL_HE");
0231     meEnergyHcal_HE = ibooker.book1D(histo, histo, 1240, -200, 6000);
0232 
0233     sprintf(histo, "CaloTowersTask_energy_ECAL_HE");
0234     meEnergyEcal_HE = ibooker.book1D(histo, histo, 1240, -200, 6000);
0235 
0236     sprintf(histo, "CaloTowersTask_number_of_fired_towers_HE");
0237     meNumFiredTowers_HE = ibooker.book1D(histo, histo, 1000, 0, 2000);
0238 
0239     sprintf(histo, "CaloTowersTask_MET_HE");
0240     MET_HE = ibooker.book1D(histo, histo, 1000, 0., 1000.);
0241 
0242     sprintf(histo, "CaloTowersTask_SET_HE");
0243     SET_HE = ibooker.book1D(histo, histo, 2000, 0., 2000.);
0244     //-------------------------------------------------------------------------------------------
0245 
0246     //Timing histos and profiles -- all six are necessary
0247     //-------------------------------------------------------------------------------------------
0248     sprintf(histo, "CaloTowersTask_EM_Timing_HE");
0249     emTiming_HE = ibooker.book1D(histo, histo, 110, -120., 100.);
0250 
0251     sprintf(histo, "CaloTowersTask_HAD_Timing_HE");
0252     hadTiming_HE = ibooker.book1D(histo, histo, 70, -48., 92.);
0253 
0254     //Energy-Timing histos are divided into low and normal to reduce memory usage
0255     //EM
0256     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_Low_HE");
0257     emEnergyTiming_Low_HE = ibooker.book2D(histo, histo, 160, 0., 160., 110, -120., 100.);
0258 
0259     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_HE");
0260     emEnergyTiming_HE = ibooker.book2D(histo, histo, 200, 0., 800., 110, -120., 100.);
0261 
0262     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_Low_HE");
0263     emEnergyTiming_profile_Low_HE = ibooker.bookProfile(histo, histo, 160, 0., 160., 110, -120., 100.);
0264 
0265     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_HE");
0266     emEnergyTiming_profile_HE = ibooker.bookProfile(histo, histo, 200, 0., 800., 110, -120., 100.);
0267 
0268     //HAD
0269     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_Low_HE");
0270     hadEnergyTiming_Low_HE = ibooker.book2D(histo, histo, 160, 0., 160., 70, -48., 92.);
0271 
0272     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_HE");
0273     hadEnergyTiming_HE = ibooker.book2D(histo, histo, 200, 0., 800., 70, -48., 92.);
0274 
0275     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_Low_HE");
0276     hadEnergyTiming_profile_Low_HE = ibooker.bookProfile(histo, histo, 160, 0., 160., 70, -48., 92.);
0277 
0278     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_HE");
0279     hadEnergyTiming_profile_HE = ibooker.bookProfile(histo, histo, 200, 0., 800., 70, -48., 92.);
0280     //-------------------------------------------------------------------------------------------
0281 
0282     //Everything else is not drawn
0283     if (useAllHistos_) {
0284       sprintf(histo, "CaloTowersTask_sum_of_energy_HCAL_vs_ECAL_HE");
0285       meEnergyHcalvsEcal_HE = ibooker.book2D(histo, histo, 500, 0., 500., 500, 0., 500.);
0286 
0287       sprintf(histo, "CaloTowersTask_energy_OUTER_HE");
0288       meEnergyHO_HE = ibooker.book1D(histo, histo, 440, -200, 2000);
0289 
0290       sprintf(histo, "CaloTowersTask_energy_of_ECAL_component_of_tower_HE");
0291       meEnergyEcalTower_HE = ibooker.book1D(histo, histo, 1100, -200, 2000);
0292 
0293       sprintf(histo, "CaloTowersTask_energy_of_HCAL_component_of_tower_HE");
0294       meEnergyHcalTower_HE = ibooker.book1D(histo, histo, 1100, -200, 2000);
0295 
0296       sprintf(histo, "CaloTowersTask_energy_HcalPlusEcalPlusHO_HE");
0297       meTotEnergy_HE = ibooker.book1D(histo, histo, 400, 0., 2000.);
0298 
0299       sprintf(histo, "CaloTowersTask_map_energy_HE");
0300       mapEnergy_HE = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0301       sprintf(histo, "CaloTowersTask_map_energy_HCAL_HE");
0302       mapEnergyHcal_HE = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0303       sprintf(histo, "CaloTowersTask_map_energy_ECAL_HE");
0304       mapEnergyEcal_HE = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0305 
0306       sprintf(histo, "CaloTowersTask_phi_MET_HE");
0307       phiMET_HE = ibooker.book1D(histo, histo, 72, -3.1415926535898, 3.1415926535898);
0308     }
0309   }
0310 
0311   if (isub == 3 || isub == 0) {
0312     //All cell histos are used
0313     sprintf(histo, "number_of_bad_cells_Hcal_HF");
0314     numBadCellsHcal_HF = ibooker.book1D(histo, histo, 5, 0, 5);
0315     sprintf(histo, "number_of_recovered_cells_Hcal_HF");
0316     numRcvCellsHcal_HF = ibooker.book1D(histo, histo, 5, 0, 5);
0317     sprintf(histo, "number_of_problematic_cells_Hcal_HF");
0318     numPrbCellsHcal_HF = ibooker.book1D(histo, histo, 5, 0, 5);
0319 
0320     //These are the five oldest CaloTower histos used: NTowers, E in HCAL/ECAL, MET and SET
0321     //-------------------------------------------------------------------------------------------
0322     sprintf(histo, "CaloTowersTask_energy_HCAL_HF");
0323     meEnergyHcal_HF = ibooker.book1D(histo, histo, 4040, -200, 20000);
0324 
0325     sprintf(histo, "CaloTowersTask_energy_ECAL_HF");
0326     meEnergyEcal_HF = ibooker.book1D(histo, histo, 2440, -200, 12000);
0327 
0328     sprintf(histo, "CaloTowersTask_number_of_fired_towers_HF");
0329     meNumFiredTowers_HF = ibooker.book1D(histo, histo, 1000, 0, 2000);
0330 
0331     sprintf(histo, "CaloTowersTask_MET_HF");
0332     MET_HF = ibooker.book1D(histo, histo, 500, 0., 500.);
0333 
0334     sprintf(histo, "CaloTowersTask_SET_HF");
0335     SET_HF = ibooker.book1D(histo, histo, 2000, 0., 2000.);
0336     //-------------------------------------------------------------------------------------------
0337 
0338     //Timing histos and profiles -- all six are necessary
0339     //-------------------------------------------------------------------------------------------
0340     sprintf(histo, "CaloTowersTask_EM_Timing_HF");
0341     emTiming_HF = ibooker.book1D(histo, histo, 110, -120., 100.);
0342 
0343     sprintf(histo, "CaloTowersTask_HAD_Timing_HF");
0344     hadTiming_HF = ibooker.book1D(histo, histo, 70, -48., 92.);
0345 
0346     //EM
0347     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_HF");
0348     emEnergyTiming_HF = ibooker.book2D(histo, histo, 150, 0., 300., 110, -120., 100.);
0349 
0350     sprintf(histo, "CaloTowersTask_EM_Energy_Timing_profile_HF");
0351     emEnergyTiming_profile_HF = ibooker.bookProfile(histo, histo, 150, 0., 300., 110, -120., 100.);
0352 
0353     //HAD (requires two different sets of histograms to lower RAM usage)
0354     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_Low_HF");
0355     hadEnergyTiming_Low_HF = ibooker.book2D(histo, histo, 40, 0., 40., 70, -48., 92.);
0356 
0357     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_HF");
0358     hadEnergyTiming_HF = ibooker.book2D(histo, histo, 200, 0., 600., 70, -48., 92.);
0359 
0360     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_Low_HF");
0361     hadEnergyTiming_profile_Low_HF = ibooker.bookProfile(histo, histo, 40, 0., 40., 70, -48., 92.);
0362 
0363     sprintf(histo, "CaloTowersTask_HAD_Energy_Timing_profile_HF");
0364     hadEnergyTiming_profile_HF = ibooker.bookProfile(histo, histo, 200, 0., 600., 70, -48., 92.);
0365     //-------------------------------------------------------------------------------------------
0366 
0367     //Everything else is not drawn
0368     if (useAllHistos_) {
0369       sprintf(histo, "CaloTowersTask_sum_of_energy_HCAL_vs_ECAL_HF");
0370       meEnergyHcalvsEcal_HF = ibooker.book2D(histo, histo, 500, 0., 500., 500, 0., 500.);
0371 
0372       sprintf(histo, "CaloTowersTask_energy_OUTER_HF");
0373       meEnergyHO_HF = ibooker.book1D(histo, histo, 440, -200, 2000);
0374 
0375       sprintf(histo, "CaloTowersTask_energy_of_ECAL_component_of_tower_HF");
0376       meEnergyEcalTower_HF = ibooker.book1D(histo, histo, 440, -200, 2000);
0377 
0378       sprintf(histo, "CaloTowersTask_energy_of_HCAL_component_of_tower_HF");
0379       meEnergyHcalTower_HF = ibooker.book1D(histo, histo, 440, -200, 2000);
0380 
0381       sprintf(histo, "CaloTowersTask_energy_HcalPlusEcalPlusHO_HF");
0382       meTotEnergy_HF = ibooker.book1D(histo, histo, 400, 0., 2000.);
0383 
0384       sprintf(histo, "CaloTowersTask_map_energy_HF");
0385       mapEnergy_HF = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0386       sprintf(histo, "CaloTowersTask_map_energy_HCAL_HF");
0387       mapEnergyHcal_HF = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0388       sprintf(histo, "CaloTowersTask_map_energy_ECAL_HF");
0389       mapEnergyEcal_HF = ibooker.book2D(histo, histo, 83, -41.5, 41.5, 72, 0., 72.);
0390 
0391       sprintf(histo, "CaloTowersTask_phi_MET_HF");
0392       phiMET_HF = ibooker.book1D(histo, histo, 72, -3.1415926535898, 3.1415926535898);
0393     }
0394   }
0395 }
0396 
0397 void CaloTowersValidation::analyze(edm::Event const& event, edm::EventSetup const& c) {
0398   nevent++;
0399 
0400   // bool     MC = false; // UNUSED
0401   double phi_MC = 9999.;
0402   double eta_MC = 9999.;
0403 
0404   if (imc != 0) {
0405     edm::Handle<edm::HepMCProduct> evtMC;
0406     event.getByToken(tok_evt_, evtMC);  // generator in late 310_preX
0407     if (!evtMC.isValid()) {
0408       std::cout << "no HepMCProduct found" << std::endl;
0409     } else {
0410       // MC=true; // UNUSED
0411       //    std::cout << "*** source HepMCProduct found"<< std::endl;
0412     }
0413 
0414     // MC particle with highest pt is taken as a direction reference
0415     double maxPt = -99999.;
0416     const HepMC::GenEvent* myGenEvent = evtMC->GetEvent();
0417     for (HepMC::GenEvent::particle_const_iterator p = myGenEvent->particles_begin(); p != myGenEvent->particles_end();
0418          ++p) {
0419       double phip = (*p)->momentum().phi();
0420       double etap = (*p)->momentum().eta();
0421       //    phi_MC = phip;
0422       //    eta_MC = etap;
0423       double pt = (*p)->momentum().perp();
0424       if (pt > maxPt) {
0425         maxPt = pt;
0426         phi_MC = phip;
0427         eta_MC = etap;
0428       }
0429     }
0430     //  std::cout << "*** Max pT = " << maxPt <<  std::endl;
0431   }
0432 
0433   edm::Handle<CaloTowerCollection> towers;
0434   event.getByToken(tok_calo_, towers);
0435   CaloTowerCollection::const_iterator cal;
0436 
0437   double met;
0438   double phimet;
0439 
0440   // ieta scan
0441   double partR = 0.3;
0442   double Rmin = 9999.;
0443   double Econe = 0.;
0444   double Hcone = 0.;
0445   double Ee1 = 0.;
0446   double Eh1 = 0.;
0447   double ieta_MC = 9999;
0448   double iphi_MC = 9999;
0449   //  double  etaM   = 9999.;
0450 
0451   // HB
0452   double sumEnergyHcal_HB = 0.;
0453   double sumEnergyEcal_HB = 0.;
0454   double sumEnergyHO_HB = 0.;
0455   Int_t numFiredTowers_HB = 0;
0456   double metx_HB = 0.;
0457   double mety_HB = 0.;
0458   double metz_HB = 0.;
0459   double sEt_HB = 0.;
0460   // HE
0461   double sumEnergyHcal_HE = 0.;
0462   double sumEnergyEcal_HE = 0.;
0463   double sumEnergyHO_HE = 0.;
0464   Int_t numFiredTowers_HE = 0;
0465   double metx_HE = 0.;
0466   double mety_HE = 0.;
0467   double metz_HE = 0.;
0468   double sEt_HE = 0.;
0469   // HF
0470   double sumEnergyHcal_HF = 0.;
0471   double sumEnergyEcal_HF = 0.;
0472   double sumEnergyHO_HF = 0.;
0473   Int_t numFiredTowers_HF = 0;
0474   double metx_HF = 0.;
0475   double mety_HF = 0.;
0476   double metz_HF = 0.;
0477   double sEt_HF = 0.;
0478 
0479   for (cal = towers->begin(); cal != towers->end(); ++cal) {
0480     double eE = cal->emEnergy();
0481     double eH = cal->hadEnergy();
0482     double eHO = cal->outerEnergy();
0483     double etaT = cal->eta();
0484     double phiT = cal->phi();
0485     double en = cal->energy();
0486     double etT = cal->et();
0487     double had_tm = cal->hcalTime();
0488     double em_tm = cal->ecalTime();
0489 
0490     int numBadEcalCells = cal->numBadEcalCells();
0491     int numRcvEcalCells = cal->numRecoveredEcalCells();
0492     int numPrbEcalCells = cal->numProblematicEcalCells();
0493 
0494     int numBadHcalCells = cal->numBadHcalCells();
0495     int numRcvHcalCells = cal->numRecoveredHcalCells();
0496     int numPrbHcalCells = cal->numProblematicHcalCells();
0497 
0498     math::RhoEtaPhiVector mom(cal->et(), cal->eta(), cal->phi());
0499     //  Vector mom  = cal->momentum();
0500 
0501     // cell properties
0502     CaloTowerDetId idT = cal->id();
0503     int ieta = idT.ieta();
0504     int iphi = idT.iphi();
0505 
0506     // ecal:  0 EcalBarrel  1 EcalEndcap
0507     // hcal:  0 hcalBarrel  1 HcalEndcap  2 HcalForward
0508     std::vector<int> inEcals(2), inHcals(3);
0509     unsigned int constitSize = cal->constituentsSize();
0510     for (unsigned int ic = 0; ic < constitSize; ic++) {
0511       DetId detId = cal->constituent(ic);
0512       if (detId.det() == DetId::Ecal) {
0513         if (detId.subdetId() == EcalBarrel)
0514           inEcals[0] = 1;
0515         else if (detId.subdetId() == EcalEndcap)
0516           inEcals[1] = 1;
0517       }
0518       if (detId.det() == DetId::Hcal) {
0519         if (HcalDetId(detId).subdet() == HcalBarrel)
0520           inHcals[0] = 1;
0521         else if (HcalDetId(detId).subdet() == HcalEndcap)
0522           inHcals[1] = 1;
0523         else if (HcalDetId(detId).subdet() == HcalForward)
0524           inHcals[2] = 1;
0525       }
0526     }
0527     //All cell histos are used
0528     if (inEcals[0]) {
0529       numBadCellsEcal_EB->Fill(numBadEcalCells);
0530       numRcvCellsEcal_EB->Fill(numRcvEcalCells);
0531       numPrbCellsEcal_EB->Fill(numPrbEcalCells);
0532     }
0533     if (inEcals[1]) {
0534       numBadCellsEcal_EE->Fill(numBadEcalCells);
0535       numRcvCellsEcal_EE->Fill(numRcvEcalCells);
0536       numPrbCellsEcal_EE->Fill(numPrbEcalCells);
0537     }
0538 
0539     if (imc != 0) {
0540       double r = dR(eta_MC, phi_MC, etaT, phiT);
0541 
0542       if (r < partR) {
0543         Econe += eE;
0544         Hcone += eH;
0545 
0546         // closest to MC
0547         if (r < Rmin) {
0548           if (fabs(eta_MC) < 3.0 && (ieta > 29 || ieta < -29)) {
0549             ;
0550           } else {
0551             Rmin = r;
0552             ieta_MC = ieta;
0553             iphi_MC = iphi;
0554             Ee1 = eE;
0555             Eh1 = eH;
0556           }
0557         }
0558       }
0559     }
0560 
0561     //Ntowers is used in EndJob, occupancy_map is used for occupancy vs ieta
0562     Ntowers_vs_ieta->Fill(double(ieta), 1.);
0563     occupancy_map->Fill(double(ieta), double(iphi));
0564 
0565     if ((isub == 0 || isub == 1) && (fabs(etaT) < etaMax[0] && fabs(etaT) >= etaMin[0])) {
0566       //All cell histos are used
0567       numBadCellsHcal_HB->Fill(numBadHcalCells);
0568       numRcvCellsHcal_HB->Fill(numRcvHcalCells);
0569       numPrbCellsHcal_HB->Fill(numPrbHcalCells);
0570 
0571       //Map energy histos are not used
0572       if (useAllHistos_) {
0573         mapEnergy_HB->Fill(double(ieta), double(iphi), en);
0574         mapEnergyHcal_HB->Fill(double(ieta), double(iphi), eH);
0575         mapEnergyEcal_HB->Fill(double(ieta), double(iphi), eE);
0576       }
0577       //      std::cout << " e_ecal = " << eE << std::endl;
0578 
0579       //  simple sums
0580       sumEnergyHcal_HB += eH;
0581       sumEnergyEcal_HB += eE;
0582       sumEnergyHO_HB += eHO;
0583 
0584       numFiredTowers_HB++;
0585 
0586       //Not used
0587       if (useAllHistos_) {
0588         meEnergyEcalTower_HB->Fill(eE);
0589         meEnergyHcalTower_HB->Fill(eH);
0590       }
0591 
0592       // MET, SET & phimet
0593       //  double  etT = cal->et();
0594       metx_HB += mom.x();
0595       mety_HB += mom.y();  //etT * sin(phiT);
0596       sEt_HB += etT;
0597 
0598       //Timing (all histos are used)
0599       emTiming_HB->Fill(em_tm);
0600       hadTiming_HB->Fill(had_tm);
0601 
0602       emEnergyTiming_Low_HB->Fill(eE, em_tm);
0603       emEnergyTiming_HB->Fill(eE, em_tm);
0604       emEnergyTiming_High_HB->Fill(eE, em_tm);
0605       emEnergyTiming_profile_Low_HB->Fill(eE, em_tm);
0606       emEnergyTiming_profile_HB->Fill(eE, em_tm);
0607       emEnergyTiming_profile_High_HB->Fill(eE, em_tm);
0608 
0609       hadEnergyTiming_Low_HB->Fill(eH, had_tm);
0610       hadEnergyTiming_HB->Fill(eH, had_tm);
0611       hadEnergyTiming_High_HB->Fill(eH, had_tm);
0612       hadEnergyTiming_profile_Low_HB->Fill(eH, had_tm);
0613       hadEnergyTiming_profile_HB->Fill(eH, had_tm);
0614       hadEnergyTiming_profile_High_HB->Fill(eH, had_tm);
0615     }
0616 
0617     if ((isub == 0 || isub == 2) && (fabs(etaT) < etaMax[1] && fabs(etaT) >= etaMin[1])) {
0618       //All cell histos are used
0619       numBadCellsHcal_HE->Fill(numBadHcalCells);
0620       numRcvCellsHcal_HE->Fill(numRcvHcalCells);
0621       numPrbCellsHcal_HE->Fill(numPrbHcalCells);
0622 
0623       //Map energy histos are not used
0624       if (useAllHistos_) {
0625         mapEnergy_HE->Fill(double(ieta), double(iphi), en);
0626         mapEnergyHcal_HE->Fill(double(ieta), double(iphi), eH);
0627         mapEnergyEcal_HE->Fill(double(ieta), double(iphi), eE);
0628       }
0629       //      std::cout << " e_ecal = " << eE << std::endl;
0630 
0631       //  simple sums
0632       sumEnergyHcal_HE += eH;
0633       sumEnergyEcal_HE += eE;
0634       sumEnergyHO_HE += eHO;
0635 
0636       numFiredTowers_HE++;
0637 
0638       //Not used
0639       if (useAllHistos_) {
0640         meEnergyEcalTower_HE->Fill(eE);
0641         meEnergyHcalTower_HE->Fill(eH);
0642       }
0643       // MET, SET & phimet
0644       //  double  etT = cal->et();
0645       metx_HE += mom.x();
0646       mety_HE += mom.y();  //etT * sin(phiT);
0647       sEt_HE += etT;
0648 
0649       //Timing (all histos are used)
0650       emTiming_HE->Fill(em_tm);
0651       hadTiming_HE->Fill(had_tm);
0652 
0653       emEnergyTiming_Low_HE->Fill(eE, em_tm);
0654       emEnergyTiming_HE->Fill(eE, em_tm);
0655       emEnergyTiming_profile_Low_HE->Fill(eE, em_tm);
0656       emEnergyTiming_profile_HE->Fill(eE, em_tm);
0657 
0658       hadEnergyTiming_Low_HE->Fill(eH, had_tm);
0659       hadEnergyTiming_HE->Fill(eH, had_tm);
0660       hadEnergyTiming_profile_Low_HE->Fill(eH, had_tm);
0661       hadEnergyTiming_profile_HE->Fill(eH, had_tm);
0662     }
0663 
0664     if ((isub == 0 || isub == 3) && (fabs(etaT) < etaMax[2] && fabs(etaT) >= etaMin[2])) {
0665       //All cell histos are used
0666       numBadCellsHcal_HF->Fill(numBadHcalCells);
0667       numRcvCellsHcal_HF->Fill(numRcvHcalCells);
0668       numPrbCellsHcal_HF->Fill(numPrbHcalCells);
0669 
0670       //Map energy histos are not used
0671       if (useAllHistos_) {
0672         mapEnergy_HF->Fill(double(ieta), double(iphi), en);
0673         mapEnergyHcal_HF->Fill(double(ieta), double(iphi), eH);
0674         mapEnergyEcal_HF->Fill(double(ieta), double(iphi), eE);
0675       }
0676       //      std::cout << " e_ecal = " << eE << std::endl;
0677 
0678       //  simple sums
0679       sumEnergyHcal_HF += eH;
0680       sumEnergyEcal_HF += eE;
0681       sumEnergyHO_HF += eHO;
0682 
0683       numFiredTowers_HF++;
0684 
0685       //Not used
0686       if (useAllHistos_) {
0687         meEnergyEcalTower_HF->Fill(eE);
0688         meEnergyHcalTower_HF->Fill(eH);
0689       }
0690       // MET, SET & phimet
0691       //  double  etT = cal->et();
0692       metx_HF += mom.x();
0693       mety_HF += mom.y();  //etT * sin(phiT);
0694       sEt_HF += etT;
0695 
0696       //Timing (all histos are used)
0697       emTiming_HF->Fill(em_tm);
0698       hadTiming_HF->Fill(had_tm);
0699       emEnergyTiming_HF->Fill(eE, em_tm);
0700       emEnergyTiming_profile_HF->Fill(eE, em_tm);
0701 
0702       hadEnergyTiming_Low_HF->Fill(eH, had_tm);
0703       hadEnergyTiming_HF->Fill(eH, had_tm);
0704       hadEnergyTiming_profile_Low_HF->Fill(eH, had_tm);
0705       hadEnergyTiming_profile_HF->Fill(eH, had_tm);
0706     }
0707 
0708   }  // end of Towers cycle
0709 
0710   //These are the six single pion histos
0711   emean_vs_ieta_E->Fill(double(ieta_MC), Econe);
0712   emean_vs_ieta_H->Fill(double(ieta_MC), Hcone);
0713   emean_vs_ieta_EH->Fill(double(ieta_MC), Econe + Hcone);
0714   emean_vs_ieta_E1->Fill(double(ieta_MC), Ee1);
0715   emean_vs_ieta_H1->Fill(double(ieta_MC), Eh1);
0716   emean_vs_ieta_EH1->Fill(double(ieta_MC), Ee1 + Eh1);
0717 
0718   //Map histos are not used except the last one in EndJob
0719   if (useAllHistos_) {
0720     mapEnergy_E->Fill(double(ieta_MC), double(iphi_MC), Ee1);
0721     mapEnergy_H->Fill(double(ieta_MC), double(iphi_MC), Eh1);
0722     mapEnergy_EH->Fill(double(ieta_MC), double(iphi_MC), Ee1 + Eh1);
0723   }
0724   mapEnergy_N->Fill(double(ieta_MC), double(iphi_MC), 1.);
0725 
0726   if (isub == 0 || isub == 1) {
0727     met = sqrt(metx_HB * metx_HB + mety_HB * mety_HB);
0728     Vector metv(metx_HB, mety_HB, metz_HB);
0729     phimet = metv.phi();
0730 
0731     //Five oldest drawn histos first; the rest are not used
0732     meEnergyHcal_HB->Fill(sumEnergyHcal_HB);
0733     meEnergyEcal_HB->Fill(sumEnergyEcal_HB);
0734     meNumFiredTowers_HB->Fill(numFiredTowers_HB);
0735     MET_HB->Fill(met);
0736     SET_HB->Fill(sEt_HB);
0737 
0738     if (useAllHistos_) {
0739       meEnergyHcalvsEcal_HB->Fill(sumEnergyEcal_HB, sumEnergyHcal_HB);
0740       meEnergyHO_HB->Fill(sumEnergyHO_HB);
0741       meTotEnergy_HB->Fill(sumEnergyHcal_HB + sumEnergyEcal_HB + sumEnergyHO_HB);
0742       phiMET_HB->Fill(phimet);
0743     }
0744   }
0745 
0746   if (isub == 0 || isub == 2) {
0747     met = sqrt(metx_HE * metx_HE + mety_HE * mety_HE);
0748     Vector metv(metx_HE, mety_HE, metz_HE);
0749     phimet = metv.phi();
0750 
0751     //Five oldest drawn histos first; the rest are not used
0752     meEnergyHcal_HE->Fill(sumEnergyHcal_HE);
0753     meEnergyEcal_HE->Fill(sumEnergyEcal_HE);
0754     meNumFiredTowers_HE->Fill(numFiredTowers_HE);
0755     MET_HE->Fill(met);
0756     SET_HE->Fill(sEt_HE);
0757 
0758     if (useAllHistos_) {
0759       meEnergyHcalvsEcal_HE->Fill(sumEnergyEcal_HE, sumEnergyHcal_HE);
0760       meEnergyHO_HE->Fill(sumEnergyHO_HE);
0761       meTotEnergy_HE->Fill(sumEnergyHcal_HE + sumEnergyEcal_HE + sumEnergyHO_HE);
0762       phiMET_HE->Fill(phimet);
0763     }
0764   }
0765 
0766   if (isub == 0 || isub == 3) {
0767     met = sqrt(metx_HF * metx_HF + mety_HF * mety_HF);
0768     Vector metv(metx_HF, mety_HF, metz_HF);
0769     phimet = metv.phi();
0770 
0771     //Five oldest drawn histos first; the rest are not used
0772     meEnergyHcal_HF->Fill(sumEnergyHcal_HF);
0773     meEnergyEcal_HF->Fill(sumEnergyEcal_HF);
0774     meNumFiredTowers_HF->Fill(numFiredTowers_HF);
0775     MET_HF->Fill(met);
0776     SET_HF->Fill(sEt_HF);
0777 
0778     if (useAllHistos_) {
0779       meEnergyHcalvsEcal_HF->Fill(sumEnergyEcal_HF, sumEnergyHcal_HF);
0780       meEnergyHO_HF->Fill(sumEnergyHO_HF);
0781       meTotEnergy_HF->Fill(sumEnergyHcal_HF + sumEnergyEcal_HF + sumEnergyHO_HF);
0782       phiMET_HF->Fill(phimet);
0783     }
0784   }
0785 }
0786 
0787 double CaloTowersValidation::dR(double eta1, double phi1, double eta2, double phi2) {
0788   double PI = 3.1415926535898;
0789   double deltaphi = phi1 - phi2;
0790   if (phi2 > phi1) {
0791     deltaphi = phi2 - phi1;
0792   }
0793   if (deltaphi > PI) {
0794     deltaphi = 2. * PI - deltaphi;
0795   }
0796   double deltaeta = eta2 - eta1;
0797   double tmp = sqrt(deltaeta * deltaeta + deltaphi * deltaphi);
0798   return tmp;
0799 }
0800 
0801 DEFINE_FWK_MODULE(CaloTowersValidation);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!