FairRoot/PandaRoot
Public Types | Public Member Functions | Protected Member Functions | Protected Attributes | Private Attributes | List of all members
DiscDIRC_Photodetector Class Reference

#include <PndDiscPhotodetector.h>

Inheritance diagram for DiscDIRC_Photodetector:
SensorGrid::SensorGridPhotodetector

Public Types

enum  DesignID { DESIGN_SIPM = 1, DESIGN_LRD }
 

Public Member Functions

 DiscDIRC_Photodetector ()
 
 DiscDIRC_Photodetector (DesignID id)
 
virtual ~DiscDIRC_Photodetector ()
 
void SetPDE (int n_entries, const double *wavelength_nm, const double *pde)
 
virtual double GetPDE (const double &wavelength_nm) const
 Derived classes should override this function and return the average pde of the sensor. More...
 
void SetPixel (int pixel_id, double efficiency, double noise_rate, double time_res_ns)
 
void UseInhomogenityFactor (bool flag)
 
int Detect (double const &hit_pos_x, double const &hit_pos_y, double const &hit_time_ns, double const &wavelength_nm, PixelInfo &pixel_info, double &smeared_time_ns) const
 Handle photon detection: More...
 
void GenerateNoise (double const &time_start_ns, double const &time_window_ns, std::vector< std::pair< int, double > > &hits) const
 Generate noise hits. More...
 
int GenerateNoise (double const &time_start_ns, double const &time_window_ns)
 Reinit noise iterator and return number of noise hits. More...
 
virtual double GetInhomegenityFactor (double const &, double const &) const
 Derived classes should override this function to account for spatial PDE deviations (e.g. due to gain inhomgenity) More...
 
double GetSmearedTime (double const &time_value, PixelInfo const &pixel_info) const
 Apply time smearing. More...
 
bool GetNoiseHit (int &pixel_number, double &hit_time_ns, double &smeared_time_ns)
 Noise hit generator function which has to be called after GenerateNoise() to retrieve the noise hits. More...
 
void SetDCR (double const &dcr_Hz)
 Set the dark count rate for all pixels to dcr_Hz. More...
 
const SensorGridBase * GetGrid ()
 

Protected Member Functions

void Init ()
 
void Init (DesignID design_id)
 
void Init (SensorGridBase *sensor_grid_, bool per_pixel_traits_, double const &efficiency_init, double const &noise_rate_init, double const &time_res_init)
 Initialization - can also be used in derived class ctors (avoid code duplication). More...
 

Protected Attributes

bool per_pixel_traits
 
bool use_inhomogenity_factor
 
PixelTraitspixel_traits
 
SensorGridBase * sensor_grid
 
int noisegen_pixel_number
 
double noisegen_time_start_ns
 
double noisegen_time_window_ns
 
int noisegen_n_hits
 
int noisegen_current_hit
 
double noisegen_current_time_sigma
 

Private Attributes

ROOT::Math::Interpolator pde_interpolator
 

Detailed Description

Definition at line 16 of file PndDiscPhotodetector.h.

Member Enumeration Documentation

Enumerator
DESIGN_SIPM 
DESIGN_LRD 

Definition at line 19 of file PndDiscPhotodetector.h.

Constructor & Destructor Documentation

DiscDIRC_Photodetector::DiscDIRC_Photodetector ( )
DiscDIRC_Photodetector::DiscDIRC_Photodetector ( DesignID  id)

Definition at line 21 of file PndDiscPhotodetector.cxx.

References Init().

21  {
22  Init(design_id);
23 }
virtual DiscDIRC_Photodetector::~DiscDIRC_Photodetector ( )
inlinevirtual

Definition at line 26 of file PndDiscPhotodetector.h.

26 {};

Member Function Documentation

int SensorGrid::SensorGridPhotodetector::Detect ( double const &  hit_pos_x,
double const &  hit_pos_y,
double const &  hit_time_ns,
double const &  wavelength_nm,
PixelInfo pixel_info,
double &  smeared_time_ns 
) const
inherited

Handle photon detection:

Detection of a hit.

Definition at line 62 of file PndDiscSensorGridPhotodetector.cxx.

References SensorGrid::PixelInfo::column_on_grid, SensorGrid::SensorGridPhotodetector::PixelTraits::efficiency, SensorGrid::SensorGridPhotodetector::GetInhomegenityFactor(), SensorGrid::SensorGridPhotodetector::GetPDE(), SensorGrid::PixelInfo::grid_id, SensorGrid::SensorGridPhotodetector::per_pixel_traits, SensorGrid::PixelInfo::pixel_number, SensorGrid::SensorGridPhotodetector::pixel_traits, SensorGrid::SensorGridBase::PositionToPixel(), SensorGrid::PixelInfo::row_on_grid, SensorGrid::SensorGridPhotodetector::sensor_grid, SensorGrid::SensorGridPhotodetector::PixelTraits::time_res_ns, and SensorGrid::SensorGridPhotodetector::use_inhomogenity_factor.

Referenced by PndDiscTaskDigitization::Exec().

62  {
63 #ifdef DETECT_DEBUG
64  ofs_detect_debug << "X: " << hit_pos_x << " Y: " << hit_pos_y ;
65 #endif
66  if(!sensor_grid->PositionToPixel(hit_pos_x, hit_pos_y, pixel_info)) {
67 #ifdef DETECT_DEBUG
68  ofs_detect_debug << " not detected\n";
69 #endif
70  return -1; // no pixel (geometric efficiency)
71  }
72 #ifdef DETECT_DEBUG
73  ofs_detect_debug << "Pixel: " << pixel_info.pixel_number << " Grid: " << pixel_info.grid_id << " Column/row: " << pixel_info.column_on_grid << " / " << pixel_info.row_on_grid << std::endl;
74 #endif
75 
76  double global_pde = GetPDE(wavelength_nm);
77  if(use_inhomogenity_factor) global_pde*=GetInhomegenityFactor(hit_pos_x, hit_pos_y);
78 
79  if(per_pixel_traits) {
80  if(gRandom->Uniform() > pixel_traits[pixel_info.pixel_number].efficiency*global_pde) return -1; // global PDE and per pixel efficiency (e.g. due to bad gain uniformity)
81  smeared_time_ns = gRandom->Gaus(hit_time_ns, pixel_traits[pixel_info.pixel_number].time_res_ns);
82  } else {
83  if(gRandom->Uniform() > pixel_traits->efficiency*global_pde) return -1;
84  smeared_time_ns = gRandom->Gaus(hit_time_ns, pixel_traits->time_res_ns);
85  }
86  return pixel_info.pixel_number;
87 }
virtual bool PositionToPixel(const double &x, const double &y, PixelInfo &pixel_info) const =0
virtual double GetPDE(const double &) const
Derived classes should override this function and return the average pde of the sensor.
virtual double GetInhomegenityFactor(double const &, double const &) const
Derived classes should override this function to account for spatial PDE deviations (e...
void SensorGrid::SensorGridPhotodetector::GenerateNoise ( double const &  time_start_ns,
double const &  time_window_ns,
std::vector< std::pair< int, double > > &  hits 
) const
inherited

Generate noise hits.

Generate noise on the sensor.

Definition at line 91 of file PndDiscSensorGridPhotodetector.cxx.

References SensorGrid::SensorGridBase::GetNumberOfPixels(), hit(), hits, i, SensorGrid::SensorGridPhotodetector::PixelTraits::noise_rate, SensorGrid::SensorGridPhotodetector::per_pixel_traits, SensorGrid::SensorGridPhotodetector::pixel_traits, and SensorGrid::SensorGridPhotodetector::sensor_grid.

91  {
92  int i, n_hits, n_pixels = sensor_grid->GetNumberOfPixels();
93  double avg_n_hits;
94  std::pair<int, double> hit;
95 
96  // loop over all pixel numbers
97  for(hit.first = 0; hit.first < n_pixels; ++hit.first) {
98  // average number of noise events per pixel:
100  avg_n_hits = pixel_traits[hit.first].noise_rate*time_window_ns*1E-9;
101  else
102  avg_n_hits = pixel_traits->noise_rate*time_window_ns*1E-9;
103 
104  n_hits = gRandom->Poisson(avg_n_hits); // compute number of hits on this pixel
105  // generate the hits in the given time window:
106  for(i=0; i<n_hits; ++i) {
107  hit.second = time_start_ns + gRandom->Uniform(time_window_ns); // distribute hit times in window
108  hits.push_back(hit); // append hits to vector
109  }
110  }
111 }
Int_t i
Definition: run_full.C:25
int hit(Int_t nEvents=0, TString inFile="sim.root", TString parFile="par.root", TString inDigi="digi.root", TString outFile="hit.root", Int_t timeBased=0)
Definition: hit.C:1
CbmHit * hits[nHits]
Definition: RiemannTest.C:19
int SensorGrid::SensorGridPhotodetector::GenerateNoise ( double const &  time_start_ns,
double const &  time_window_ns 
)
inherited
const SensorGridBase* SensorGrid::SensorGridPhotodetector::GetGrid ( )
inlineinherited
virtual double SensorGrid::SensorGridPhotodetector::GetInhomegenityFactor ( double const &  ,
double const &   
) const
inlinevirtualinherited

Derived classes should override this function to account for spatial PDE deviations (e.g. due to gain inhomgenity)

This function is used to implement inhomogenitites over the whole sensor, which cannot be treated as single pixel efficiencies, e.g. if the pixels are large compared to the spatial deviation of the inhomogenity.

Definition at line 56 of file PndDiscSensorGridPhotodetector.h.

Referenced by SensorGrid::SensorGridPhotodetector::Detect().

56 {return 1.0;} //FIXME//hit_pos_x hit_pos_y//[R.K.03/2017] unused variable(s)
bool SensorGrid::SensorGridPhotodetector::GetNoiseHit ( int &  pixel_number,
double &  hit_time_ns,
double &  smeared_time_ns 
)
inherited

Noise hit generator function which has to be called after GenerateNoise() to retrieve the noise hits.

Definition at line 127 of file PndDiscSensorGridPhotodetector.cxx.

References SensorGrid::SensorGridBase::GetNumberOfPixels(), SensorGrid::SensorGridPhotodetector::PixelTraits::noise_rate, SensorGrid::SensorGridPhotodetector::noisegen_current_hit, SensorGrid::SensorGridPhotodetector::noisegen_current_time_sigma, SensorGrid::SensorGridPhotodetector::noisegen_n_hits, SensorGrid::SensorGridPhotodetector::noisegen_pixel_number, SensorGrid::SensorGridPhotodetector::noisegen_time_start_ns, SensorGrid::SensorGridPhotodetector::noisegen_time_window_ns, SensorGrid::SensorGridPhotodetector::per_pixel_traits, SensorGrid::SensorGridPhotodetector::pixel_traits, SensorGrid::SensorGridPhotodetector::sensor_grid, and SensorGrid::SensorGridPhotodetector::PixelTraits::time_res_ns.

127  {
128  if(noisegen_pixel_number >= sensor_grid->GetNumberOfPixels()) return false;
129 
130  // after iterating over all hits, switch to next pixel.
132  do{
133  // next pixel:
135  if(noisegen_pixel_number >= sensor_grid->GetNumberOfPixels()) return false;
136 
137  // average number of noise events on this pixel:
138  double avg_n_hits = 0.0;
139  if(per_pixel_traits) {
140  avg_n_hits = pixel_traits[pixel_number].noise_rate*noisegen_time_window_ns*1E-9;
142  }
143  else {
146  }
147 
148  // number of hits on this pixel:
149  noisegen_n_hits = gRandom->Poisson(avg_n_hits); // compute number of hits on this pixel
150  }
151  while(noisegen_n_hits == 0);
152 
154  LOG(INFO) << Form("Number of generated hits on pixel %d are: %d\n", noisegen_pixel_number, noisegen_n_hits);
155  }
156 
157  pixel_number = noisegen_pixel_number;
158  hit_time_ns = noisegen_time_start_ns + gRandom->Uniform(noisegen_time_window_ns); // distribute hit times in window
159  smeared_time_ns = gRandom->Gaus(hit_time_ns, noisegen_current_time_sigma);
161 
162  return true;
163 }
virtual double DiscDIRC_Photodetector::GetPDE ( const double &  ) const
inlinevirtual

Derived classes should override this function and return the average pde of the sensor.

Per pixel deviations of the PDE can be taken into account by including them in a individual pixel efficiency.

Reimplemented from SensorGrid::SensorGridPhotodetector.

Definition at line 33 of file PndDiscPhotodetector.h.

References pde_interpolator.

34  {
35  return pde_interpolator.Eval(wavelength_nm);
36  }
ROOT::Math::Interpolator pde_interpolator
double SensorGrid::SensorGridPhotodetector::GetSmearedTime ( double const &  hit_time_ns,
PixelInfo const &  pixel_info 
) const
inherited

Apply time smearing.

Can be used to smear the time of noise entries. Noise is returned with its true time stamp for dead-time handling and an additional smeared time would increase the array, so it appears to me that it is more handy to smear afterwards. A better way might be a noise iterator which generates one hit at a time (-> implemented, will be tested).

Definition at line 171 of file PndDiscSensorGridPhotodetector.cxx.

References SensorGrid::SensorGridPhotodetector::per_pixel_traits, SensorGrid::PixelInfo::pixel_number, SensorGrid::SensorGridPhotodetector::pixel_traits, and SensorGrid::SensorGridPhotodetector::PixelTraits::time_res_ns.

171  {
172  double smeared_time_ns;
173  if(per_pixel_traits) {
174  smeared_time_ns = gRandom->Gaus(hit_time_ns, pixel_traits[pixel_info.pixel_number].time_res_ns);
175  } else {
176  smeared_time_ns = gRandom->Gaus(hit_time_ns, pixel_traits->time_res_ns);
177  }
178  return smeared_time_ns;
179 }
void SensorGrid::SensorGridPhotodetector::Init ( SensorGridBase sensor_grid_,
bool  per_pixel_traits_,
double const &  efficiency_init,
double const &  noise_rate_init,
double const &  time_res_init 
)
protectedinherited

Initialization - can also be used in derived class ctors (avoid code duplication).

Definition at line 34 of file PndDiscSensorGridPhotodetector.cxx.

References SensorGrid::SensorGridPhotodetector::PixelTraits::efficiency, SensorGrid::SensorGridBase::GetNumberOfPixels(), i, SensorGrid::SensorGridBase::LockGrid(), SensorGrid::SensorGridPhotodetector::PixelTraits::noise_rate, SensorGrid::SensorGridPhotodetector::per_pixel_traits, SensorGrid::SensorGridPhotodetector::pixel_traits, SensorGrid::SensorGridPhotodetector::sensor_grid, and SensorGrid::SensorGridPhotodetector::PixelTraits::time_res_ns.

Referenced by Init(), and SensorGrid::SensorGridPhotodetector::SensorGridPhotodetector().

34  {
35  sensor_grid = sensor_grid_;
36  sensor_grid->LockGrid(true);
37  per_pixel_traits = per_pixel_traits_;
38  if(per_pixel_traits) {
39  pixel_traits = new PixelTraits[sensor_grid->GetNumberOfPixels()];
40  //if(pixel_traits == NULL)
41  // throw ...
42  for(int i = 0; i< sensor_grid->GetNumberOfPixels(); ++i) {
43  pixel_traits[i].efficiency = efficiency_init;
44  pixel_traits[i].noise_rate = noise_rate_init;
45  pixel_traits[i].time_res_ns = time_res_init;
46  }
47  } else {
48  pixel_traits = new PixelTraits();
49  pixel_traits->efficiency = efficiency_init;
50  pixel_traits->noise_rate = noise_rate_init;
51  pixel_traits->time_res_ns = time_res_init;
52  }
53 }
Int_t i
Definition: run_full.C:25
void LockGrid(bool lock)
Lock the grid:
void DiscDIRC_Photodetector::Init ( )
protected

Referenced by DiscDIRC_Photodetector().

void DiscDIRC_Photodetector::Init ( DesignID  design_id)
protected

Definition at line 26 of file PndDiscPhotodetector.cxx.

References SensorGrid::MultipleGrids::AddGrid(), DESIGN_LRD, DESIGN_SIPM, SensorGrid::SensorGridPhotodetector::Init(), SensorGrid::SensorGridBase::SetUserColumnOffset(), and SensorGrid::SensorGridBase::SetUserRowOffset().

26  {
27  using namespace SensorGrid;
28 
29  if(design_id == DESIGN_SIPM) {
30  // 210 stripes SiPM sensor (acceptance gap due to bonding is covered with stripes !).
31  MultipleGrids * new_sensor_grid = new MultipleGrids();
32  if(new_sensor_grid == 0) LOG(FATAL) << "Allocation of sensor failed.";
33 
34  BasicGrid * grid = NULL;
35 
36  grid = new BasicGrid( -1.5900/2.0, 1.5900/2.0, 1.0, 1 ,
37  -0.672, 0.0064, 0.0064, 105 );
38  if(grid == 0) LOG(FATAL) << "Allocation of sensor failed.";
39 
40  new_sensor_grid->AddGrid(grid);
41 
42 
43  grid = new BasicGrid( 0.0, 1.5900/2.0, 1.0, 1 ,
44  -0.672, 0.0064, 0.0064, 105 );
45  if(grid == 0) LOG(FATAL) << "Allocation of sensor failed.";
46 
47  grid->SetUserColumnOffset(1);
48  new_sensor_grid->AddGrid(grid);
49 
50 
51  grid = new BasicGrid( -1.5900/2.0, 1.5900/2.0, 1.0, 1 ,
52  0.0, 0.0064, 0.0064, 105 );
53  if(grid == 0) LOG(FATAL) << "Allocation of sensor failed.";
54 
55  grid->SetUserRowOffset(105);
56  new_sensor_grid->AddGrid(grid);
57 
58 
59  grid = new BasicGrid( 0.0, 1.5900/2.0, 1.0, 1 ,
60  0.0, 0.0064, 0.0064, 105 );
61  if(grid == 0) LOG(FATAL) << "Allocation of sensor failed.";
62 
63  grid->SetUserColumnOffset(1);
64  grid->SetUserRowOffset(105);
65  new_sensor_grid->AddGrid(grid);
66 
67  double active_cell_fraction = 0.9;
68  double SPAD_DCR = 25.0; // Hz
69  SensorGrid::SensorGridPhotodetector::Init(new_sensor_grid, false, active_cell_fraction, 256. * SPAD_DCR, 0.060);
70  }
71  else if(design_id == DESIGN_LRD) {
72  // Implement Low Rate Design ...
73  BasicGrid * grid = NULL;
74 
75  // digitization grid 3 x 100 for 50 x 50 mm^2 active area (Planacon)
76  grid = new BasicGrid( -2.55, 1.6, 1.7, 3 , // X: 1 mm gap between the 16 mm wide stripes
77  -2.5, 0.05, 0.05, 100 ); // Y: 100 stripes without gap
78  if(grid == 0) LOG(FATAL) << "Allocation of sensor failed.";
79 
80  SensorGrid::SensorGridPhotodetector::Init(grid, false, 0.7, 1.0, 0.050);
81  }
82 }
A generic regular pixel grid with dead space between cells.
A grid to group other grids or to create nested grids.
void Init(SensorGridBase *sensor_grid_, bool per_pixel_traits_, double const &efficiency_init, double const &noise_rate_init, double const &time_res_init)
Initialization - can also be used in derived class ctors (avoid code duplication).
void AddGrid(SensorGridBase *grid)
void SensorGrid::SensorGridPhotodetector::SetDCR ( double const &  dcr_Hz)
inherited

Set the dark count rate for all pixels to dcr_Hz.

Will set the DCR for all pixels to dcr_Hz. The current implementation does not provide the functionality to model statistical fluctuations in the DCR.

Definition at line 185 of file PndDiscSensorGridPhotodetector.cxx.

References SensorGrid::SensorGridBase::GetNumberOfPixels(), i, SensorGrid::SensorGridPhotodetector::PixelTraits::noise_rate, SensorGrid::SensorGridPhotodetector::per_pixel_traits, SensorGrid::SensorGridPhotodetector::pixel_traits, and SensorGrid::SensorGridPhotodetector::sensor_grid.

185  {
186  if(per_pixel_traits) {
187  int n_pixels = sensor_grid->GetNumberOfPixels();
188  for(int i=0; i<n_pixels; ++i) {
189  pixel_traits[i].noise_rate = dcr_Hz;
190  }
191  }
192  else
193  pixel_traits->noise_rate = dcr_Hz;
194 }
Int_t i
Definition: run_full.C:25
void DiscDIRC_Photodetector::SetPDE ( int  n_entries,
const double *  wavelength_nm,
const double *  pde 
)
inline

Definition at line 28 of file PndDiscPhotodetector.h.

References pde_interpolator.

Referenced by PndDiscTaskDigitization::Init().

29  {
30  pde_interpolator.SetData(n_entries, wavelength_nm, pde);
31  }
ROOT::Math::Interpolator pde_interpolator
void SensorGrid::SensorGridPhotodetector::SetPixel ( int  pixel_id,
double  efficiency,
double  noise_rate,
double  time_res_ns 
)
inherited
void SensorGrid::SensorGridPhotodetector::UseInhomogenityFactor ( bool  flag)
inlineinherited

Member Data Documentation

int SensorGrid::SensorGridPhotodetector::noisegen_current_hit
protectedinherited
double SensorGrid::SensorGridPhotodetector::noisegen_current_time_sigma
protectedinherited
int SensorGrid::SensorGridPhotodetector::noisegen_n_hits
protectedinherited
int SensorGrid::SensorGridPhotodetector::noisegen_pixel_number
protectedinherited
double SensorGrid::SensorGridPhotodetector::noisegen_time_start_ns
protectedinherited
double SensorGrid::SensorGridPhotodetector::noisegen_time_window_ns
protectedinherited
ROOT::Math::Interpolator DiscDIRC_Photodetector::pde_interpolator
private

Definition at line 43 of file PndDiscPhotodetector.h.

Referenced by GetPDE(), and SetPDE().

bool SensorGrid::SensorGridPhotodetector::per_pixel_traits
protectedinherited
PixelTraits* SensorGrid::SensorGridPhotodetector::pixel_traits
protectedinherited
SensorGridBase* SensorGrid::SensorGridPhotodetector::sensor_grid
protectedinherited
bool SensorGrid::SensorGridPhotodetector::use_inhomogenity_factor
protectedinherited

The documentation for this class was generated from the following files: