#include <EventFilter/SiPixelRawToDigi/plugins/SiPixelRawToDigi.h>
Public Member Functions | |
virtual void | beginJob (const edm::EventSetup &) |
initialisation. Retrieves cabling map from EventSetup. | |
virtual void | endJob () |
dummy end of job | |
virtual void | produce (edm::Event &, const edm::EventSetup &) |
get data, convert to digis attach againe to Event | |
SiPixelRawToDigi (const edm::ParameterSet &) | |
ctor | |
virtual | ~SiPixelRawToDigi () |
dtor | |
Private Attributes | |
bool | checkOrder |
edm::ParameterSet | config_ |
unsigned long | eventCounter_ |
const SiPixelFedCablingMap * | fedCablingMap_ |
std::vector< int > | fedList_ |
TH1D * | hCPU |
TH1D * | hDigi |
bool | includeErrors |
TFile * | rootFile |
R2DTimerObserver * | theTimer |
Definition at line 20 of file SiPixelRawToDigi.h.
SiPixelRawToDigi::SiPixelRawToDigi | ( | const edm::ParameterSet & | conf | ) | [explicit] |
ctor
Definition at line 33 of file SiPixelRawToDigi.cc.
References checkOrder, config_, lat::endl(), edm::ParameterSet::getUntrackedParameter(), hCPU, hDigi, includeErrors, rootFile, and theTimer.
00034 : eventCounter_(0), 00035 config_(conf), 00036 fedCablingMap_(0), 00037 hCPU(0), hDigi(0), rootFile(0), 00038 theTimer(0) 00039 { 00040 edm::LogInfo("SiPixelRawToDigi")<< " HERE ** constructor!" << endl; 00041 bool timing = config_.getUntrackedParameter<bool>("Timing",false); 00042 includeErrors = config_.getUntrackedParameter<bool>("IncludeErrors",false); 00043 checkOrder = config_.getUntrackedParameter<bool>("CheckPixelOrder",false); 00044 // Products 00045 produces< edm::DetSetVector<PixelDigi> >(); 00046 if(includeErrors) produces< edm::DetSetVector<SiPixelRawDataError> >(); 00047 // Timing 00048 if (timing) { 00049 theTimer = new R2DTimerObserver("**** MY TIMING REPORT ***"); 00050 rootFile = new TFile("analysis.root", "RECREATE", "my histograms"); 00051 hCPU = new TH1D ("hCPU","hCPU",60,0.,0.030); 00052 hDigi = new TH1D("hDigi","hDigi",50,0.,15000.); 00053 } 00054 }
SiPixelRawToDigi::~SiPixelRawToDigi | ( | ) | [virtual] |
dtor
Definition at line 58 of file SiPixelRawToDigi.cc.
References rootFile, and theTimer.
00058 { 00059 edm::LogInfo("SiPixelRawToDigi") << " HERE ** SiPixelRawToDigi destructor!"; 00060 00061 if (theTimer) { 00062 rootFile->Write(); 00063 delete theTimer; 00064 } 00065 00066 }
void SiPixelRawToDigi::beginJob | ( | const edm::EventSetup & | c | ) | [virtual] |
initialisation. Retrieves cabling map from EventSetup.
Reimplemented from edm::EDProducer.
Definition at line 70 of file SiPixelRawToDigi.cc.
dummy end of job
Reimplemented from edm::EDProducer.
Definition at line 33 of file SiPixelRawToDigi.h.
void SiPixelRawToDigi::produce | ( | edm::Event & | ev, | |
const edm::EventSetup & | es | |||
) | [virtual] |
get data, convert to digis attach againe to Event
Implements edm::EDProducer.
Definition at line 75 of file SiPixelRawToDigi.cc.
References edm::ESWatcher< T >::check(), checkOrder, collection, config_, GenMuonPlsPt100GeV_cfg::cout, edm::DetSet< T >::data, lat::endl(), HLT_VtxMuL3::errors, fedCablingMap_, SiPixelFedCablingMap::fedList(), fedList_, edm::EventSetup::get(), edm::Event::getByLabel(), edm::ParameterSet::getUntrackedParameter(), hCPU, hDigi, includeErrors, PixelDataFormatter::interpretRawData(), it, label, R2DTimerObserver::lastMeasurement(), LogDebug, PixelDataFormatter::nDigis(), PixelDataFormatter::nWords(), edm::ESHandle< T >::product(), edm::Event::put(), PixelDataFormatter::setErrorStatus(), R2DTimerObserver::start(), R2DTimerObserver::stop(), and theTimer.
00077 { 00078 00079 // initialize cabling map or update if necessary 00080 static edm::ESWatcher<SiPixelFedCablingMapRcd> recordWatcher; 00081 if (recordWatcher.check( es )) { 00082 edm::ESHandle<SiPixelFedCablingMap> map; 00083 es.get<SiPixelFedCablingMapRcd>().get( map ); 00084 LogDebug("map version:")<< map->version(); 00085 fedCablingMap_ = map.product(); 00086 typedef std::vector<const sipixelobjects::PixelFEDCabling *>::iterator FLI; 00087 std::vector<const sipixelobjects::PixelFEDCabling *> feds = map.product()->fedList(); 00088 00089 for (FLI fedIds = feds.begin(); fedIds != feds.end(); fedIds++) { 00090 int fedId = (*fedIds)->id(); 00091 fedList_.push_back( fedId ); 00092 } 00093 } 00094 00095 edm::Handle<FEDRawDataCollection> buffers; 00096 static string label = config_.getUntrackedParameter<string>("InputLabel","source"); 00097 static string instance = config_.getUntrackedParameter<string>("InputInstance",""); 00098 ev.getByLabel( label, instance, buffers); 00099 00100 // create product (digis & errors) 00101 std::auto_ptr< edm::DetSetVector<PixelDigi> > collection( new edm::DetSetVector<PixelDigi> ); 00102 std::auto_ptr< edm::DetSetVector<SiPixelRawDataError> > errorcollection( new edm::DetSetVector<SiPixelRawDataError> ); 00103 static int ndigis = 0; 00104 static int nwords = 0; 00105 00106 PixelDataFormatter formatter(fedCablingMap_); 00107 formatter.setErrorStatus(includeErrors, checkOrder); 00108 00109 if (theTimer) theTimer->start(); 00110 bool errorsInEvent = false; 00111 00112 typedef std::vector<int>::iterator IF; 00113 for (IF theFed = fedList_.begin(); theFed != fedList_.end(); theFed++) { 00114 int fedId = *theFed; 00115 LogDebug("SiPixelRawToDigi")<< " PRODUCE DIGI FOR FED: " << fedId << endl; 00116 PixelDataFormatter::Digis digis; 00117 PixelDataFormatter::Errors errors; 00118 00119 //get event data for this fed 00120 const FEDRawData& fedRawData = buffers->FEDData( fedId ); 00121 00122 //convert data to digi and strip off errors 00123 formatter.interpretRawData( errorsInEvent, fedId, fedRawData, digis, errors); 00124 00125 //pack digi into collection 00126 typedef PixelDataFormatter::Digis::iterator ID; 00127 for (ID it = digis.begin(); it != digis.end(); it++) { 00128 // uint32_t detid = it->id; 00129 uint32_t detid = it->first; 00130 edm::DetSet<PixelDigi>& detSet = collection->find_or_insert(detid); 00131 // detSet.data = it->data; 00132 detSet.data = it->second; 00133 } // end digi loop over detIds 00134 //pack errors into collection 00135 if(includeErrors) { 00136 typedef PixelDataFormatter::Errors::iterator IE; 00137 for (IE is = errors.begin(); is != errors.end(); is++) { 00138 // uint32_t errordetid = is->id; 00139 uint32_t errordetid = is->first; 00140 edm::DetSet<SiPixelRawDataError>& errorDetSet = errorcollection->find_or_insert(errordetid); 00141 // detSet.data = is->data; 00142 errorDetSet.data = is->second; 00143 } // end error loop over detIds 00144 } // end if(includeErrors) 00145 } // end loop over feds 00146 if (errorsInEvent) LogDebug("SiPixelRawToDigi") << "Error words were stored in this event"; 00147 00148 if (theTimer) { 00149 theTimer->stop(); 00150 cout << "TIMING IS: (real)" << theTimer->lastMeasurement().real() << endl; 00151 ndigis += formatter.nDigis(); 00152 nwords += formatter.nWords(); 00153 cout << " (Words/Digis) this ev: "<<formatter.nWords()<<"/"<<formatter.nDigis() 00154 << "--- all :"<<nwords<<"/"<<ndigis<<endl; 00155 hCPU->Fill( theTimer->lastMeasurement().real() ); 00156 hDigi->Fill(formatter.nDigis()); 00157 } 00158 00159 //send digis and errors back to framework 00160 ev.put( collection ); 00161 if(includeErrors) ev.put( errorcollection ); 00162 }
bool SiPixelRawToDigi::checkOrder [private] |
edm::ParameterSet SiPixelRawToDigi::config_ [private] |
unsigned long SiPixelRawToDigi::eventCounter_ [private] |
Definition at line 40 of file SiPixelRawToDigi.h.
const SiPixelFedCablingMap* SiPixelRawToDigi::fedCablingMap_ [private] |
std::vector<int> SiPixelRawToDigi::fedList_ [private] |
TH1D* SiPixelRawToDigi::hCPU [private] |
TH1D * SiPixelRawToDigi::hDigi [private] |
bool SiPixelRawToDigi::includeErrors [private] |
TFile* SiPixelRawToDigi::rootFile [private] |
Definition at line 46 of file SiPixelRawToDigi.h.
Referenced by SiPixelRawToDigi(), and ~SiPixelRawToDigi().
R2DTimerObserver* SiPixelRawToDigi::theTimer [private] |
Definition at line 47 of file SiPixelRawToDigi.h.
Referenced by produce(), SiPixelRawToDigi(), and ~SiPixelRawToDigi().