Belle II Software  release-08-01-10
Belle2::PXD Namespace Reference

Namespace to encapsulate code needed for simulation and reconstrucion of the PXD. More...

Classes

class  GeoPXDCreator
 The creator for the PXD geometry of the Belle II detector. More...
 
class  SensorInfo
 Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific information. More...
 
class  PXDBackgroundModule
 PXD Background module. More...
 
class  PXDBeamBackHitFilterModule
 The PXDBeamBackHitFilter module. More...
 
class  PXDBgTupleProducerModule
 PXD Background Tuple Producer. More...
 
class  PXDMCBgTupleProducerModule
 PXD MC Background Tuple Producer. More...
 
class  PXDDAQDQMModule
 The PXD DAQ DQM module. More...
 
class  PXDGatedDHCDQMModule
 The PXD Gatint after Injection DQM module. More...
 
class  PXDGatedModeDQMModule
 The PXD for GatedMode DQM module. More...
 
class  PXDInjectionDQMModule
 The PXD Occupancy after Injection DQM module. More...
 
class  PXDRawDQMChipsModule
 The raw PXD DQM module. More...
 
class  PXDRawDQMModule
 The raw PXD DQM module. More...
 
class  PXDROIDQMModule
 The raw PXD DQM module. More...
 
class  PXDBadSensorTagModule
 The PXD bad sensor tagger module. More...
 
class  PXDEventPlotModule
 Plot each event with ROI and Pixels. More...
 
class  PXDRawDumperModule
 Dump Raw PXD/ ONSEN event data. More...
 
class  PXDRawHitMaskingModule
 The PXDRawHitMasking module. More...
 
class  PXDROIPlotModule
 Plot each event with ROI and Pixels. More...
 
class  ActivatePXDClusterPositionEstimatorModule
 The ActivatePXDClusterPositionEstimator module. More...
 
class  ActivatePXDGainCalibratorModule
 The ActivatePXDGainCalibrator module. More...
 
class  ActivatePXDPixelMaskerModule
 The ActivatePXDPixelMasker module. More...
 
class  PXDClusterCheckModule
 The PXDClusterCheck module. More...
 
class  PXDClusterizerModule
 The PXDClusterizer module. More...
 
class  PXDDigitSorterModule
 The PXDDigitSorter module. More...
 
class  PXDRawHitSorterModule
 The PXDRawHitSorter module. More...
 
class  Digit
 Class to represent the coordinates of one pixel. More...
 
class  DigitValue
 Class representing the charge and particle contributions for one pixel. More...
 
class  PXDDigitizerModule
 The PXD Digitizer module. More...
 
class  PXDGatedInfoFillerModule
 PXD Gates Mode infromation on readout gate basis. More...
 
class  PXDPackerErrModule
 The PXDPackerErr module. More...
 
class  PXDPackerModule
 The PXDPacker module. More...
 
class  PXDPostErrorCheckerModule
 The PXD DAQ Post Unpacking Error Check. More...
 
class  PXDReadRawBonnDAQModule
 Module to Load Raw PXD Data from DHH network-dump file and store it as RawPXD in Data Store This is meant for lab use (standalone testing, debugging) without an event builder. More...
 
class  PXDReadRawBonnDAQMatchedModule
 Module to Load BonnDAQ file and store it as RawPXD in Data Store This is meant for lab use (standalone testing, debugging) without an event builder. More...
 
class  PXDReadRawONSENModule
 A class definition of an input module for Sequential ROOT I/O. More...
 
class  PXDUnpackerModule
 The PXDUnpacker module. More...
 
class  PXDUnpackerOldModule
 The PXDUnpacker module. More...
 
class  PXDUnpackerOTModule
 The PXDUnpackerOT module. More...
 
class  ClusterCache
 Class to remember recently assigned clusters This class will remember the current and the last pixel row to allow fast finding of the correct cluster a pixel belongs to. More...
 
class  ClusterCandidate
 Class representing a possible cluster during clustering of the PXD It supports merging of different clusters and keeps track of the highest charge inside the cluster. More...
 
class  ClusterProjection
 Helper struct to collect information about the 1D projection of a Pixel cluster. More...
 
class  NoiseMap
 Base Class to represent pixel dependent Noise Map. More...
 
class  Pixel
 Class to represent one pixel, used in clustering for fast access. More...
 
class  PXDClusterPositionEstimator
 Singleton class that estimates cluster positions taking into account the estimated track incidence angles into the sensor. More...
 
class  PXDClusterShape
 Class to correct estimation of cluster error and position base on its shape. More...
 
class  PXDGainCalibrator
 Singleton class for managing gain corrections for the PXD. More...
 
class  PXDPixelMasker
 Singleton class for managing pixel masking for the PXD. More...
 
class  PXDMappingLookup
 Class to make the mapping between u/v cell ID of pixels back to DCD drain lines, pixel row/col, DCD and Switcher IDs Details: Belle Note: BELLE2-NOTE-TE-2015-01 "The vertex detector numbering scheme" https://docs.belle2.org/record/243/files/Belle%20II%20note%200010.pdf PXD WhiteBook 3.1.3 Sensor Design and Appendix #3 https://confluence.desy.de/display/BI/PXD+WebHome?preview=/34029260/56330158/PXDwb.pdf. More...
 
struct  dhc_frame_header_word0
 DHC frame header word data struct. More...
 
struct  dhc_start_frame
 DHC start frame data struct. More...
 
struct  dhc_dhe_start_frame
 DHH start frame data struct. More...
 
struct  dhc_commode_frame
 DHH common mode frame data struct. More...
 
struct  dhc_direct_readout_frame
 DHC direct readout frame data struct. More...
 
struct  dhc_direct_readout_frame_raw
 DHC RAW direct readout frame data struct. More...
 
struct  dhc_direct_readout_frame_zsd
 DHC Zero supressed direct readout frame data struct. More...
 
struct  dhc_onsen_trigger_frame
 ONSEN Trigger frame data struct. More...
 
struct  dhc_onsen_roi_frame
 ONSEN (debug) ROI frame data struct. More...
 
struct  dhc_ghost_frame
 DHC Ghost frame data struct. More...
 
struct  dhc_end_frame
 DHC End frame data struct. More...
 
struct  dhc_dhe_end_frame
 DHE End frame data struct. More...
 
class  dhc_frames
 DHC frame wrapper class. More...
 
struct  Cluster_t
 Struct to hold variables for PXD clusters. More...
 
struct  TrackPoint_t
 Struct to hold variables for intersection points. More...
 
struct  TrackCluster_t
 Struct to hold variables for track clusters. More...
 
struct  TrackBase_t
 Struct to hold variables from a track which contains a vector of data type like TrackCluster. More...
 

Typedefs

typedef std::map< Digit, DigitValueSensor
 Map of all hits in one Sensor.
 
typedef std::map< VxdID, SensorSensors
 Map of all hits in all Sensors.
 
typedef std::map< pxdClusterShapeType, std::string > pxdClusterShapeDescr
 Type specifies cluster shape type description.
 
typedef VXD::SensitiveDetector< PXDSimHit, PXDTrueHitSensitiveDetector
 The PXD Sensitive Detector class.
 
using ulittle16_t = boost::endian::little_uint16_t
 define alias ulittle16_t
 
using ulittle32_t = boost::endian::little_uint32_t
 define alias ulittle32_t
 
using ubig16_t = boost::endian::big_uint16_t
 define alias ubig16_t
 
using ubig32_t = boost::endian::big_uint32_t
 define alias ubig32_t
 
typedef boost::crc_optimal< 32, 0x04C11DB7, 0, 0, false, false > dhc_crc_32_type
 define our CRC function
 
typedef TrackBase_t< TrackCluster_tTrack_t
 Typedef TrackBase_t<TrackCluster_t> Track_t.
 
typedef genfit::MeasuredStateOnPlane TrackState
 Typedef TrackState (genfit::MeasuredStateOnPlane)
 

Enumerations

enum class  pxdClusterShapeType {
  no_shape_set = 0 ,
  shape_1 ,
  shape_2_u ,
  shape_2_v ,
  shape_2_uv_diag ,
  shape_2_uv_antidiag ,
  shape_N1 ,
  shape_1M ,
  shape_N2 ,
  shape_2M ,
  shape_4 ,
  shape_3_L ,
  shape_3_L_mirr_u ,
  shape_3_L_mirr_v ,
  shape_3_L_mirr_uv ,
  shape_large
}
 Type specifies cluster shape type.
 
enum  EDHPFrameHeaderDataType {
  c_RAW = 0x0 ,
  c_ZSD = 0x5
}
 Enums for DHP data modes in the DHP header. More...
 
enum  EDHCFrameHeaderDataType {
  c_DHP_RAW = 0x0 ,
  c_DHP_ZSD = 0x5 ,
  c_FCE_RAW = 0x1 ,
  c_COMMODE = 0x6 ,
  c_GHOST = 0x2 ,
  c_DHE_START = 0x3 ,
  c_DHE_END = 0x4 ,
  c_DHC_START = 0xB ,
  c_DHC_END = 0xC ,
  c_ONSEN_DHP = 0xD ,
  c_ONSEN_FCE = 0x9 ,
  c_ONSEN_ROI = 0xF ,
  c_ONSEN_TRG = 0xE ,
  c_UNUSED_7 = 0x7 ,
  c_UNUSED_8 = 0x8 ,
  c_UNUSED_A = 0xA
}
 Enums for DHC data frame types. More...
 
enum  EDHEStateMachineError {
  c_DHESM_NO_ERROR = 0x0 ,
  c_DHESM_MISS_DHP_FRM = 0x1 ,
  c_DHESM_TIMEOUT = 0x2 ,
  c_DHESM_DHP_LINKDOWN = 0x3 ,
  c_DHESM_DHP_MASKED = 0x4 ,
  c_DHESM_EVTNR_MM = 0x5 ,
  c_DHESM_DHP_SIZE_OVERFLOW = 0x6
}
 Enums for DHE DHP StateMachine Error States. More...
 

Functions

void getNumberOfBins (const std::shared_ptr< TH1I > &histo_ptr, unsigned short &nBinsU, unsigned short &nBinsV)
 Helper function to extract number of bins along u side and v side from counter histogram labels.
 
unsigned short getNumberOfSensors (const std::shared_ptr< TH1I > &histo_ptr)
 Helper function to extract number of sensors from counter histogram labels.
 
double CalculateMedian (std::vector< double > &signals)
 Helper function to calculate a median from unsorted signal vector. More...
 
double CalculateMedian (TH1 *hist)
 Helper function to calculate a median from 1D histogram.
 
double FitLandau (TH1 *hist)
 Helper function to estimate MPV from 1D histogram.
 
double FitLandau (std::vector< double > &signals)
 Helper function to calculate MPV from a vector. More...
 
 TEST (ClusterCache, FindNeighbours)
 Check that we cluster to hits next to each other but not if one is in between. More...
 
 TEST (ClusterCache, Merging)
 Test that clusters get merged if they are found to have a common pixel. More...
 
 TEST (ClusterCache, Empty)
 Check if the cluster cache is empty.
 
 TEST (ClusterCache, OutOfRange)
 Check that out_of_range exceptions are raised if the pixel is out of range.
 
double xiBeta2_L (const int Z=Z_Si, const double A=A_Si, const double rho=rho_Si, const int z=1)
 xi = (K/2)*(Z/A)*z*z*(rho*L)/beta2 in MeV More...
 
double hbarWp (const int Z=Z_Si, const double A=A_Si, const double rho=rho_Si)
 hbarWp = sqrt(rho*Z/A)*28.816 in eV More...
 
double getDeltaP (const double mom, const double length, const double mass=Const::electronMass)
 helper function to estimate the most probable energy loss for a given track length. More...
 
unsigned short getPXDModuleID (const VxdID &sensorID)
 Helper function to get DHE id like module id from VxdID.
 
VxdID getVxdIDFromPXDModuleID (const unsigned short &id)
 Helper function to get VxdID from DHE id like module iid.
 
std::shared_ptr< TrackStategetTrackStateOnModule (const VXD::SensorInfoBase &pxdSensorInfo, RecoTrack &recoTrack, double lambda=0.0)
 Helper function to get a track state on a module. More...
 
bool isCloseToBorder (int u, int v, int checkDistance)
 Helper function to check if a pixel is close to the border. More...
 
bool isDefectivePixelClose (int u, int v, int checkDistance, const VxdID &moduleID)
 Helper function to check if a defective (hot/dead) pixel is close. More...
 
bool isClusterAtUEdge (VxdID id, unsigned int umin, unsigned int umax)
 Helper function to check if one of the end pixels are at the edge of the sensor. More...
 
bool isClusterAtVEdge (VxdID id, unsigned int vmin, unsigned int vmax)
 Helper function to check if one of the end pixels are at the edge of the sensor. More...
 
bool isClusterAtLadderJoint (VxdID id, unsigned int vmin, unsigned int vmax)
 Helper function to check if one of the end pixels are at the ladder joint. More...
 
unsigned short getBinU (VxdID id, unsigned int uid, unsigned int vid, unsigned short nBinsU)
 Function to return a bin number for equal sized binning in U. More...
 
unsigned short getBinV (VxdID id, unsigned int vid, unsigned short nBinsV)
 Function to return a bin number for equal sized binning in V. More...
 

Variables

geometry::CreatorFactory< GeoPXDCreatorGeoPXDFactory ("PXDCreator")
 Register the creator.
 
const int Z_Si = 14
 Const and Const expressions Only valid when g_mol is the default unit. More...
 
const double A_Si = 28.085
 Atomic mass of silicon in g mol^-1.
 
const double rho_Si = 2.3290 * Unit::g_cm3
 Silicon density in g cm^-3.
 

Detailed Description

Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.

Enumeration Type Documentation

◆ EDHCFrameHeaderDataType

Enums for DHC data frame types.

4 bits value; found in the first header word of each frame. See Data format definitions [BELLE2-NOTE-TE-2016-009] on https://docs.belle2.org/

Definition at line 28 of file PXDRawDataDefinitions.h.

28  {
29  c_DHP_RAW = 0x0, // DHP memory dump ("pedestals")
30  c_DHP_ZSD = 0x5, // DHP zero supressed data
31  c_FCE_RAW = 0x1, // Clustered data
32  c_COMMODE = 0x6, // Common mode data
33  c_GHOST = 0x2, // Ghost frame, no data
34  // DHE envelope
35  c_DHE_START = 0x3, // DHE Start
36  c_DHE_END = 0x4, // DHE End
37  // DHC envelope
38  c_DHC_START = 0xB, // DHC Start
39  c_DHC_END = 0xC, // DHC End
40  // Onsen processed data, new
41  c_ONSEN_DHP = 0xD, // Onsen processed zero supressed DHP
42  c_ONSEN_FCE = 0x9, // Onsen processed clustered
43  c_ONSEN_ROI = 0xF, // Onsen ROIs (HLT+DATCON)
44  c_ONSEN_TRG = 0xE, // Trigger frame (the 1st frame)
45  // Free IDs are 0x7 0x8 0xA
46  c_UNUSED_7 = 0x7,
47  c_UNUSED_8 = 0x8,
48  c_UNUSED_A = 0xA,
49  };

◆ EDHEStateMachineError

Enums for DHE DHP StateMachine Error States.

4 bits value; Currently the same definitions as the ones encoded in ghost frame (new! but 3 bit only) See Data format definitions [BELLE2-NOTE-TE-2016-009] on https://docs.belle2.org/

Definition at line 56 of file PXDRawDataDefinitions.h.

◆ EDHPFrameHeaderDataType

Enums for DHP data modes in the DHP header.

DHP modes have the same value as for the DHC/DHE frame See Data format definitions [BELLE2-NOTE-TE-2016-009] on https://docs.belle2.org/

Definition at line 22 of file PXDRawDataDefinitions.h.

Function Documentation

◆ CalculateMedian()

double CalculateMedian ( std::vector< double > &  signals)

Helper function to calculate a median from unsorted signal vector.

The input vector gets sorted.

Definition at line 88 of file PXDCalibrationUtilities.cc.

89  {
90  auto size = signals.size();
91 
92  if (size == 0) {
93  return 0.0; // Undefined, really.
94  } else if (size <= 100) {
95  // sort() or partial_sort is in O(NlogN)
96  sort(signals.begin(), signals.end());
97  if (size % 2 == 0) {
98  return (signals[size / 2 - 1] + signals[size / 2]) / 2;
99  } else {
100  return signals[size / 2];
101  }
102  } else {
103  // nth_element or max_element in O(N) only
104  // All elements before the nth are guanranteed smaller
105  auto n = size / 2;
106  nth_element(signals.begin(), signals.begin() + n, signals.end());
107  auto med = signals[n];
108  if (!(size & 1)) { // if size is even
109  auto max_it = max_element(signals.begin(), signals.begin() + n);
110  med = (*max_it + med) / 2.0;
111  }
112  return med;
113  }
114  }

◆ FitLandau()

double FitLandau ( std::vector< double > &  signals)

Helper function to calculate MPV from a vector.

The input vector gets sorted.

Definition at line 157 of file PXDCalibrationUtilities.cc.

◆ getBinU()

unsigned short Belle2::PXD::getBinU ( VxdID  id,
unsigned int  uid,
unsigned int  vid,
unsigned short  nBinsU 
)
inline

Function to return a bin number for equal sized binning in U.

Parameters
idVxdID of the sensor
uiduID of the pixel
vidvID of the pixel
nBinsUnumber of bins in U within a sensor

Definition at line 162 of file PXDUtilities.h.

163  {
164  unsigned int drainsPerBin = 4 * Belle2::VXD::GeoCache::getInstance().get(id).getUCells() / nBinsU;
165  return (uid * 4 + vid % 4) / drainsPerBin;
166  }
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:214
static const SensorInfoBase & get(Belle2::VxdID id)
Return a reference to the SensorInfo of a given SensorID.
Definition: GeoCache.h:139
int getUCells() const
Return number of pixel/strips in u direction.

◆ getBinV()

unsigned short Belle2::PXD::getBinV ( VxdID  id,
unsigned int  vid,
unsigned short  nBinsV 
)
inline

Function to return a bin number for equal sized binning in V.

Parameters
idVxdID of the sensor
vidvID of the pixel
nBinsVnumber of bins in V within a sensor

Definition at line 172 of file PXDUtilities.h.

◆ getDeltaP()

double Belle2::PXD::getDeltaP ( const double  mom,
const double  length,
const double  mass = Const::electronMass 
)
inline

helper function to estimate the most probable energy loss for a given track length.

Parameters
momMagnitude of the momentum
lengthTrack path length
massMass of the incident particle, using e- as default
Returns
the most probable energy

Definition at line 69 of file PXDUtilities.h.

◆ getTrackStateOnModule()

std::shared_ptr< TrackState > getTrackStateOnModule ( const VXD::SensorInfoBase pxdSensorInfo,
RecoTrack recoTrack,
double  lambda = 0.0 
)

Helper function to get a track state on a module.

Parameters
pxdSensorInfoof the PXD module intersecting with the track.
recoTrackthe recoTrack to be extrapolated.
lambdathe extrapolation length from track POCA.
Returns
the shared pointer of the intersection track state on the module.

Definition at line 21 of file PXDUtilities.cc.

24  {
25  // get sensor plane, always enable alignment.
26  auto centerP = pxdSensorInfo.pointToGlobal(ROOT::Math::XYZVector(0, 0, 0), true);
27  auto normalV = pxdSensorInfo.vectorToGlobal(ROOT::Math::XYZVector(0, 0, 1), true);
28  genfit::SharedPlanePtr sensorPlaneSptr(new genfit::DetPlane(XYZToTVector(centerP), XYZToTVector(normalV)));
29 
30  // genfit track and measured state on plane
31  const genfit::Track& gfTrack = RecoTrackGenfitAccess::getGenfitTrack(recoTrack);
32  auto statePtr = std::make_shared<TrackState>();
33 
34  try {
35  *statePtr = gfTrack.getFittedState();
36  lambda = statePtr->extrapolateToPlane(sensorPlaneSptr);
37  } catch (...) {
38  B2DEBUG(20, "extrapolation to plane failed! Lambda = " << lambda);
39  return std::shared_ptr<TrackState>(nullptr);
40  }
41  auto intersec = pxdSensorInfo.pointToLocal(ROOT::Math::XYZVector(statePtr->getPos()), true);
42 
43  // check if the intersection is inside (no tolerance).
44  double tolerance = 0.0;
45  bool inside = pxdSensorInfo.inside(intersec.X(), intersec.Y(), tolerance, tolerance);
46  if (!inside) return std::shared_ptr<TrackState>(nullptr);
47 
48  return statePtr;
49  }
Detector plane.
Definition: DetPlane.h:59
Collection of TrackPoint objects, AbsTrackRep objects and FitStatus objects.
Definition: Track.h:71
const MeasuredStateOnPlane & getFittedState(int id=0, const AbsTrackRep *rep=nullptr, bool biased=true) const
Shortcut to get FittedStates.
Definition: Track.cc:294
static constexpr auto XYZToTVector
Helper function to convert XYZVector to TVector3.
Definition: VectorUtil.h:24
std::shared_ptr< genfit::DetPlane > SharedPlanePtr
Shared Pointer to a DetPlane.

◆ hbarWp()

double Belle2::PXD::hbarWp ( const int  Z = Z_Si,
const double  A = A_Si,
const double  rho = rho_Si 
)
inline

hbarWp = sqrt(rho*Z/A)*28.816 in eV

Parameters
ZAtomic number of absorber
AAtomic mass of absorber in g*mol^{-1}
rhoDensity of the absorber in g*cm^{-3}
Returns
plasma energy

Definition at line 56 of file PXDUtilities.h.

◆ isCloseToBorder()

bool isCloseToBorder ( int  u,
int  v,
int  checkDistance 
)

Helper function to check if a pixel is close to the border.

Parameters
uuID of the pixel of interest
vvID of the pixel of interest
checkDistancethe distance along u/v
Returns
true if the pixel is close to border

Definition at line 51 of file PXDUtilities.cc.

◆ isClusterAtLadderJoint()

bool Belle2::PXD::isClusterAtLadderJoint ( VxdID  id,
unsigned int  vmin,
unsigned int  vmax 
)
inline

Helper function to check if one of the end pixels are at the ladder joint.

Parameters
idVxdID of the sensor
vminminimum vID within the cluster
vmaxmaximum vID within the cluster
Returns
true if one of the end pixels are at the ladder joint

Definition at line 149 of file PXDUtilities.h.

◆ isClusterAtUEdge()

bool Belle2::PXD::isClusterAtUEdge ( VxdID  id,
unsigned int  umin,
unsigned int  umax 
)
inline

Helper function to check if one of the end pixels are at the edge of the sensor.

Parameters
idVxdID of the sensor
uminminimum uID within the cluster
umaxmaximum uID within the cluster
Returns
true if one of the end pixels are at the edge of the sensor

Definition at line 126 of file PXDUtilities.h.

◆ isClusterAtVEdge()

bool Belle2::PXD::isClusterAtVEdge ( VxdID  id,
unsigned int  vmin,
unsigned int  vmax 
)
inline

Helper function to check if one of the end pixels are at the edge of the sensor.

Parameters
idVxdID of the sensor
vminminimum vID within the cluster
vmaxmaximum vID within the cluster
Returns
true if one of the end pixels are at the edge of the sensor

Definition at line 137 of file PXDUtilities.h.

◆ isDefectivePixelClose()

bool isDefectivePixelClose ( int  u,
int  v,
int  checkDistance,
const VxdID moduleID 
)

Helper function to check if a defective (hot/dead) pixel is close.

Parameters
uuID of the pixel of interest
vvID of the pixel of interest
checkDistancethe distance along u/v
moduleIDVxdID of the module
Returns
true if a defective pixel is found in matrix

Definition at line 61 of file PXDUtilities.cc.

◆ TEST() [1/2]

Belle2::PXD::TEST ( ClusterCache  ,
FindNeighbours   
)

Check that we cluster to hits next to each other but not if one is in between.

So we set a cluster at position (2,0) marked c and the findCluster should return this cluster for all positions marked 1. All other pixels should return NULL. We do not check (0,0) and (1,0) or a row above C since the cluster cache is meant to be used with ordered data so these should already be done once we set (2,0).

* u 0 1 2 3 4 5 6 7 8 9 →
* v┌─┬─┬─┬─┬─┬─┬─┬─┬─┬─┬─┐
* 0│ │ │C│1│0│ │ │ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* 1│0│1│1│1│0│ │ │ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* 3│0│0│0│0│0│ │ │ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* 4│ │ │ │ │ │ │ │ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* ↓│ │ │ │ │ │ │ │ │ │ │ │
*  └─┴─┴─┴─┴─┴─┴─┴─┴─┴─┴─┘
* 

Definition at line 45 of file ClusterCache.cc.

46  {
47  ClusterCache cache;
48  for (int v = 0; v < 4; ++v) {
49  for (int u = 0; u < 5; ++u) {
50  if (v == 0 && u <= 2) continue;
51  cache.clear();
52  ClusterCandidate& cls1 = cache.findCluster(2, 0);
53  if ((v == 0 && u == 3) || (v == 1 && u >= 1 && u <= 3)) {
54  //Check that neighboring pixels return the same cluster
55  EXPECT_EQ(&cls1, &cache.findCluster(u, v)) << "u: " << u << " v: " << v;
56  } else {
57  //And all other pixels return another cluster
58  EXPECT_NE(&cls1, &cache.findCluster(u, v)) << "u: " << u << " v: " << v;
59  }
60  }
61  }
62  }

◆ TEST() [2/2]

Belle2::PXD::TEST ( ClusterCache  ,
Merging   
)

Test that clusters get merged if they are found to have a common pixel.

If we have to clusters called 1 and 2 like shown below and we add a pixel at X than all occurences of 1 and 2 should be merged to one of those two. It does not matter who is merged to whom but after adding X the result must be

  • one of the two clusters is empty
  • the other cluster contains all pixels
  • all pointers for 1 and 2 should point to the same pixel so that when we add pixels Y and Z we get the same pointer as for X
* u 0 1 2 3 4 5 6 7 8 9 →
* v┌─┬─┬─┬─┬─┬─┬─┬─┬─┬─┬─┐
* 0│ │ │1│ │2│ │3│ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* 1│4│ │ │X│ │ │5│ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* 2│ │6│Y│ │7│Z│ │ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* 3│ │ │ │ │ │ │ │ │ │ │ │
*  ├─┼─┼─┼─┼─┼─┼─┼─┼─┼─┼─┤
* ↓│ │ │ │ │ │ │ │ │ │ │ │
*  └─┴─┴─┴─┴─┴─┴─┴─┴─┴─┴─┘
* 

Definition at line 90 of file ClusterCache.cc.

◆ xiBeta2_L()

double Belle2::PXD::xiBeta2_L ( const int  Z = Z_Si,
const double  A = A_Si,
const double  rho = rho_Si,
const int  z = 1 
)
inline

xi = (K/2)*(Z/A)*z*z*(rho*L)/beta2 in MeV

Parameters
ZAtomic number of absorber
AAtomic mass of absorber in g*mol^{-1}
rhoDensity of the absorber in g*cm^{-3}
zCharge number of incident particle
Returns
xi*beta^2/L in MeV/cm where L is track length

Definition at line 41 of file PXDUtilities.h.

Variable Documentation

◆ Z_Si

const int Z_Si = 14

Const and Const expressions Only valid when g_mol is the default unit.

Atomic number of silicon

Definition at line 30 of file PXDUtilities.h.