Belle II Software  release-05-01-25
TOPreco Class Reference

TOP reconstruction: this class provides interface to fortran code. More...

#include <TOPreco.h>

Public Types

enum  PDFoption {
  c_Rough = 0,
  c_Fine,
  c_Optimal
}
 Options for PDF: rough: no dependence on y fine: y dependent PDF everywhere optimal: y dependent PDF only where necessary (default)
 
enum  StoreOption {
  c_Reduced,
  c_Full
}
 Options for storing PDF parameters in Fortran common TOP_PIK reduced: only position, width, nphot and fic (default) full: also number of reflections etc.
 

Public Member Functions

 TOPreco (int NumHyp, double Masses[], int pdgCodes[], double BkgPerModule=0, double ScaleN0=1)
 Constructor. More...
 
void setMass (double mass, int pdg)
 Set mass of the particle hypothesis (overrides settings in the constructor) More...
 
void setTimeWindow (double Tmin, double Tmax)
 Set time window for photons. More...
 
void getTimeWindow (double &Tmin, double &Tmax)
 Returns time window for photons. More...
 
void setPDFoption (PDFoption opt, int NP=0, int NC=0)
 Sets PDF option. More...
 
void setStoreOption (StoreOption opt)
 Sets option for storing PDF parameters in Fortran common TOP_PIK. More...
 
void clearData ()
 Clear data list.
 
int addData (int moduleID, int pixelID, double time, double timeError)
 Add data. More...
 
int getDataSize ()
 Return size of data list.
 
void reconstruct (const TOPtrack &trk, int pdg=0)
 Run reconstruction for a given track. More...
 
int getFlag ()
 Return status. More...
 
double getExpectedPhotons (int i=0)
 Return expected number of photons (signal + background) for i-th mass hypothesis. More...
 
double getExpectedBG ()
 Return expected number of background photons. More...
 
int getNumOfPhotons ()
 Return number of measured photons. More...
 
double getLogL (int i=0)
 Return log likelihood for i-th mass hypothesis. More...
 
void getLogL (int Size, double LogL[], double ExpNphot[], int &Nphot)
 Get all reconstruction results at once. More...
 
double getLogL (double timeShift, double timeMin, double timeMax, double sigma=0.0)
 Return log likelihood for the last mass hypothesis using time-shifted PDF If timeMax <= timeMin use those set by setTimeWindow(double Tmin, double Tmax) More...
 
void getLogL (double timeShift, double timeMin, double timeMax, double sigma, float *logL)
 Return pixel log likelihoods for the last mass hypothesis using time-shifted PDF If timeMax <= timeMin use those set by setTimeWindow(double Tmin, double Tmax) More...
 
void getTrackHit (int LocGlob, double R[3], double Dir[3], double &Len, double &Tlen, double &Mom, int &moduleID)
 Return track hit at the bar in Local or Global frame. More...
 
void dumpLogL (int NumHyp)
 Print log likelihoods to std output. More...
 
void dumpTrackHit (int LocGlob)
 Print track to std output. More...
 
int getPullSize ()
 Get pulls: size. More...
 
void getPull (int K, int &pixelID, float &T, float &T0, float &Wid, float &PhiCer, float &Wt)
 Get pulls: K-th pull. More...
 
double getPDF (int pixelID, double t, double mass, int PDG, double jitter=0)
 Return PDF for pixel pixelID at time t for mass hypothesis mass. More...
 
void setBeta (double beta)
 Set track beta (for beta resolution studies) if beta>0 this value is used instead of beta from momentum and mass. More...
 
double getBeta () const
 Return track beta. More...
 
void redoPDF (double mass, int PDG)
 Re-calculate PDF for a given particle mass using option c_Fine. More...
 
int getNumofPDFPeaks (int pixelID) const
 Returns number of peaks for given pixel describing signal PDF. More...
 
void getPDFPeak (int pixelID, int k, float &position, float &width, float &numPhotons) const
 Returns k-th PDF peak for given pixel describing signal PDF. More...
 
float getBkgLevel (int pixelID) const
 Returns estimated background level for given pixel. More...
 
int getPDFPeakType (int pixelID, int k) const
 Returns type of the k-th PDF peak for given pixel. More...
 
float getPDFPeakFic (int pixelID, int k) const
 Returns Cerenkov azimuthal angle of PDF peak. More...
 
float getPDFPeakE (int pixelID, int k) const
 Returns photon energy of PDF peak. More...
 
float getPDFPeakSigE (int pixelID, int k) const
 Returns photon energy spread of PDF peak. More...
 
int getPDFPeakNx (int pixelID, int k) const
 Returns total number of reflections in x of PDF peak. More...
 
int getPDFPeakNy (int pixelID, int k) const
 Returns total number of reflections in y of PDF peak. More...
 
int getPDFPeakNxm (int pixelID, int k) const
 Returns number of reflections in x before mirror. More...
 
int getPDFPeakNym (int pixelID, int k) const
 Returns number of reflections in y before mirror. More...
 
int getPDFPeakNxe (int pixelID, int k) const
 Returns number of reflections in x in prism. More...
 
int getPDFPeakNye (int pixelID, int k) const
 Returns number of reflections in y in prism. More...
 
float getPDFPeakXD (int pixelID, int k) const
 Returns unfolded x position of pixel. More...
 
float getPDFPeakYD (int pixelID, int k) const
 Returns unfolded y position of pixel. More...
 
float getPDFPeakKxe (int pixelID, int k) const
 Returns photon reconstructed direction in x at emission. More...
 
float getPDFPeakKye (int pixelID, int k) const
 Returns photon reconstructed direction in y at emission. More...
 
float getPDFPeakKze (int pixelID, int k) const
 Returns photon reconstructed direction in z at emission. More...
 
float getPDFPeakKxd (int pixelID, int k) const
 Returns photon reconstructed direction in x at detection. More...
 
float getPDFPeakKyd (int pixelID, int k) const
 Returns photon reconstructed direction in y at detection. More...
 
float getPDFPeakKzd (int pixelID, int k) const
 Returns photon reconstructed direction in z at detection. More...
 

Static Public Member Functions

static void setChannelMask (const DBObjPtr< TOPCalChannelMask > &mask, const TOPAsicMask &asicMask)
 Set channel mask. More...
 
static void setUncalibratedChannelsOff (const DBObjPtr< TOPCalChannelT0 > &channelT0)
 Set uncalibrated channels off. More...
 
static void setUncalibratedChannelsOff (const DBObjPtr< TOPCalTimebase > &timebase)
 Set uncalibrated channels off. More...
 
static void printChannelMask ()
 Print channel mask.
 
static void setChannelEffi ()
 Set relative efficiencies of pixels.
 

Private Member Functions

void reconstruct (double X, double Y, double Z, double Tlen, double Px, double Py, double Pz, int Q, int pdg=0, int moduleID=0)
 Run reconstruction. More...
 

Private Attributes

int m_hypID = 0
 true hypothesis ID
 
double m_beta = 0
 beta value, if set
 

Detailed Description

TOP reconstruction: this class provides interface to fortran code.

Definition at line 36 of file TOPreco.h.

Constructor & Destructor Documentation

◆ TOPreco()

TOPreco ( int  NumHyp,
double  Masses[],
int  pdgCodes[],
double  BkgPerModule = 0,
double  ScaleN0 = 1 
)

Constructor.

Parameters
NumHypnumber of mass hypotheses
Massesmasses
pdgCodesPDG codes
BkgPerModuleestimation for minimal number of background hits
ScaleN0scale factor to scale N0

Definition at line 65 of file TOPreco.cc.

67  {
68  data_clear_();
69  rtra_clear_();
70 
71  std::vector<float> masses;
72  for (int i = 0; i < Num; i++) {
73  masses.push_back((float) Masses[i]);
74  }
75  rtra_set_hypo_(&Num, masses.data());
76  rtra_set_hypid_(&Num, pdgCodes);
77 
78  float b = (float) BkgPerModule;
79  float s = (float) ScaleN0;
80  set_top_par_(&b, &s);
81 
82  setPDFoption(c_Optimal); // default option
83  setTimeWindow(0.0, 0.0); // use default (TDC range)
84  setBeta(0); // use default: beta from momentum and mass
85  }

Member Function Documentation

◆ addData()

int addData ( int  moduleID,
int  pixelID,
double  time,
double  timeError 
)

Add data.

Parameters
moduleIDmodule ID
pixelIDpixel ID (e.g. software channel, 1-based)
timeTBC and t0-corrected time in [ns]
timeErrortime uncertainty in [ns]

Definition at line 214 of file TOPreco.cc.

◆ dumpLogL()

void dumpLogL ( int  NumHyp)

Print log likelihoods to std output.

Parameters
NumHypnumber of hypotheses

Definition at line 356 of file TOPreco.cc.

◆ dumpTrackHit()

void dumpTrackHit ( int  LocGlob)

Print track to std output.

Parameters
LocGlobin Local or Global frame

Definition at line 386 of file TOPreco.cc.

◆ getBeta()

double getBeta ( ) const
inline

Return track beta.

Returns
beta value

Definition at line 290 of file TOPreco.h.

290 {return m_beta;};

◆ getBkgLevel()

float getBkgLevel ( int  pixelID) const

Returns estimated background level for given pixel.

Parameters
pixelIDpixel ID (1-based)
Returns
number of background hits per nano second

Definition at line 464 of file TOPreco.cc.

◆ getExpectedBG()

double getExpectedBG ( )

Return expected number of background photons.

Returns
expected number of background

Definition at line 291 of file TOPreco.cc.

◆ getExpectedPhotons()

double getExpectedPhotons ( int  i = 0)

Return expected number of photons (signal + background) for i-th mass hypothesis.

Parameters
iindex of mass hypothesis
Returns
expected number of photons

Definition at line 284 of file TOPreco.cc.

◆ getFlag()

int getFlag ( )

Return status.

Returns
status: 1=OK, 0=out of acceptance, -1=inside gap

Definition at line 278 of file TOPreco.cc.

◆ getLogL() [1/4]

void getLogL ( double  timeShift,
double  timeMin,
double  timeMax,
double  sigma,
float *  logL 
)

Return pixel log likelihoods for the last mass hypothesis using time-shifted PDF If timeMax <= timeMin use those set by setTimeWindow(double Tmin, double Tmax)

Parameters
timeShifttime shift of PDF
timeMinlower edge of time window within which the photons are accepted
timeMaxupper edge of time window within which the photons are accepted
sigmaadditional time smearing sigma
logLreturn array of pixel log likelihood values (must be zeroed on input)

Definition at line 332 of file TOPreco.cc.

◆ getLogL() [2/4]

double getLogL ( double  timeShift,
double  timeMin,
double  timeMax,
double  sigma = 0.0 
)

Return log likelihood for the last mass hypothesis using time-shifted PDF If timeMax <= timeMin use those set by setTimeWindow(double Tmin, double Tmax)

Parameters
timeShifttime shift of PDF
timeMinlower edge of time window within which the photons are accepted
timeMaxupper edge of time window within which the photons are accepted
sigmaadditional time smearing sigma
Returns
log likelihood

Definition at line 322 of file TOPreco.cc.

◆ getLogL() [3/4]

double getLogL ( int  i = 0)

Return log likelihood for i-th mass hypothesis.

Parameters
iindex of mass hypothesis
Returns
log likelihood

Definition at line 303 of file TOPreco.cc.

◆ getLogL() [4/4]

void getLogL ( int  Size,
double  LogL[],
double  ExpNphot[],
int &  Nphot 
)

Get all reconstruction results at once.

Parameters
Sizesize of arrays
LogLlog likelihoods for Masses[]
ExpNphotexpected number of photons for Masses[]
Nphotmeasured number of photons

Definition at line 310 of file TOPreco.cc.

◆ getNumofPDFPeaks()

int getNumofPDFPeaks ( int  pixelID) const

Returns number of peaks for given pixel describing signal PDF.

Parameters
pixelIDpixel ID (1-based)
Returns
number of peaks

Definition at line 450 of file TOPreco.cc.

◆ getNumOfPhotons()

int getNumOfPhotons ( )

Return number of measured photons.

Returns
number of photons

Definition at line 297 of file TOPreco.cc.

◆ getPDF()

double getPDF ( int  pixelID,
double  t,
double  mass,
int  PDG,
double  jitter = 0 
)

Return PDF for pixel pixelID at time t for mass hypothesis mass.

Parameters
pixelIDpixel ID (e.g. software channel, 1-based)
ttime
massparticle mass
PDGparticle code
jitteradditional time jitter, like electronic jitter

Definition at line 428 of file TOPreco.cc.

◆ getPDFPeak()

void getPDFPeak ( int  pixelID,
int  k,
float &  position,
float &  width,
float &  numPhotons 
) const

Returns k-th PDF peak for given pixel describing signal PDF.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)

Definition at line 456 of file TOPreco.cc.

◆ getPDFPeakE()

float getPDFPeakE ( int  pixelID,
int  k 
) const

Returns photon energy of PDF peak.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
photon energy [eV]

Definition at line 484 of file TOPreco.cc.

◆ getPDFPeakFic()

float getPDFPeakFic ( int  pixelID,
int  k 
) const

Returns Cerenkov azimuthal angle of PDF peak.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
azimuthal angle

Definition at line 477 of file TOPreco.cc.

◆ getPDFPeakKxd()

float getPDFPeakKxd ( int  pixelID,
int  k 
) const

Returns photon reconstructed direction in x at detection.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
direction in x

Definition at line 575 of file TOPreco.cc.

◆ getPDFPeakKxe()

float getPDFPeakKxe ( int  pixelID,
int  k 
) const

Returns photon reconstructed direction in x at emission.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
direction in x

Definition at line 554 of file TOPreco.cc.

◆ getPDFPeakKyd()

float getPDFPeakKyd ( int  pixelID,
int  k 
) const

Returns photon reconstructed direction in y at detection.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
direction in y

Definition at line 582 of file TOPreco.cc.

◆ getPDFPeakKye()

float getPDFPeakKye ( int  pixelID,
int  k 
) const

Returns photon reconstructed direction in y at emission.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
direction in y

Definition at line 561 of file TOPreco.cc.

◆ getPDFPeakKzd()

float getPDFPeakKzd ( int  pixelID,
int  k 
) const

Returns photon reconstructed direction in z at detection.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
direction in z

Definition at line 589 of file TOPreco.cc.

◆ getPDFPeakKze()

float getPDFPeakKze ( int  pixelID,
int  k 
) const

Returns photon reconstructed direction in z at emission.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
direction in z

Definition at line 568 of file TOPreco.cc.

◆ getPDFPeakNx()

int getPDFPeakNx ( int  pixelID,
int  k 
) const

Returns total number of reflections in x of PDF peak.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
total number of reflections

Definition at line 498 of file TOPreco.cc.

◆ getPDFPeakNxe()

int getPDFPeakNxe ( int  pixelID,
int  k 
) const

Returns number of reflections in x in prism.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
total number of reflections before mirror

Definition at line 526 of file TOPreco.cc.

◆ getPDFPeakNxm()

int getPDFPeakNxm ( int  pixelID,
int  k 
) const

Returns number of reflections in x before mirror.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
total number of reflections before mirror

Definition at line 512 of file TOPreco.cc.

◆ getPDFPeakNy()

int getPDFPeakNy ( int  pixelID,
int  k 
) const

Returns total number of reflections in y of PDF peak.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
total number of reflections

Definition at line 505 of file TOPreco.cc.

◆ getPDFPeakNye()

int getPDFPeakNye ( int  pixelID,
int  k 
) const

Returns number of reflections in y in prism.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
total number of reflections before mirror

Definition at line 533 of file TOPreco.cc.

◆ getPDFPeakNym()

int getPDFPeakNym ( int  pixelID,
int  k 
) const

Returns number of reflections in y before mirror.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
total number of reflections before mirror

Definition at line 519 of file TOPreco.cc.

◆ getPDFPeakSigE()

float getPDFPeakSigE ( int  pixelID,
int  k 
) const

Returns photon energy spread of PDF peak.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
photon energy sigma [eV]

Definition at line 491 of file TOPreco.cc.

◆ getPDFPeakType()

int getPDFPeakType ( int  pixelID,
int  k 
) const

Returns type of the k-th PDF peak for given pixel.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
0 unknown, 1 direct, 2 reflected

Definition at line 470 of file TOPreco.cc.

◆ getPDFPeakXD()

float getPDFPeakXD ( int  pixelID,
int  k 
) const

Returns unfolded x position of pixel.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
unfolded x

Definition at line 540 of file TOPreco.cc.

◆ getPDFPeakYD()

float getPDFPeakYD ( int  pixelID,
int  k 
) const

Returns unfolded y position of pixel.

Parameters
pixelIDpixel ID (1-based)
kpeak counter (in C++ sense - starts with 0)
Returns
unfolded y

Definition at line 547 of file TOPreco.cc.

◆ getPull()

void getPull ( int  K,
int &  pixelID,
float &  T,
float &  T0,
float &  Wid,
float &  PhiCer,
float &  Wt 
)

Get pulls: K-th pull.

Parameters
Kcounter
pixelIDpixel ID (e.g. software channel, 1-based)
Tphoton time
T0PDF mean time
WidPDF sigma
PhiCerazimuthal Cerenkov angle
Wtweight

Definition at line 420 of file TOPreco.cc.

◆ getPullSize()

int getPullSize ( )

Get pulls: size.

Returns
: size

Definition at line 413 of file TOPreco.cc.

◆ getTimeWindow()

void getTimeWindow ( double &  Tmin,
double &  Tmax 
)

Returns time window for photons.

Parameters
Tminminimum time [ns]
Tmaxmaximum time [ns]

Definition at line 187 of file TOPreco.cc.

◆ getTrackHit()

void getTrackHit ( int  LocGlob,
double  R[3],
double  Dir[3],
double &  Len,
double &  Tlen,
double &  Mom,
int &  moduleID 
)

Return track hit at the bar in Local or Global frame.

Parameters
LocGlobselect Local or Global frame
Rtrack spatial position
Dirtrack direction (unit vector)
Lentrack length inside bar
Tlentrack length from IP
Momtrack momentum
moduleIDmodule ID if hit, else -1

Definition at line 342 of file TOPreco.cc.

◆ reconstruct() [1/2]

void reconstruct ( const TOPtrack trk,
int  pdg = 0 
)

Run reconstruction for a given track.

Parameters
trktrack
pdgPDG code for which to compute pulls. If 0 use MC true PDG.

Definition at line 268 of file TOPreco.cc.

◆ reconstruct() [2/2]

void reconstruct ( double  X,
double  Y,
double  Z,
double  Tlen,
double  Px,
double  Py,
double  Pz,
int  Q,
int  pdg = 0,
int  moduleID = 0 
)
private

Run reconstruction.

Parameters
Xtrack spatial position x
Ytrack spatial position y
Ztrack spatial position z
Tlentrack length from IP
Pxtrack momentum component x
Pytrack momentum component y
Pztrack momentum component z
Qtrack charge
pdgPDG code for which to compute pulls
moduleIDmodule ID (optional)

Definition at line 250 of file TOPreco.cc.

◆ redoPDF()

void redoPDF ( double  mass,
int  PDG 
)

Re-calculate PDF for a given particle mass using option c_Fine.

Parameters
massparticle mass
PDGparticle code

Definition at line 444 of file TOPreco.cc.

◆ setBeta()

void setBeta ( double  beta)

Set track beta (for beta resolution studies) if beta>0 this value is used instead of beta from momentum and mass.

Parameters
betabeta value

Definition at line 437 of file TOPreco.cc.

◆ setChannelMask()

void setChannelMask ( const DBObjPtr< TOPCalChannelMask > &  mask,
const TOPAsicMask asicMask 
)
static

Set channel mask.

Parameters
maskchannel mask
asicMaskmasked asics

Definition at line 88 of file TOPreco.cc.

◆ setMass()

void setMass ( double  mass,
int  pdg 
)

Set mass of the particle hypothesis (overrides settings in the constructor)

Parameters
massmass
pdgPDG code

Definition at line 172 of file TOPreco.cc.

◆ setPDFoption()

void setPDFoption ( PDFoption  opt,
int  NP = 0,
int  NC = 0 
)

Sets PDF option.

Parameters
optoption - see definition of PDFoption
NPnumber of emission positions along track segment (equidistant)
NCnumber of Cerenkov angles (equdistant in photon energies)

Definition at line 196 of file TOPreco.cc.

◆ setStoreOption()

void setStoreOption ( StoreOption  opt)

Sets option for storing PDF parameters in Fortran common TOP_PIK.

Parameters
optoption - see definition of StoreOption

Definition at line 202 of file TOPreco.cc.

◆ setTimeWindow()

void setTimeWindow ( double  Tmin,
double  Tmax 
)

Set time window for photons.

Allows to set window different than that defined by parameters of TOPNominalTDC.

If Tmax <= Tmin the window is set to default from TOPNominalTDC.

Window edges must not exceed those used during data taking or in simulation

Parameters
Tminminimum time [ns]
Tmaxmaximum time [ns]

Definition at line 180 of file TOPreco.cc.

◆ setUncalibratedChannelsOff() [1/2]

void setUncalibratedChannelsOff ( const DBObjPtr< TOPCalChannelT0 > &  channelT0)
static

Set uncalibrated channels off.

Parameters
channelT0channel T0 calibration

Definition at line 107 of file TOPreco.cc.

◆ setUncalibratedChannelsOff() [2/2]

void setUncalibratedChannelsOff ( const DBObjPtr< TOPCalTimebase > &  timebase)
static

Set uncalibrated channels off.

Parameters
timebasetimebase calibration

Definition at line 125 of file TOPreco.cc.


The documentation for this class was generated from the following files:
Belle2::TOP::TOPreco::m_beta
double m_beta
beta value, if set
Definition: TOPreco.h:485
Belle2::TOP::TOPreco::setBeta
void setBeta(double beta)
Set track beta (for beta resolution studies) if beta>0 this value is used instead of beta from moment...
Definition: TOPreco.cc:437
Belle2::TOP::TOPreco::setPDFoption
void setPDFoption(PDFoption opt, int NP=0, int NC=0)
Sets PDF option.
Definition: TOPreco.cc:196
Belle2::TOP::TOPreco::setTimeWindow
void setTimeWindow(double Tmin, double Tmax)
Set time window for photons.
Definition: TOPreco.cc:180