CMS 3D CMS Logo

Public Member Functions | Private Attributes

HcalSimpleReconstructor Class Reference

#include <HcalSimpleReconstructor.h>

Inheritance diagram for HcalSimpleReconstructor:
edm::EDProducer edm::ProducerBase edm::ProductRegistryHelper

List of all members.

Public Member Functions

 HcalSimpleReconstructor (const edm::ParameterSet &ps)
virtual void produce (edm::Event &e, const edm::EventSetup &c)
virtual ~HcalSimpleReconstructor ()

Private Attributes

DetId::Detector det_
bool dropZSmarkedPassed_
edm::InputTag inputLabel_
HcalSimpleRecAlgo reco_
int subdet_
HcalOtherSubdetector subdetOther_

Detailed Description

Date:
2009/03/27 17:38:31
Revision:
1.3
Author:
J. Mans - Minnesota
Date:
2010/01/21 14:37:40
Revision:
1.1
Author:
E. Garcia - CSU Based on HcalSimpleReconstructor.h by J. Mans

Definition at line 21 of file HcalSimpleReconstructor.h.


Constructor & Destructor Documentation

HcalSimpleReconstructor::HcalSimpleReconstructor ( const edm::ParameterSet ps) [explicit]

Definition at line 18 of file HcalSimpleReconstructor.cc.

References gather_cfg::cout, edm::ParameterSet::getParameter(), HcalBarrel, HcalForward, HcalOuter, and subdet_.

                                                                           :
  reco_(conf.getParameter<int>("firstSample"),conf.getParameter<int>("samplesToAdd"),conf.getParameter<bool>("correctForTimeslew"),
        conf.getParameter<bool>("correctForPhaseContainment"),conf.getParameter<double>("correctionPhaseNS")),
  det_(DetId::Hcal),
  inputLabel_(conf.getParameter<edm::InputTag>("digiLabel")),
  dropZSmarkedPassed_(conf.getParameter<bool>("dropZSmarkedPassed"))
{
  std::string subd=conf.getParameter<std::string>("Subdetector");
  if (!strcasecmp(subd.c_str(),"HBHE")) {
    subdet_=HcalBarrel;
    produces<HBHERecHitCollection>();
  } else if (!strcasecmp(subd.c_str(),"HO")) {
    subdet_=HcalOuter;
    produces<HORecHitCollection>();
  } else if (!strcasecmp(subd.c_str(),"HF")) {
    subdet_=HcalForward;
    produces<HFRecHitCollection>();
  } else {
    std::cout << "HcalSimpleReconstructor is not associated with a specific subdetector!" << std::endl;
  }       
  
}
HcalSimpleReconstructor::~HcalSimpleReconstructor ( ) [virtual]

Definition at line 41 of file HcalSimpleReconstructor.cc.

                                                  {
}

Member Function Documentation

void HcalSimpleReconstructor::produce ( edm::Event e,
const edm::EventSetup c 
) [virtual]

Implements edm::EDProducer.

Definition at line 44 of file HcalSimpleReconstructor.cc.

References det_, dropZSmarkedPassed_, edm::EventSetup::get(), edm::Event::getByLabel(), DetId::Hcal, HcalBarrel, HcalCalibration, HcalEndcap, HcalForward, HcalOther, HcalOuter, i, inputLabel_, edm::Event::put(), reco_, HcalSimpleRecAlgo::reconstruct(), subdet_, and subdetOther_.

{
  // get conditions
  edm::ESHandle<HcalDbService> conditions;
  eventSetup.get<HcalDbRecord>().get(conditions);
  const HcalQIEShape* shape = conditions->getHcalShape (); // this one is generic
  
  
  if (det_==DetId::Hcal) {
    if (subdet_==HcalBarrel || subdet_==HcalEndcap) {
      edm::Handle<HBHEDigiCollection> digi;
      
      e.getByLabel(inputLabel_,digi);
      
      // create empty output
      std::auto_ptr<HBHERecHitCollection> rec(new HBHERecHitCollection);
      rec->reserve(digi->size());
      // run the algorithm
      HBHEDigiCollection::const_iterator i;
      for (i=digi->begin(); i!=digi->end(); i++) {
        HcalDetId cell = i->id();
        // rof 27.03.09: drop ZS marked and passed digis:
        if (dropZSmarkedPassed_)
          if (i->zsMarkAndPass()) continue;

        const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
        const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
        HcalCoderDb coder (*channelCoder, *shape);
        rec->push_back(reco_.reconstruct(*i,coder,calibrations));
      }
      // return result
      e.put(rec);
    } else if (subdet_==HcalOuter) {
      edm::Handle<HODigiCollection> digi;
      e.getByLabel(inputLabel_,digi);
      
      // create empty output
      std::auto_ptr<HORecHitCollection> rec(new HORecHitCollection);
      rec->reserve(digi->size());
      // run the algorithm
      HODigiCollection::const_iterator i;
      for (i=digi->begin(); i!=digi->end(); i++) {
        HcalDetId cell = i->id();
        // rof 27.03.09: drop ZS marked and passed digis:
        if (dropZSmarkedPassed_)
          if (i->zsMarkAndPass()) continue;

        const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
        const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
        HcalCoderDb coder (*channelCoder, *shape);
        rec->push_back(reco_.reconstruct(*i,coder,calibrations));
      }
      // return result
      e.put(rec);    
    } else if (subdet_==HcalForward) {
      edm::Handle<HFDigiCollection> digi;
      e.getByLabel(inputLabel_,digi);
      
      // create empty output
      std::auto_ptr<HFRecHitCollection> rec(new HFRecHitCollection);
      rec->reserve(digi->size());
      // run the algorithm
      HFDigiCollection::const_iterator i;
      for (i=digi->begin(); i!=digi->end(); i++) {
        HcalDetId cell = i->id();        
        // rof 27.03.09: drop ZS marked and passed digis:
        if (dropZSmarkedPassed_)
          if (i->zsMarkAndPass()) continue;
 
        const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
        const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
        HcalCoderDb coder (*channelCoder, *shape);
        rec->push_back(reco_.reconstruct(*i,coder,calibrations));
      }
      // return result
      e.put(rec);     
    } else if (subdet_==HcalOther && subdetOther_==HcalCalibration) {
      edm::Handle<HcalCalibDigiCollection> digi;
      e.getByLabel(inputLabel_,digi);
      
      // create empty output
      std::auto_ptr<HcalCalibRecHitCollection> rec(new HcalCalibRecHitCollection);
      rec->reserve(digi->size());
      // run the algorithm
      HcalCalibDigiCollection::const_iterator i;
      for (i=digi->begin(); i!=digi->end(); i++) {
        HcalCalibDetId cell = i->id();    
        // rof 27.03.09: drop ZS marked and passed digis:
        if (dropZSmarkedPassed_)
          if (i->zsMarkAndPass()) continue;

        const HcalCalibrations& calibrations=conditions->getHcalCalibrations(cell);
        const HcalQIECoder* channelCoder = conditions->getHcalCoder (cell);
        HcalCoderDb coder (*channelCoder, *shape);
        rec->push_back(reco_.reconstruct(*i,coder,calibrations));
      }
      // return result
      e.put(rec);     
    }
  } 
}

Member Data Documentation

Definition at line 28 of file HcalSimpleReconstructor.h.

Referenced by produce().

Definition at line 33 of file HcalSimpleReconstructor.h.

Referenced by produce().

Definition at line 31 of file HcalSimpleReconstructor.h.

Referenced by produce().

Definition at line 27 of file HcalSimpleReconstructor.h.

Referenced by produce().

Definition at line 29 of file HcalSimpleReconstructor.h.

Referenced by HcalSimpleReconstructor(), and produce().

Definition at line 30 of file HcalSimpleReconstructor.h.

Referenced by produce().