CMS 3D CMS Logo

List of all members | Public Member Functions | Private Attributes | Static Private Attributes
SiPixelCPEGenericDBErrorParametrization Class Reference

#include <SiPixelCPEGenericDBErrorParametrization.h>

Public Member Functions

std::pair< float, float > getError (const SiPixelCPEGenericErrorParm *parmErrors, GeomDetType::SubDetector pixelPart, int sizex, int sizey, float alpha, float beta, bool bigInX=false, bool bigInY=false)
 
std::pair< float, float > getError (GeomDetType::SubDetector pixelPart, int sizex, int sizey, float alpha, float beta, bool bigInX=false, bool bigInY=false)
 
float index (int ind_subpart, int size, float alpha, float beta, bool big)
 
void setDBAccess (const edm::EventSetup &es)
 
 SiPixelCPEGenericDBErrorParametrization ()
 
 ~SiPixelCPEGenericDBErrorParametrization ()
 

Private Attributes

edm::ESHandle< SiPixelCPEGenericErrorParmerrorsH
 

Static Private Attributes

static const float * a_max [4] = {by_a_max, bx_a_max, fy_a_max, fx_a_max}
 
static const float * a_min [4] = {by_a_min, bx_a_min, fy_a_min, fx_a_min}
 
static const int alpha_bin_size [4] = {10, 1, 10, 1}
 
static const float * b_max [4] = {by_b_max, garbage, fy_b_max, fx_b_max}
 
static const float * b_min [4] = {by_b_min, garbage, fy_b_min, fx_b_min}
 
static const int beta_bin_size [4] = {1, 10, 1, 10}
 
static const float bx_a_max [3] = {1.725, 1.675, 1.625}
 
static const float bx_a_min [3] = {1.525, 1.475, 1.425}
 
static const float by_a_max [6]
 
static const float by_a_min [6]
 
static const float by_b_max [6] = {0.50, 0.90, 1.05, 1.15, 1.20, 1.40}
 
static const float by_b_min [6] = {0.05, 0.15, 0.70, 0.95, 1.15, 1.20}
 
static const float errors_big_pix [4] = {0.0070, 0.0030, 0.0068, 0.0040}
 
static const float fx_a_max [2] = {0.285, 0.465}
 
static const float fx_a_min [2] = {0.165, 0.185}
 
static const float fx_b_max [2] = {999, 999}
 
static const float fx_b_min [2] = {998, 998}
 
static const float fy_a_max [2] = {999, 999}
 
static const float fy_a_min [2] = {998, 998}
 
static const float fy_b_max [2] = {0.39, 0.39}
 
static const float fy_b_min [2] = {0.31, 0.31}
 
static const int part_bin_size [4] = {0, 240, 360, 380}
 
static const int size_bin_size [4] = {40, 40, 40, 40}
 
static const int size_max [4] = {5, 2, 0, 0}
 

Detailed Description

Definition at line 11 of file SiPixelCPEGenericDBErrorParametrization.h.

Constructor & Destructor Documentation

SiPixelCPEGenericDBErrorParametrization::SiPixelCPEGenericDBErrorParametrization ( )

Definition at line 48 of file SiPixelCPEGenericDBErrorParametrization.cc.

48 {}
SiPixelCPEGenericDBErrorParametrization::~SiPixelCPEGenericDBErrorParametrization ( )

Definition at line 50 of file SiPixelCPEGenericDBErrorParametrization.cc.

50 {}

Member Function Documentation

std::pair< float, float > SiPixelCPEGenericDBErrorParametrization::getError ( const SiPixelCPEGenericErrorParm parmErrors,
GeomDetType::SubDetector  pixelPart,
int  sizex,
int  sizey,
float  alpha,
float  beta,
bool  bigInX = false,
bool  bigInY = false 
)

Definition at line 57 of file SiPixelCPEGenericDBErrorParametrization.cc.

References MessageLogger_cfi::errors, SiPixelCPEGenericErrorParm::errors(), Exception, index(), createfilelist::int, GeomDetEnumerators::isBarrel(), and GeomDetEnumerators::isTrackerPixel().

64  {
65  std::pair<float, float> element;
66  std::pair<float, float> errors;
67 
68  if (!GeomDetEnumerators::isTrackerPixel(pixelPart))
69  throw cms::Exception("PixelCPEGenericDBErrorParametrization::getError") << "Non-pixel detector type !!!";
70  if (GeomDetEnumerators::isBarrel(pixelPart)) {
71  element = std::pair<float, float>(index(1, sizex, alpha, beta, bigInX), //1 -- Bx
72  index(0, sizey, alpha, beta, bigInY)); //0 -- By
73  } else {
74  element = std::pair<float, float>(index(3, sizex, alpha, beta, bigInX), //3 -- Fx
75  index(2, sizey, alpha, beta, bigInY)); //2 -- Fy
76  }
77 
78  if (bigInX && sizex == 1)
79  errors.first = element.first;
80  else
81  errors.first = parmErrors->errors()[(int)element.first].sigma;
82  if (bigInY && sizey == 1)
83  errors.second = element.second;
84  else
85  errors.second = parmErrors->errors()[(int)element.second].sigma;
86 
87  return errors;
88 }
bool isBarrel(GeomDetEnumerators::SubDetector m)
float index(int ind_subpart, int size, float alpha, float beta, bool big)
DbVector & errors()
Accessors for the vectors – non-const version.
alpha
zGenParticlesMatch = cms.InputTag(""),
bool isTrackerPixel(GeomDetEnumerators::SubDetector m)
std::pair< float, float > SiPixelCPEGenericDBErrorParametrization::getError ( GeomDetType::SubDetector  pixelPart,
int  sizex,
int  sizey,
float  alpha,
float  beta,
bool  bigInX = false,
bool  bigInY = false 
)

Definition at line 91 of file SiPixelCPEGenericDBErrorParametrization.cc.

References MessageLogger_cfi::errors, SiPixelCPEGenericErrorParm::errors(), errorsH, Exception, index(), createfilelist::int, GeomDetEnumerators::isBarrel(), and GeomDetEnumerators::isTrackerPixel().

92  {
93  std::pair<float, float> element;
94  std::pair<float, float> errors;
95 
96  if (!GeomDetEnumerators::isTrackerPixel(pixelPart))
97  throw cms::Exception("PixelCPEGenericDBErrorParametrization::getError") << "Non-pixel detector type !!!";
98  if (GeomDetEnumerators::isBarrel(pixelPart)) {
99  element = std::pair<float, float>(index(1, sizex, alpha, beta, bigInX), //1 -- Bx
100  index(0, sizey, alpha, beta, bigInY)); //0 -- By
101  } else {
102  element = std::pair<float, float>(index(3, sizex, alpha, beta, bigInX), //3 -- Fx
103  index(2, sizey, alpha, beta, bigInY)); //2 -- Fy
104  }
105 
106  if (bigInX && sizex == 1)
107  errors.first = element.first;
108  else
109  errors.first = errorsH->errors()[(int)element.first].sigma;
110  if (bigInY && sizey == 1)
111  errors.second = element.second;
112  else
113  errors.second = errorsH->errors()[(int)element.second].sigma;
114 
115  return errors;
116 }
bool isBarrel(GeomDetEnumerators::SubDetector m)
float index(int ind_subpart, int size, float alpha, float beta, bool big)
DbVector & errors()
Accessors for the vectors – non-const version.
alpha
zGenParticlesMatch = cms.InputTag(""),
edm::ESHandle< SiPixelCPEGenericErrorParm > errorsH
bool isTrackerPixel(GeomDetEnumerators::SubDetector m)
float SiPixelCPEGenericDBErrorParametrization::index ( int  ind_subpart,
int  size,
float  alpha,
float  beta,
bool  big 
)

Definition at line 118 of file SiPixelCPEGenericDBErrorParametrization.cc.

References a_max, a_min, alpha_bin_size, b_max, b_min, beta_bin_size, errors_big_pix, createfilelist::int, math_pi, min(), part_bin_size, size_bin_size, and size_max.

Referenced by getError(), and BeautifulSoup.PageElement::insert().

118  {
119  //This is a check for big pixels. If it passes, the code returns a given error and the function ends.
120  if (big && size == 1)
121  return errors_big_pix[ind_subpart];
122 
123  int ind_size = std::min(size - 1, size_max[ind_subpart]);
124 
125  float alpha_rad = -999.9;
126  float betap_rad = -999.9;
127 
128  int ind_alpha = -99;
129  int ind_beta = -99;
130 
131  float binw_a = -999.9;
132  float binw_b = -999.9;
133  int maxbin_a = -99;
134  int maxbin_b = -99;
135 
136  betap_rad = fabs(math_pi / 2.0 - beta);
137  //We must take into account that Fx(subpart=3) has different alpha parametrization
138  if (ind_subpart == 3)
139  alpha_rad = fabs(math_pi / 2.0 - alpha);
140  else
141  alpha_rad = fabs(alpha);
142 
143  //Sets the correct binning for alpha and beta based on whether in x or y
144  if (ind_subpart == 0 || ind_subpart == 2) {
145  binw_a = (a_max[ind_subpart][ind_size] - a_min[ind_subpart][ind_size]) / 2.0;
146  binw_b = (b_max[ind_subpart][ind_size] - b_min[ind_subpart][ind_size]) / 8.0;
147  maxbin_a = 3;
148  maxbin_b = 9;
149  } else {
150  binw_a = (a_max[ind_subpart][ind_size] - a_min[ind_subpart][ind_size]) / 8.0;
151  binw_b = (b_max[ind_subpart][ind_size] - b_min[ind_subpart][ind_size]) / 2.0;
152  maxbin_a = 3;
153  maxbin_b = 9;
154  }
155 
156  //Binning for alpha
157  if (alpha_rad < a_min[ind_subpart][ind_size])
158  ind_alpha = 0;
159  else if (alpha_rad >= a_max[ind_subpart][ind_size])
160  ind_alpha = maxbin_a;
161  else
162  ind_alpha = 1 + (int)((alpha_rad - a_min[ind_subpart][ind_size]) / binw_a);
163 
164  //Binning for beta -- we need to account for Bx(subpart=1) having uneven binning
165  if (ind_subpart == 1) {
166  if (betap_rad <= 0.7)
167  ind_beta = 0;
168  else if (0.7 < betap_rad && betap_rad <= 1.0)
169  ind_beta = 1;
170  else if (1.0 < betap_rad && betap_rad <= 1.2)
171  ind_beta = 2;
172  else if (1.2 <= betap_rad)
173  ind_beta = 3;
174  } else if (betap_rad < b_min[ind_subpart][ind_size])
175  ind_beta = 0;
176  else if (betap_rad >= b_max[ind_subpart][ind_size])
177  ind_beta = maxbin_b;
178  else
179  ind_beta = 1 + (int)((betap_rad - b_min[ind_subpart][ind_size]) / binw_b);
180 
181  //Index to be used to find error in database
182  int index = part_bin_size[ind_subpart] + size_bin_size[ind_subpart] * ind_size +
183  alpha_bin_size[ind_subpart] * ind_alpha + beta_bin_size[ind_subpart] * ind_beta;
184 
185  return index;
186 }
size
Write out results.
T min(T a, T b)
Definition: MathUtil.h:58
float index(int ind_subpart, int size, float alpha, float beta, bool big)
Definition: big.h:8
alpha
zGenParticlesMatch = cms.InputTag(""),
void SiPixelCPEGenericDBErrorParametrization::setDBAccess ( const edm::EventSetup es)

Member Data Documentation

const float * SiPixelCPEGenericDBErrorParametrization::a_max = {by_a_max, bx_a_max, fy_a_max, fx_a_max}
staticprivate

Definition at line 67 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const float * SiPixelCPEGenericDBErrorParametrization::a_min = {by_a_min, bx_a_min, fy_a_min, fx_a_min}
staticprivate

Definition at line 66 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const int SiPixelCPEGenericDBErrorParametrization::alpha_bin_size = {10, 1, 10, 1}
staticprivate

Definition at line 63 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const float * SiPixelCPEGenericDBErrorParametrization::b_max = {by_b_max, garbage, fy_b_max, fx_b_max}
staticprivate

Definition at line 69 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const float * SiPixelCPEGenericDBErrorParametrization::b_min = {by_b_min, garbage, fy_b_min, fx_b_min}
staticprivate

Definition at line 68 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const int SiPixelCPEGenericDBErrorParametrization::beta_bin_size = {1, 10, 1, 10}
staticprivate

Definition at line 64 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const float SiPixelCPEGenericDBErrorParametrization::bx_a_max = {1.725, 1.675, 1.625}
staticprivate

Definition at line 41 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::bx_a_min = {1.525, 1.475, 1.425}
staticprivate

Definition at line 40 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::by_a_max
staticprivate
Initial value:
= {
1.67078, 1.67078, 1.67078, 1.67078, 1.67078, 1.67078}

Definition at line 49 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::by_a_min
staticprivate
Initial value:
= {
1.47078, 1.47078, 1.47078, 1.47078, 1.47078, 1.47078}

Definition at line 48 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::by_b_max = {0.50, 0.90, 1.05, 1.15, 1.20, 1.40}
staticprivate

Definition at line 51 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::by_b_min = {0.05, 0.15, 0.70, 0.95, 1.15, 1.20}
staticprivate

Definition at line 50 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::errors_big_pix = {0.0070, 0.0030, 0.0068, 0.0040}
staticprivate

Definition at line 58 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

edm::ESHandle<SiPixelCPEGenericErrorParm> SiPixelCPEGenericDBErrorParametrization::errorsH
private

Definition at line 38 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by getError(), and setDBAccess().

const float SiPixelCPEGenericDBErrorParametrization::fx_a_max = {0.285, 0.465}
staticprivate

Definition at line 44 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fx_a_min = {0.165, 0.185}
staticprivate

Definition at line 43 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fx_b_max = {999, 999}
staticprivate

Definition at line 46 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fx_b_min = {998, 998}
staticprivate

Definition at line 45 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fy_a_max = {999, 999}
staticprivate

Definition at line 54 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fy_a_min = {998, 998}
staticprivate

Definition at line 53 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fy_b_max = {0.39, 0.39}
staticprivate

Definition at line 56 of file SiPixelCPEGenericDBErrorParametrization.h.

const float SiPixelCPEGenericDBErrorParametrization::fy_b_min = {0.31, 0.31}
staticprivate

Definition at line 55 of file SiPixelCPEGenericDBErrorParametrization.h.

const int SiPixelCPEGenericDBErrorParametrization::part_bin_size = {0, 240, 360, 380}
staticprivate

Definition at line 61 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const int SiPixelCPEGenericDBErrorParametrization::size_bin_size = {40, 40, 40, 40}
staticprivate

Definition at line 62 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().

const int SiPixelCPEGenericDBErrorParametrization::size_max = {5, 2, 0, 0}
staticprivate

Definition at line 59 of file SiPixelCPEGenericDBErrorParametrization.h.

Referenced by index().