Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-08-10 22:50:01

0001 // user include files
0002 
0003 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0004 #include "FWCore/Framework/interface/Frameworkfwd.h"
0005 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0006 
0007 #include "FWCore/Framework/interface/ESHandle.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/MakerMacros.h"
0010 
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 
0013 #include "DQMServices/Core/interface/DQMStore.h"
0014 
0015 #include "DataFormats/ForwardDetId/interface/HGCalDetId.h"
0016 #include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
0017 
0018 #include "DataFormats/CaloRecHit/interface/CaloClusterFwd.h"
0019 #include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
0020 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0021 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0022 #include "DataFormats/EgammaCandidates/interface/GsfElectron.h"
0023 #include "DataFormats/EgammaCandidates/interface/Photon.h"
0024 
0025 #include "SimDataFormats/CaloAnalysis/interface/CaloParticle.h"
0026 #include "SimDataFormats/CaloAnalysis/interface/SimCluster.h"
0027 #include "SimDataFormats/TrackingAnalysis/interface/TrackingParticle.h"
0028 #include "SimDataFormats/TrackingAnalysis/interface/TrackingVertex.h"
0029 
0030 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0031 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0032 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0033 #include "Geometry/CaloGeometry/interface/TruncatedPyramid.h"
0034 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0035 
0036 #include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"
0037 
0038 #include <map>
0039 #include <array>
0040 #include <string>
0041 #include <numeric>
0042 
0043 class HGCalHitCalibration : public DQMEDAnalyzer {
0044 public:
0045   explicit HGCalHitCalibration(const edm::ParameterSet&);
0046   ~HGCalHitCalibration() override;
0047 
0048   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0049 
0050 private:
0051   void bookHistograms(DQMStore::IBooker&, edm::Run const&, edm::EventSetup const&) override;
0052   void analyze(const edm::Event&, const edm::EventSetup&) override;
0053 
0054   void fillWithRecHits(std::map<DetId, const HGCRecHit*>&, DetId, unsigned int, float, int&, float&);
0055 
0056   edm::EDGetTokenT<HGCRecHitCollection> recHitsEE_;
0057   edm::EDGetTokenT<HGCRecHitCollection> recHitsFH_;
0058   edm::EDGetTokenT<HGCRecHitCollection> recHitsBH_;
0059   edm::EDGetTokenT<std::vector<CaloParticle> > caloParticles_;
0060   edm::EDGetTokenT<std::vector<reco::PFCluster> > hgcalMultiClusters_;
0061   edm::EDGetTokenT<std::vector<reco::GsfElectron> > electrons_;
0062   edm::EDGetTokenT<std::vector<reco::Photon> > photons_;
0063 
0064   int algo_, depletion0_;
0065   edm::ESGetToken<CaloGeometry, CaloGeometryRecord> caloGeomToken_;
0066   bool rawRecHits_;
0067   int debug_;
0068   hgcal::RecHitTools recHitTools_;
0069   static constexpr int depletion1_ = 200;
0070   static constexpr int depletion2_ = 300;
0071   static constexpr int scint_ = 400;
0072 
0073   std::map<int, MonitorElement*> h_EoP_CPene_calib_fraction_;
0074   std::map<int, MonitorElement*> hgcal_EoP_CPene_calib_fraction_;
0075   std::map<int, MonitorElement*> hgcal_ele_EoP_CPene_calib_fraction_;
0076   std::map<int, MonitorElement*> hgcal_photon_EoP_CPene_calib_fraction_;
0077   MonitorElement* LayerOccupancy_;
0078 
0079   static constexpr int layers_ = 60;
0080   std::array<float, layers_> Energy_layer_calib_;
0081   std::array<float, layers_> Energy_layer_calib_fraction_;
0082 };
0083 
0084 HGCalHitCalibration::HGCalHitCalibration(const edm::ParameterSet& iConfig)
0085     : caloGeomToken_(esConsumes<CaloGeometry, CaloGeometryRecord>()),
0086       rawRecHits_(iConfig.getParameter<bool>("rawRecHits")),
0087       debug_(iConfig.getParameter<int>("debug")) {
0088   auto detector = iConfig.getParameter<std::string>("detector");
0089   auto recHitsEE = iConfig.getParameter<edm::InputTag>("recHitsEE");
0090   auto recHitsFH = iConfig.getParameter<edm::InputTag>("recHitsFH");
0091   auto recHitsBH = iConfig.getParameter<edm::InputTag>("recHitsBH");
0092   auto caloParticles = iConfig.getParameter<edm::InputTag>("caloParticles");
0093   auto hgcalMultiClusters = iConfig.getParameter<edm::InputTag>("hgcalMultiClusters");
0094   auto electrons = iConfig.getParameter<edm::InputTag>("electrons");
0095   auto photons = iConfig.getParameter<edm::InputTag>("photons");
0096   if (detector == "all") {
0097     recHitsEE_ = consumes<HGCRecHitCollection>(recHitsEE);
0098     recHitsFH_ = consumes<HGCRecHitCollection>(recHitsFH);
0099     recHitsBH_ = consumes<HGCRecHitCollection>(recHitsBH);
0100     algo_ = 1;
0101   } else if (detector == "EM") {
0102     recHitsEE_ = consumes<HGCRecHitCollection>(recHitsEE);
0103     algo_ = 2;
0104   } else if (detector == "HAD") {
0105     recHitsFH_ = consumes<HGCRecHitCollection>(recHitsFH);
0106     recHitsBH_ = consumes<HGCRecHitCollection>(recHitsBH);
0107     algo_ = 3;
0108   }
0109   depletion0_ = iConfig.getParameter<int>("depletionFine");
0110   caloParticles_ = consumes<std::vector<CaloParticle> >(caloParticles);
0111   hgcalMultiClusters_ = consumes<std::vector<reco::PFCluster> >(hgcalMultiClusters);
0112   electrons_ = consumes<std::vector<reco::GsfElectron> >(electrons);
0113   photons_ = consumes<std::vector<reco::Photon> >(photons);
0114 
0115   // size should be HGC layers 52 is enough
0116   Energy_layer_calib_.fill(0.);
0117   Energy_layer_calib_fraction_.fill(0.);
0118 }
0119 
0120 HGCalHitCalibration::~HGCalHitCalibration() {
0121   // do anything here that needs to be done at desctruction time
0122   // (e.g. close files, deallocate resources etc.)
0123 }
0124 
0125 void HGCalHitCalibration::bookHistograms(DQMStore::IBooker& ibooker,
0126                                          edm::Run const& iRun,
0127                                          edm::EventSetup const& /* iSetup */) {
0128   ibooker.cd();
0129   ibooker.setCurrentFolder("HGCalHitCalibration");
0130   h_EoP_CPene_calib_fraction_[depletion0_] = ibooker.book1D("h_EoP_CPene_100_calib_fraction", "", 1000, -0.5, 2.5);
0131   h_EoP_CPene_calib_fraction_[depletion1_] = ibooker.book1D("h_EoP_CPene_200_calib_fraction", "", 1000, -0.5, 2.5);
0132   h_EoP_CPene_calib_fraction_[depletion2_] = ibooker.book1D("h_EoP_CPene_300_calib_fraction", "", 1000, -0.5, 2.5);
0133   h_EoP_CPene_calib_fraction_[scint_] = ibooker.book1D("h_EoP_CPene_scint_calib_fraction", "", 1000, -0.5, 2.5);
0134   hgcal_EoP_CPene_calib_fraction_[depletion0_] =
0135       ibooker.book1D("hgcal_EoP_CPene_100_calib_fraction", "", 1000, -0.5, 2.5);
0136   hgcal_EoP_CPene_calib_fraction_[depletion1_] =
0137       ibooker.book1D("hgcal_EoP_CPene_200_calib_fraction", "", 1000, -0.5, 2.5);
0138   hgcal_EoP_CPene_calib_fraction_[depletion2_] =
0139       ibooker.book1D("hgcal_EoP_CPene_300_calib_fraction", "", 1000, -0.5, 2.5);
0140   hgcal_EoP_CPene_calib_fraction_[scint_] = ibooker.book1D("hgcal_EoP_CPene_scint_calib_fraction", "", 1000, -0.5, 2.5);
0141   hgcal_ele_EoP_CPene_calib_fraction_[depletion0_] =
0142       ibooker.book1D("hgcal_ele_EoP_CPene_100_calib_fraction", "", 1000, -0.5, 2.5);
0143   hgcal_ele_EoP_CPene_calib_fraction_[depletion1_] =
0144       ibooker.book1D("hgcal_ele_EoP_CPene_200_calib_fraction", "", 1000, -0.5, 2.5);
0145   hgcal_ele_EoP_CPene_calib_fraction_[depletion2_] =
0146       ibooker.book1D("hgcal_ele_EoP_CPene_300_calib_fraction", "", 1000, -0.5, 2.5);
0147   hgcal_ele_EoP_CPene_calib_fraction_[scint_] =
0148       ibooker.book1D("hgcal_ele_EoP_CPene_scint_calib_fraction", "", 1000, -0.5, 2.5);
0149   hgcal_photon_EoP_CPene_calib_fraction_[depletion0_] =
0150       ibooker.book1D("hgcal_photon_EoP_CPene_100_calib_fraction", "", 1000, -0.5, 2.5);
0151   hgcal_photon_EoP_CPene_calib_fraction_[depletion1_] =
0152       ibooker.book1D("hgcal_photon_EoP_CPene_200_calib_fraction", "", 1000, -0.5, 2.5);
0153   hgcal_photon_EoP_CPene_calib_fraction_[depletion2_] =
0154       ibooker.book1D("hgcal_photon_EoP_CPene_300_calib_fraction", "", 1000, -0.5, 2.5);
0155   hgcal_photon_EoP_CPene_calib_fraction_[scint_] =
0156       ibooker.book1D("hgcal_photon_EoP_CPene_scint_calib_fraction", "", 1000, -0.5, 2.5);
0157   LayerOccupancy_ = ibooker.book1D("LayerOccupancy", "", layers_, 0., (float)layers_);
0158 }
0159 
0160 void HGCalHitCalibration::fillWithRecHits(std::map<DetId, const HGCRecHit*>& hitmap,
0161                                           const DetId hitid,
0162                                           const unsigned int hitlayer,
0163                                           const float fraction,
0164                                           int& seedDet,
0165                                           float& seedEnergy) {
0166   if (hitmap.find(hitid) == hitmap.end()) {
0167     // Hit was not reconstructed
0168     IfLogTrace(debug_ > 0, "HGCalHitCalibration")
0169         << ">>> Failed to find detid " << std::hex << hitid.rawId() << std::dec << " Det " << hitid.det() << " Subdet "
0170         << hitid.subdetId() << std::endl;
0171     return;
0172   }
0173   if ((hitid.det() != DetId::Forward) && (hitid.det() != DetId::HGCalEE) && (hitid.det() != DetId::HGCalHSi) &&
0174       (hitid.det() != DetId::HGCalHSc)) {
0175     IfLogTrace(debug_ > 0, "HGCalHitCalibration")
0176         << ">>> Wrong type of detid " << std::hex << hitid.rawId() << std::dec << " Det " << hitid.det() << " Subdet "
0177         << hitid.subdetId() << std::endl;
0178     return;
0179   }
0180 
0181   unsigned int layer = recHitTools_.getLayerWithOffset(hitid);
0182   assert(hitlayer == layer);
0183   Energy_layer_calib_fraction_[layer] += hitmap[hitid]->energy() * fraction;
0184   LayerOccupancy_->Fill(layer);
0185   if (seedEnergy < hitmap[hitid]->energy()) {
0186     seedEnergy = hitmap[hitid]->energy();
0187     seedDet = recHitTools_.getSiThickness(hitid);
0188     if (hitid.det() == DetId::HGCalHSc) {
0189       seedDet = scint_;
0190     }
0191   }
0192 }
0193 
0194 void HGCalHitCalibration::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0195   using namespace edm;
0196 
0197   edm::ESHandle<CaloGeometry> geom = iSetup.getHandle(caloGeomToken_);
0198   recHitTools_.setGeometry(*geom);
0199 
0200   Handle<HGCRecHitCollection> recHitHandleEE;
0201   Handle<HGCRecHitCollection> recHitHandleFH;
0202   Handle<HGCRecHitCollection> recHitHandleBH;
0203 
0204   Handle<std::vector<CaloParticle> > caloParticleHandle;
0205   iEvent.getByToken(caloParticles_, caloParticleHandle);
0206   const std::vector<CaloParticle>& caloParticles = *caloParticleHandle;
0207 
0208   Handle<std::vector<reco::PFCluster> > hgcalMultiClustersHandle;
0209   iEvent.getByToken(hgcalMultiClusters_, hgcalMultiClustersHandle);
0210 
0211   Handle<std::vector<reco::GsfElectron> > PFElectronHandle;
0212   iEvent.getByToken(electrons_, PFElectronHandle);
0213 
0214   Handle<std::vector<reco::Photon> > PFPhotonHandle;
0215   iEvent.getByToken(photons_, PFPhotonHandle);
0216 
0217   // make a map detid-rechit
0218   std::map<DetId, const HGCRecHit*> hitmap;
0219   switch (algo_) {
0220     case 1: {
0221       iEvent.getByToken(recHitsEE_, recHitHandleEE);
0222       iEvent.getByToken(recHitsFH_, recHitHandleFH);
0223       iEvent.getByToken(recHitsBH_, recHitHandleBH);
0224       const auto& rechitsEE = *recHitHandleEE;
0225       const auto& rechitsFH = *recHitHandleFH;
0226       const auto& rechitsBH = *recHitHandleBH;
0227       for (unsigned int i = 0; i < rechitsEE.size(); ++i) {
0228         hitmap[rechitsEE[i].detid()] = &rechitsEE[i];
0229       }
0230       for (unsigned int i = 0; i < rechitsFH.size(); ++i) {
0231         hitmap[rechitsFH[i].detid()] = &rechitsFH[i];
0232       }
0233       for (unsigned int i = 0; i < rechitsBH.size(); ++i) {
0234         hitmap[rechitsBH[i].detid()] = &rechitsBH[i];
0235       }
0236       break;
0237     }
0238     case 2: {
0239       iEvent.getByToken(recHitsEE_, recHitHandleEE);
0240       const HGCRecHitCollection& rechitsEE = *recHitHandleEE;
0241       for (unsigned int i = 0; i < rechitsEE.size(); i++) {
0242         hitmap[rechitsEE[i].detid()] = &rechitsEE[i];
0243       }
0244       break;
0245     }
0246     case 3: {
0247       iEvent.getByToken(recHitsFH_, recHitHandleFH);
0248       iEvent.getByToken(recHitsBH_, recHitHandleBH);
0249       const auto& rechitsFH = *recHitHandleFH;
0250       const auto& rechitsBH = *recHitHandleBH;
0251       for (unsigned int i = 0; i < rechitsFH.size(); i++) {
0252         hitmap[rechitsFH[i].detid()] = &rechitsFH[i];
0253       }
0254       for (unsigned int i = 0; i < rechitsBH.size(); i++) {
0255         hitmap[rechitsBH[i].detid()] = &rechitsBH[i];
0256       }
0257       break;
0258     }
0259     default:
0260       assert(false);
0261       break;
0262   }
0263 
0264   // loop over caloParticles
0265   int seedDet = 0;
0266   float seedEnergy = 0.;
0267   IfLogTrace(debug_ > 0, "HGCalHitCalibration") << "Number of caloParticles: " << caloParticles.size() << std::endl;
0268   for (const auto& it_caloPart : caloParticles) {
0269     if (it_caloPart.g4Tracks()[0].eventId().event() != 0 or it_caloPart.g4Tracks()[0].eventId().bunchCrossing() != 0) {
0270       LogDebug("HGCalHitCalibration") << "Excluding CaloParticles from event: "
0271                                       << it_caloPart.g4Tracks()[0].eventId().event()
0272                                       << " with BX: " << it_caloPart.g4Tracks()[0].eventId().bunchCrossing()
0273                                       << std::endl;
0274       continue;
0275     }
0276 
0277     const SimClusterRefVector& simClusterRefVector = it_caloPart.simClusters();
0278     Energy_layer_calib_.fill(0.);
0279     Energy_layer_calib_fraction_.fill(0.);
0280 
0281     seedDet = 0;
0282     seedEnergy = 0.;
0283     for (const auto& it_sc : simClusterRefVector) {
0284       const SimCluster& simCluster = (*(it_sc));
0285       IfLogTrace(debug_ > 1, "HGCalHitCalibration")
0286           << ">>> SC.energy(): " << simCluster.energy() << " SC.simEnergy(): " << simCluster.simEnergy() << std::endl;
0287       const std::vector<std::pair<uint32_t, float> >& hits_and_fractions = simCluster.hits_and_fractions();
0288 
0289       // loop over hits
0290       for (const auto& it_haf : hits_and_fractions) {
0291         unsigned int hitlayer = recHitTools_.getLayerWithOffset(it_haf.first);
0292         DetId hitid = (it_haf.first);
0293         // dump raw RecHits and match
0294         if (rawRecHits_) {
0295           if ((hitid.det() == DetId::Forward &&
0296                (hitid.subdetId() == HGCEE || hitid.subdetId() == HGCHEF || hitid.subdetId() == HGCHEB)) ||
0297               (hitid.det() == DetId::HGCalEE) || (hitid.det() == DetId::HGCalHSi) || (hitid.det() == DetId::HGCalHSc) ||
0298               (hitid.det() == DetId::Hcal && hitid.subdetId() == HcalEndcap))
0299             fillWithRecHits(hitmap, hitid, hitlayer, it_haf.second, seedDet, seedEnergy);
0300         }
0301       }  // end simHit
0302     }    // end simCluster
0303 
0304     auto sumCalibRecHitCalib_fraction =
0305         std::accumulate(Energy_layer_calib_fraction_.begin(), Energy_layer_calib_fraction_.end(), 0.);
0306     IfLogTrace(debug_ > 0, "HGCalHitCalibration")
0307         << ">>> MC Energy: " << it_caloPart.energy() << " reco energy: " << sumCalibRecHitCalib_fraction << std::endl;
0308     if (h_EoP_CPene_calib_fraction_.count(seedDet))
0309       h_EoP_CPene_calib_fraction_[seedDet]->Fill(sumCalibRecHitCalib_fraction / it_caloPart.energy());
0310 
0311     // Loop on reconstructed SC.
0312     const auto& clusters = *hgcalMultiClustersHandle;
0313     float total_energy = 0.;
0314     float max_dR2 = 0.0025;
0315     auto closest =
0316         std::min_element(clusters.begin(), clusters.end(), [&](const reco::PFCluster& a, const reco::PFCluster& b) {
0317           auto dR2_a = reco::deltaR2(it_caloPart, a);
0318           auto dR2_b = reco::deltaR2(it_caloPart, b);
0319           auto ERatio_a = a.correctedEnergy() / it_caloPart.energy();
0320           auto ERatio_b = b.correctedEnergy() / it_caloPart.energy();
0321           // If both clusters are within 0.0025, mark as first (min) the
0322           // element with the highest ratio against the SimCluster
0323           if (dR2_a < max_dR2 && dR2_b < max_dR2)
0324             return ERatio_a > ERatio_b;
0325           return dR2_a < dR2_b;
0326         });
0327     if (closest != clusters.end() && reco::deltaR2(*closest, it_caloPart) < 0.01) {
0328       total_energy = closest->correctedEnergy();
0329       seedDet = recHitTools_.getSiThickness(closest->seed());
0330       if (closest->seed().det() == DetId::HGCalHSc) {
0331         seedDet = scint_;
0332       }
0333       if (hgcal_EoP_CPene_calib_fraction_.count(seedDet)) {
0334         hgcal_EoP_CPene_calib_fraction_[seedDet]->Fill(total_energy / it_caloPart.energy());
0335       }
0336     }
0337 
0338     auto closest_fcn = [&](auto const& a, auto const& b) {
0339       auto dR2_a = reco::deltaR2(it_caloPart, a);
0340       auto dR2_b = reco::deltaR2(it_caloPart, b);
0341       auto ERatio_a = a.energy() / it_caloPart.energy();
0342       auto ERatio_b = b.energy() / it_caloPart.energy();
0343       // If both clusters are within 0.0025, mark as first (min) the
0344       // element with the highest ratio against the SimCluster
0345       if (dR2_a < max_dR2 && dR2_b < max_dR2)
0346         return ERatio_a > ERatio_b;
0347       return dR2_a < dR2_b;
0348     };
0349     // ELECTRONS in HGCAL
0350     if (PFElectronHandle.isValid()) {
0351       auto const& ele = (*PFElectronHandle);
0352       auto closest = std::min_element(ele.begin(), ele.end(), closest_fcn);
0353       if (closest != ele.end() &&
0354           (closest->superCluster()->seed()->seed().det() == DetId::Forward ||
0355            closest->superCluster()->seed()->seed().det() == DetId::HGCalEE) &&
0356           reco::deltaR2(*closest, it_caloPart) < 0.01) {
0357         seedDet = recHitTools_.getSiThickness(closest->superCluster()->seed()->seed());
0358         if (closest->superCluster()->seed()->seed().det() == DetId::HGCalHSc) {
0359           seedDet = scint_;
0360         }
0361         if (hgcal_ele_EoP_CPene_calib_fraction_.count(seedDet)) {
0362           hgcal_ele_EoP_CPene_calib_fraction_[seedDet]->Fill(closest->energy() / it_caloPart.energy());
0363         }
0364       }
0365     }
0366 
0367     // PHOTONS in HGCAL
0368     if (PFPhotonHandle.isValid()) {
0369       auto const& photon = (*PFPhotonHandle);
0370       auto closest = std::min_element(photon.begin(), photon.end(), closest_fcn);
0371       if (closest != photon.end() &&
0372           (closest->superCluster()->seed()->seed().det() == DetId::Forward ||
0373            closest->superCluster()->seed()->seed().det() == DetId::HGCalEE) &&
0374           reco::deltaR2(*closest, it_caloPart) < 0.01) {
0375         seedDet = recHitTools_.getSiThickness(closest->superCluster()->seed()->seed());
0376         if (hgcal_photon_EoP_CPene_calib_fraction_.count(seedDet)) {
0377           hgcal_photon_EoP_CPene_calib_fraction_[seedDet]->Fill(closest->energy() / it_caloPart.energy());
0378         }
0379       }
0380     }
0381   }  // end caloparticle
0382 }
0383 
0384 // ------------ method fills 'descriptions' with the allowed parameters for the
0385 // module  ------------
0386 void HGCalHitCalibration::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0387   edm::ParameterSetDescription desc;
0388   desc.add<int>("debug", 0);
0389   desc.add<bool>("rawRecHits", true);
0390   desc.add<std::string>("detector", "all");
0391   desc.add<int>("depletionFine", 120);
0392   desc.add<edm::InputTag>("caloParticles", edm::InputTag("mix", "MergedCaloTruth"));
0393   desc.add<edm::InputTag>("recHitsEE", edm::InputTag("HGCalRecHit", "HGCEERecHits"));
0394   desc.add<edm::InputTag>("recHitsFH", edm::InputTag("HGCalRecHit", "HGCHEFRecHits"));
0395   desc.add<edm::InputTag>("recHitsBH", edm::InputTag("HGCalRecHit", "HGCHEBRecHits"));
0396   desc.add<edm::InputTag>("hgcalMultiClusters", edm::InputTag("particleFlowClusterHGCal"));
0397   desc.add<edm::InputTag>("electrons", edm::InputTag("ecalDrivenGsfElectronsHGC"));
0398   desc.add<edm::InputTag>("photons", edm::InputTag("photonsHGC"));
0399   descriptions.add("hgcalHitCalibrationDefault", desc);
0400 }
0401 
0402 // define this as a plug-in
0403 DEFINE_FWK_MODULE(HGCalHitCalibration);