CMS 3D CMS Logo

/afs/cern.ch/work/a/aaltunda/public/www/CMSSW_6_2_5/src/RecoLocalCalo/HcalRecProducers/src/HcalSimpleReconstructor.cc

Go to the documentation of this file.
00001 #include "HcalSimpleReconstructor.h"
00002 #include "DataFormats/HcalDigi/interface/HcalDigiCollections.h"
00003 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
00004 #include "DataFormats/Common/interface/EDCollection.h"
00005 #include "DataFormats/Common/interface/Handle.h"
00006 #include "FWCore/Framework/interface/ESHandle.h"
00007 #include "FWCore/Framework/interface/EventSetup.h"
00008 #include "CalibFormats/HcalObjects/interface/HcalCoderDb.h"
00009 #include "CalibFormats/HcalObjects/interface/HcalCalibrations.h"
00010 #include "CalibFormats/HcalObjects/interface/HcalDbService.h"
00011 #include "CalibFormats/HcalObjects/interface/HcalDbRecord.h"
00012 #include "Geometry/CaloTopology/interface/HcalTopology.h"
00013 
00014 #include <iostream>
00015     
00016 HcalSimpleReconstructor::HcalSimpleReconstructor(edm::ParameterSet const& conf):
00017   reco_(conf.getParameter<bool>("correctForTimeslew"),
00018         conf.getParameter<bool>("correctForPhaseContainment"),conf.getParameter<double>("correctionPhaseNS")),
00019   det_(DetId::Hcal),
00020   inputLabel_(conf.getParameter<edm::InputTag>("digiLabel")),
00021   dropZSmarkedPassed_(conf.getParameter<bool>("dropZSmarkedPassed")),
00022   firstSample_(conf.getParameter<int>("firstSample")),
00023   samplesToAdd_(conf.getParameter<int>("samplesToAdd")),
00024   tsFromDB_(conf.getParameter<bool>("tsFromDB")),
00025   upgradeHBHE_(false),
00026   upgradeHF_(false),
00027   paramTS(0),
00028   theTopology(0)
00029 {
00030 
00031   std::string subd=conf.getParameter<std::string>("Subdetector");
00032   if(!strcasecmp(subd.c_str(),"upgradeHBHE")) {
00033      upgradeHBHE_ = true;
00034      produces<HBHERecHitCollection>();
00035   }
00036   else if (!strcasecmp(subd.c_str(),"upgradeHF")) {
00037      upgradeHF_ = true;
00038      produces<HFRecHitCollection>();
00039   }
00040   else if (!strcasecmp(subd.c_str(),"HO")) {
00041     subdet_=HcalOuter;
00042     produces<HORecHitCollection>();
00043   }  
00044   else if (!strcasecmp(subd.c_str(),"HBHE")) {
00045     if( !upgradeHBHE_) {
00046       subdet_=HcalBarrel;
00047       produces<HBHERecHitCollection>();
00048     }
00049   } 
00050   else if (!strcasecmp(subd.c_str(),"HF")) {
00051     if( !upgradeHF_) {
00052     subdet_=HcalForward;
00053     produces<HFRecHitCollection>();
00054     }
00055   } 
00056   else {
00057     std::cout << "HcalSimpleReconstructor is not associated with a specific subdetector!" << std::endl;
00058   }       
00059   
00060 }
00061 
00062 HcalSimpleReconstructor::~HcalSimpleReconstructor() { 
00063   delete paramTS;
00064   delete theTopology;
00065 }
00066 
00067 void HcalSimpleReconstructor::beginRun(edm::Run const&r, edm::EventSetup const & es){
00068   if(tsFromDB_) {
00069     edm::ESHandle<HcalRecoParams> p;
00070     es.get<HcalRecoParamsRcd>().get(p);
00071     paramTS = new HcalRecoParams(*p.product());
00072 
00073     edm::ESHandle<HcalTopology> htopo;
00074     es.get<IdealGeometryRecord>().get(htopo);
00075     theTopology=new HcalTopology(*htopo);
00076     paramTS->setTopo(theTopology);
00077 
00078   }
00079   reco_.beginRun(es);
00080 }
00081 
00082 void HcalSimpleReconstructor::endRun(edm::Run const&r, edm::EventSetup const & es){
00083   if(tsFromDB_ && paramTS) {
00084     delete paramTS;
00085     paramTS = 0;
00086     reco_.endRun();
00087   }
00088 }
00089 
00090 
00091 template<class DIGICOLL, class RECHITCOLL> 
00092 void HcalSimpleReconstructor::process(edm::Event& e, const edm::EventSetup& eventSetup)
00093 {
00094   // get conditions
00095   edm::ESHandle<HcalDbService> conditions;
00096   eventSetup.get<HcalDbRecord>().get(conditions);
00097 
00098   edm::Handle<DIGICOLL> digi;
00099   e.getByLabel(inputLabel_,digi);
00100 
00101   // create empty output
00102   std::auto_ptr<RECHITCOLL> rec(new RECHITCOLL);
00103   rec->reserve(digi->size());
00104   // run the algorithm
00105   int first = firstSample_;
00106   int toadd = samplesToAdd_;
00107   typename DIGICOLL::const_iterator i;
00108   for (i=digi->begin(); i!=digi->end(); i++) {
00109     HcalDetId cell = i->id();
00110     DetId detcell=(DetId)cell;
00111     // rof 27.03.09: drop ZS marked and passed digis:
00112     if (dropZSmarkedPassed_)
00113       if (i->zsMarkAndPass()) continue;
00114 
00115     const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00116     const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00117     const HcalQIEShape* shape = conditions->getHcalShape (channelCoder); 
00118     HcalCoderDb coder (*channelCoder, *shape);
00119 
00120     //>>> firstSample & samplesToAdd
00121     if(tsFromDB_) {
00122       const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00123       first = param_ts->firstSample();
00124       toadd = param_ts->samplesToAdd();
00125     }
00126     rec->push_back(reco_.reconstruct(*i,first,toadd,coder,calibrations));   
00127   }
00128   // return result
00129   e.put(rec);
00130 }
00131 
00132 
00133 void HcalSimpleReconstructor::processUpgrade(edm::Event& e, const edm::EventSetup& eventSetup)
00134 {
00135   // get conditions
00136   edm::ESHandle<HcalDbService> conditions;
00137   eventSetup.get<HcalDbRecord>().get(conditions);
00138 
00139   if(upgradeHBHE_){
00140    
00141     edm::Handle<HBHEUpgradeDigiCollection> digi;
00142     e.getByLabel(inputLabel_, digi);
00143 
00144     // create empty output
00145     std::auto_ptr<HBHERecHitCollection> rec(new HBHERecHitCollection);
00146     rec->reserve(digi->size()); 
00147 
00148     // run the algorithm
00149     int first = firstSample_;
00150     int toadd = samplesToAdd_;
00151     HBHEUpgradeDigiCollection::const_iterator i;
00152     for (i=digi->begin(); i!=digi->end(); i++) {
00153       HcalDetId cell = i->id();
00154       DetId detcell=(DetId)cell;
00155       // rof 27.03.09: drop ZS marked and passed digis:
00156       if (dropZSmarkedPassed_)
00157       if (i->zsMarkAndPass()) continue;
00158       
00159       const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00160       const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00161       const HcalQIEShape* shape = conditions->getHcalShape (channelCoder); 
00162       HcalCoderDb coder (*channelCoder, *shape);
00163 
00164       //>>> firstSample & samplesToAdd
00165       if(tsFromDB_) {
00166         const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00167         first = param_ts->firstSample();
00168         toadd = param_ts->samplesToAdd();
00169       }
00170       rec->push_back(reco_.reconstructHBHEUpgrade(*i,first,toadd,coder,calibrations));
00171 
00172     }
00173 
00174     e.put(rec); // put results
00175   }// End of upgradeHBHE
00176 
00177   if(upgradeHF_){
00178 
00179     edm::Handle<HFUpgradeDigiCollection> digi;
00180     e.getByLabel(inputLabel_, digi);
00181 
00182     // create empty output
00183     std::auto_ptr<HFRecHitCollection> rec(new HFRecHitCollection);
00184     rec->reserve(digi->size()); 
00185 
00186     // run the algorithm
00187     int first = firstSample_;
00188     int toadd = samplesToAdd_;
00189     HFUpgradeDigiCollection::const_iterator i;
00190     for (i=digi->begin(); i!=digi->end(); i++) {
00191       HcalDetId cell = i->id();
00192       DetId detcell=(DetId)cell;
00193       // rof 27.03.09: drop ZS marked and passed digis:
00194       if (dropZSmarkedPassed_)
00195       if (i->zsMarkAndPass()) continue;
00196       
00197       const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
00198       const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
00199       const HcalQIEShape* shape = conditions->getHcalShape (channelCoder); 
00200       HcalCoderDb coder (*channelCoder, *shape);
00201 
00202       //>>> firstSample & samplesToAdd
00203       if(tsFromDB_) {
00204         const HcalRecoParam* param_ts = paramTS->getValues(detcell.rawId());
00205         first = param_ts->firstSample();
00206         toadd = param_ts->samplesToAdd();
00207       }
00208       rec->push_back(reco_.reconstructHFUpgrade(*i,first,toadd,coder,calibrations));
00209 
00210     }  
00211     e.put(rec); // put results
00212   }// End of upgradeHF
00213 
00214 }
00215 
00216 
00217 
00218 void HcalSimpleReconstructor::produce(edm::Event& e, const edm::EventSetup& eventSetup)
00219 {
00220   // HACK related to HB- corrections
00221   if(e.isRealData()) reco_.setForData();
00222  
00223   // What to produce, better to avoid the same subdet Upgrade and regular 
00224   // rechits "clashes"
00225   if(upgradeHBHE_ || upgradeHF_) {
00226       processUpgrade(e, eventSetup);
00227   } else if (det_==DetId::Hcal) {
00228     if ((subdet_==HcalBarrel || subdet_==HcalEndcap) && !upgradeHBHE_) {
00229       process<HBHEDigiCollection, HBHERecHitCollection>(e, eventSetup);
00230     } else if (subdet_==HcalForward && !upgradeHF_) {
00231       process<HFDigiCollection, HFRecHitCollection>(e, eventSetup);
00232     } else if (subdet_==HcalOuter) {
00233       process<HODigiCollection, HORecHitCollection>(e, eventSetup);
00234     } else if (subdet_==HcalOther && subdetOther_==HcalCalibration) {
00235       process<HcalCalibDigiCollection, HcalCalibRecHitCollection>(e, eventSetup);
00236     }
00237   } 
00238 }