Belle II Software  release-05-01-25
PXDRawHitMaskingModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Bjoern.Spruck@belle2.org *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <pxd/modules/pxdHelper/PXDRawHitMaskingModule.h>
12 
13 #include <framework/datastore/DataStore.h>
14 #include <framework/datastore/RelationArray.h>
15 #include <framework/logging/Logger.h>
16 
17 #include <vxd/geometry/GeoCache.h>
18 
19 #include <pxd/reconstruction/Pixel.h>
20 #include <pxd/reconstruction/PXDPixelMasker.h>
21 
22 #include <set>
23 
24 using namespace std;
25 using namespace Belle2;
26 using namespace Belle2::PXD;
27 
28 //-----------------------------------------------------------------
29 // Register the Module
30 //-----------------------------------------------------------------
31 REG_MODULE(PXDRawHitMasking);
32 
33 //-----------------------------------------------------------------
34 // Implementation
35 //-----------------------------------------------------------------
36 
37 PXDRawHitMaskingModule::PXDRawHitMaskingModule() : Module()
38 {
39  //Set module properties
40  setDescription("This module masks the input collection of PXDRawHits into "
41  "a PXDRawHit output collection");
43 
44  addParam("zeroSuppressionCut", m_0cut, "Minimum charge for a digit to carry", 0);
45  addParam("trimOutOfRange", m_trimOutOfRange, "Discard rawhits whit out-of-range coordinates", true);
46  addParam("rawHits", m_storeRawHitsName, "PXDRawHit collection name", string(""));
47  addParam("rawHitsOut", m_storeRawHitsNameOut, "PXDRawHit Out collection name", string("PXDRawHitsOut"));
48 }
49 
50 
52 {
53  //Register collections
54  m_pxdRawHit.isRequired(m_storeRawHitsName);
55  m_pxdRawHitOut.registerInDataStore(m_storeRawHitsNameOut);
56 }
57 
59 {
60 
61  // if no input, nothing to do
62  if (!m_pxdRawHit || !m_pxdRawHit.getEntries()) return;
63 
65 
66  for (auto& it : m_pxdRawHit) {
67  VxdID sensorID = it.getSensorID();
68  if (!geo.validSensorID(sensorID)) {
69  B2WARNING("Malformed PXDRawHit, VxdID $" << hex << sensorID.getID() << ", dropping. (" << sensorID << ")");
70  continue;
71  }
72  if (m_trimOutOfRange && !goodHit(it))
73  continue;
74  // Zero-suppression cut
75  if (it.getCharge() < m_0cut) continue;
76 
77  // We need some protection against crap data
78  if (sensorID.getLayerNumber() && sensorID.getLadderNumber() && sensorID.getSensorNumber()) {
79  if (PXDPixelMasker::getInstance().pixelOK(sensorID, it.getColumn(), it.getRow())) {
80  m_pxdRawHitOut.appendNew(it);
81  }
82  }
83  }
84 }
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::PXD::PXDRawHitMaskingModule::goodHit
bool goodHit(const PXDRawHit &rawhit) const
Utility function to check pixel coordinates.
Definition: PXDRawHitMaskingModule.h:59
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Module::c_ParallelProcessingCertified
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition: Module.h:82
Belle2::VxdID::getLadderNumber
baseType getLadderNumber() const
Get the ladder id.
Definition: VxdID.h:108
Belle2::PXD::PXDRawHitMaskingModule::m_0cut
int m_0cut
Minimum charge for a digit to carry.
Definition: PXDRawHitMaskingModule.h:72
Belle2::VxdID::getID
baseType getID() const
Get the unique id.
Definition: VxdID.h:104
Belle2::VXD::GeoCache::validSensorID
bool validSensorID(Belle2::VxdID id) const
Check that id is a valid sensor number.
Definition: GeoCache.cc:53
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::Module::setPropertyFlags
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:210
Belle2::PXD::PXDRawHitMaskingModule::event
virtual void event() override
do the filtering
Definition: PXDRawHitMaskingModule.cc:58
Belle2::PXD
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
Definition: PXDCalibrationUtilities.h:28
Belle2::VXD::GeoCache::getInstance
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:215
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::PXD::PXDRawHitMaskingModule::m_storeRawHitsName
std::string m_storeRawHitsName
Name of the collection to use for PXDRawHits.
Definition: PXDRawHitMaskingModule.h:68
Belle2::PXD::PXDRawHitMaskingModule::m_storeRawHitsNameOut
std::string m_storeRawHitsNameOut
Name of the collection to use for Output PXDRawHits.
Definition: PXDRawHitMaskingModule.h:70
Belle2::PXD::PXDRawHitMaskingModule::m_pxdRawHitOut
StoreArray< PXDRawHit > m_pxdRawHitOut
Required output for PXDRawHit.
Definition: PXDRawHitMaskingModule.h:57
Belle2::VxdID::getSensorNumber
baseType getSensorNumber() const
Get the sensor id.
Definition: VxdID.h:110
Belle2::Module::addParam
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition: Module.h:562
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::PXD::PXDRawHitMaskingModule::initialize
virtual void initialize() override
Initialize the module.
Definition: PXDRawHitMaskingModule.cc:51
Belle2::VxdID::getLayerNumber
baseType getLayerNumber() const
Get the layer id.
Definition: VxdID.h:106
Belle2::PXD::PXDRawHitMaskingModule::m_pxdRawHit
StoreArray< PXDRawHit > m_pxdRawHit
Required input for PXDRawHit.
Definition: PXDRawHitMaskingModule.h:56
Belle2::PXD::PXDRawHitMaskingModule::m_trimOutOfRange
bool m_trimOutOfRange
Discard out-of-range hits.
Definition: PXDRawHitMaskingModule.h:74
Belle2::PXD::PXDPixelMasker::getInstance
static PXDPixelMasker & getInstance()
Main (and only) way to access the PXDPixelMasker.
Definition: PXDPixelMasker.cc:37