CMS 3D CMS Logo

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

#include <TSFit.h>

Inheritance diagram for TSFit:

Public Member Functions

double fit_third_degree_polynomial (double *, double *)
 
double fpol3dg (int, double *, double *, double *)
 
void init_errmat (double)
 
double inverms (int, double xx[matdim][matdim], double yy[matdim][matdim])
 
void set_params (int, int, int, int, int, double, double, int, int)
 
 TSFit (int size=SDIM, int size_sh=PLSHDIM)
 
 ~TSFit () override
 

Public Attributes

int plshdim
 
int sdim
 

Static Public Attributes

static constexpr unsigned int diminpar = 10
 
static constexpr unsigned int dimoutpar = 10
 
static constexpr int matdim = 5
 
static constexpr unsigned int npar_moni = 4
 
static constexpr unsigned int PLSHDIM = 650
 
static constexpr unsigned int SDIM = 14
 

Private Attributes

double acc [SDIM]
 
double adcp [SDIM]
 
double adfmx [SDIM]
 
double al [matdim][matdim]
 
double alpha_th
 
double avtm
 
double be [matdim][matdim]
 
double beta_th
 
double corel [SDIM]
 
double cov [matdim][matdim]
 
double der [SDIM][5]
 
double errmat [SDIM][SDIM]
 
double f [SDIM]
 
double ff [SDIM][4]
 
int iinf
 
double invcov [matdim][matdim]
 
int isup
 
double maskp3 [SDIM]
 
int n_presamples
 
int n_samples_aft_max
 
int n_samples_bef_max
 
double nbcor [SDIM]
 
int nbr_iter_fit
 
int nbs
 
int nmxu_sto
 
double norme
 
double parfp3 [dimoutpar]
 
int sample_flag [SDIM]
 
double t [SDIM]
 
double xki2_max
 
double z [SDIM]
 

Detailed Description

Definition at line 6 of file TSFit.h.

Constructor & Destructor Documentation

◆ TSFit()

TSFit::TSFit ( int  size = SDIM,
int  size_sh = PLSHDIM 
)

Definition at line 19 of file TSFit.cc.

References plshdim, sdim, and findQualityFiles::size.

19  {
20  sdim = size;
21  plshdim = size_sh;
22  // sample_flag = new int[size];
23  // t = new double[sdim];
24  // z = new double[sdim];
25  // f = new double[sdim];
26  // acc = new double[sdim];
27  // adfmx = new double[sdim];
28  // adcp = new double[sdim];
29  // maskp3 = new double[sdim];
30  // corel = new double[sdim];
31  // nbcor = new double[sdim];
32  // int k;
33  // ff = new double *[sdim];
34  // for(k=0;k<sdim;k++)ff[k] = new double[4];
35  // der = new double *[sdim];
36  // for(k=0;k<sdim;k++)der[k] = new double[5];
37 }
size
Write out results.
int plshdim
Definition: TSFit.h:64
int sdim
Definition: TSFit.h:63

◆ ~TSFit()

TSFit::~TSFit ( )
inlineoverride

Definition at line 68 of file TSFit.h.

68 {}

Member Function Documentation

◆ fit_third_degree_polynomial()

double TSFit::fit_third_degree_polynomial ( double *  bdc,
double *  ret_dat 
)

Definition at line 255 of file TSFit.cc.

References adfmx, errmat, dqmMemoryStats::float, fpol3dg(), mps_fire::i, iinf, isup, dqmdumpme::k, maskp3, nbs, norme, parfp3, sample_flag, createJobs::tmp, and xki2_max.

Referenced by EcalTestPulseAnalyzer::analyze().

255  {
256  // third degree polynomial fit of the pulse summit.
257  // samples are contained in array bdc and must be pedestal
258  // substracted.
259  // only samples having sample_flag >= 1 are used.
260  // the unit of time is one clock unit, that is to say 25 ns.
261  // output: ret_dat[0] = pulse height
262  // ret_dat[1] position of maximum in the sample frame in clock units
263  // ret_dat[2] adc value of the highest sample
264  // ret_dat[3] number of the highest sample
265  // ret_dat[4] lower sample number used for fitting
266  // ret_dat[5] upper sample number used for fitting
267  // errmat inverse of the error matrix
268 
269  int i;
270  int nus;
271  double xki2;
272  double tm, tmp, amp;
273 
274  ret_dat[0] = -999.;
275  ret_dat[1] = -999.;
276 
277  // search the maximum
278  double val_max = 0.;
279  int imax = 0;
280  for (i = 0; i < nbs; i++) {
281  if (sample_flag[i] == 0)
282  continue;
283  if (bdc[i] > val_max) {
284  val_max = bdc[i];
285  imax = i;
286  }
287  }
288 
289  if ((val_max * val_max) * errmat[imax][imax] < 16.)
290  return -118;
291 
292  // if( imax != 9 )printf( "imax : %d !!!!!!!!!!!!!!!!!!!!!!!!!!!\n", imax );
293 
294  if (norme == 0.)
295  norme = val_max;
296 
297  // look for samples above 1/3 of maximum before and 1/2 after
298  double val2 = val_max / 2.;
299  double val3 = val_max / 2.;
300  int ilow = iinf;
301  int ihig = 0;
302 
303  for (i = iinf; i <= isup; i++) {
304  if (sample_flag[i] >= 1) {
305  if ((bdc[i] < val3) && (i < imax))
306  ilow = i;
307  if (bdc[i] > val2)
308  ihig = i;
309  }
310  }
311 
312  ilow++;
313 
314  //ilow = imax - 1;
315 
316  /* le test suivant, apparemment idiot, est mis a cause des sequences 0. 2048. qui apparaissent dans certains mauvais evts JPP 11/09/00 */
317 
318  if (ihig == ilow)
319  return -105;
320  if (ilow == imax)
321  ilow = ilow - 1;
322  //if( ihig - ilow < 3 )ihig = ilow + 3;
323  ihig = ilow + 3;
324 
325  /* printf(" third degree: ilow %d ihig %d \n",ilow,ihig); */
326  nus = 0;
327  int number_of_good_samples = 0;
328  for (i = ilow; i <= ihig; i++) {
329  maskp3[nus] = 0;
330  adfmx[nus] = 0.;
331  /* printf(" adc %f sample_flag %d number_of good_samples %d \n",bdc[i],sample_flag[i],number_of_good_samples); */
332  if (sample_flag[i] >= 1) {
333  adfmx[nus] = bdc[i];
334  maskp3[nus] = 1.;
335  number_of_good_samples++;
336  }
337  nus++;
338  }
339 
340  if (number_of_good_samples < 4) {
341  return (-106);
342  }
343 
344  xki2 = fpol3dg(nus, &parfp3[0], &maskp3[0], &adfmx[0]);
345 
346  /* printf( "fpol3dg-----------------------------------> %f %f %f %f %f\n",
347  parfp3[0], parfp3[1], parfp3[2], parfp3[3], parfp3[4] ); */
348 
349  tm = parfp3[4] + (float)ilow;
350  amp = parfp3[5];
351 
352  if (amp * amp * errmat[0][0] < 2.)
353  return -101.;
354  tmp = parfp3[6] + (float)ilow;
355 
356  /*
357  validation of fit quality. Most of the time the fit is done with
358  four samples, therefore there is no possible ki2 check. When more than
359  4 samples are used the ki2 is often bad. So, in order to suppress some
360  events with bad samples, a consistency check on the position of the
361  maximum and minimum of the 3rd degree polynomial is used.
362  */
363 
364  if (xki2 > xki2_max) {
365  return -102.;
366  }
367  if ((tm < (double)ilow) || (tm > (double)ihig)) {
368  return -103.;
369  }
370 
371  if ((tmp > (double)ilow) && (tmp < (double)ihig - 1.)) {
372  return -104.;
373  }
374 
375  ret_dat[0] = amp;
376  ret_dat[1] = tm;
377  ret_dat[2] = val_max;
378  ret_dat[3] = (double)imax;
379  ret_dat[4] = (double)ilow;
380  ret_dat[5] = (double)ihig;
381  ret_dat[6] = (double)tmp;
382  int k;
383  for (i = 0; i < 4; i++) {
384  k = i + 7;
385  ret_dat[k] = parfp3[i];
386  }
387 
388  return xki2;
389 }
double parfp3[dimoutpar]
Definition: TSFit.h:40
int nbs
Definition: TSFit.h:21
double xki2_max
Definition: TSFit.h:30
double fpol3dg(int, double *, double *, double *)
Definition: TSFit.cc:105
double norme
Definition: TSFit.h:30
double errmat[SDIM][SDIM]
Definition: TSFit.h:47
double maskp3[SDIM]
Definition: TSFit.h:55
double adfmx[SDIM]
Definition: TSFit.h:53
int sample_flag[SDIM]
Definition: TSFit.h:48
int isup
Definition: TSFit.h:23
int iinf
Definition: TSFit.h:23
tmp
align.sh
Definition: createJobs.py:716

◆ fpol3dg()

double TSFit::fpol3dg ( int  nmxul,
double *  parom,
double *  mask,
double *  adc 
)

Definition at line 105 of file TSFit.cc.

References gpuClustering::adc, cov, dumpMFGeometry_cfg::delta, errmat, ff, h, mps_fire::i, invcov, inverms(), dqmdumpme::k, cmsLHEtoEOSManager::l, gpuClustering::pixelStatus::mask, alignCSCRings::s, mathSSE::sqrt(), t, RandomServiceHelper::t2, and createJobs::tmp.

Referenced by fit_third_degree_polynomial().

105  {
106  // fit third degree polynomial
107  // nmxul array adc[] length
108  // parom return parameters (a0,a1,a2,a3,pos max,height)
109  // fplo3dg uses only the diagonal terms of errmat[][]
110  // errmat inverse of the error matrix
111 
112  int i, k, l;
113  double h, t2, tm, delta, tmp;
114  double xki2, dif, difmx, deglib;
115  double bv[4], s;
116 
117  deglib = (double)nmxul - 4.;
118  for (i = 0; i < nmxul; i++) {
119  t[i] = i;
120  ff[i][0] = 1.;
121  ff[i][1] = t[i];
122  ff[i][2] = t[i] * t[i];
123  ff[i][3] = ff[i][2] * t[i];
124  }
125  /* computation of covariance matrix */
126  for (k = 0; k < 4; k++) {
127  for (l = 0; l < 4; l++) {
128  s = 0.;
129  for (i = 0; i < nmxul; i++) {
130  s = s + ff[i][k] * ff[i][l] * errmat[i][i] * mask[i];
131  }
132  cov[k][l] = s;
133  }
134  s = 0.;
135  for (i = 0; i < nmxul; i++) {
136  s = s + ff[i][k] * adc[i] * errmat[i][i] * mask[i];
137  }
138  bv[k] = s;
139  }
140  /* parameters */
141  inverms(4, cov, invcov);
142  for (k = 0; k < 4; k++) {
143  s = 0.;
144  for (l = 0; l < 4; l++) {
145  s = s + bv[l] * invcov[l][k];
146  }
147  parom[k] = s;
148  }
149 
150  if (parom[3] == 0.) {
151  parom[4] = -1000.;
152  parom[5] = -1000.;
153  parom[6] = -1000.;
154  return 1000000.;
155  }
156  /* worst hit and ki2 */
157  xki2 = 0.;
158  difmx = 0.;
159  for (i = 0; i < nmxul; i++) {
160  t2 = t[i] * t[i];
161  h = parom[0] + parom[1] * t[i] + parom[2] * t2 + parom[3] * t2 * t[i];
162  dif = (adc[i] - h) * mask[i];
163  xki2 = xki2 + dif * dif * errmat[i][i];
164  if (dif > difmx) {
165  difmx = dif;
166  }
167  }
168  if (deglib > 0.5)
169  xki2 = xki2 / deglib;
170  /* amplitude and maximum position */
171  delta = parom[2] * parom[2] - 3. * parom[3] * parom[1];
172  if (delta > 0.) {
173  delta = sqrt(delta);
174  tm = -(delta + parom[2]) / (3. * parom[3]);
175  tmp = (delta - parom[2]) / (3. * parom[3]);
176  } else {
177  parom[4] = -1000.;
178  parom[5] = -1000.;
179  parom[6] = -1000.;
180  return xki2;
181  }
182  parom[4] = tm;
183  parom[5] = parom[0] + parom[1] * tm + parom[2] * tm * tm + parom[3] * tm * tm * tm;
184  parom[6] = tmp;
185  return xki2;
186 }
double invcov[matdim][matdim]
Definition: TSFit.h:35
double t[SDIM]
Definition: TSFit.h:49
constexpr uint32_t mask
Definition: gpuClustering.h:26
double cov[matdim][matdim]
Definition: TSFit.h:35
double inverms(int, double xx[matdim][matdim], double yy[matdim][matdim])
Definition: TSFit.cc:187
T sqrt(T t)
Definition: SSEVec.h:19
double errmat[SDIM][SDIM]
Definition: TSFit.h:47
double ff[SDIM][4]
Definition: TSFit.h:59
The Signals That Services Can Subscribe To This is based on ActivityRegistry h
Helper function to determine trigger accepts.
Definition: Activities.doc:4
tmp
align.sh
Definition: createJobs.py:716
uint16_t *__restrict__ uint16_t const *__restrict__ adc

◆ init_errmat()

void TSFit::init_errmat ( double  noise_initialvalue)

Definition at line 87 of file TSFit.cc.

References errmat, mps_fire::i, dqmiolumiharvest::j, and sdim.

Referenced by EcalTestPulseAnalyzer::analyze().

87  {
88  // Input: noise_initial value noise (in adc channels) as read in the
89  // data base.
90  /*------------------------------------------------------------------------*/
91 
92  int i, j;
93  //double one_over_noisesq;
94 
95  //one_over_noisesq = 1. / ( noise_initialvalue * noise_initialvalue );
96  for (i = 0; i < sdim; i++) {
97  for (j = 0; j < sdim; j++) {
98  errmat[i][j] = noise_initialvalue;
99  //errmat[i][j] = 0.;
100  }
101  //errmat[i][i] = one_over_noisesq;
102  }
103 }
int sdim
Definition: TSFit.h:63
double errmat[SDIM][SDIM]
Definition: TSFit.h:47

◆ inverms()

double TSFit::inverms ( int  n,
double  xx[matdim][matdim],
double  yy[matdim][matdim] 
)

Definition at line 187 of file TSFit.cc.

References al, be, g, mps_fire::i, dqmiolumiharvest::j, findQualityFiles::jj, dqmdumpme::k, matdim, dqmiodumpmetadata::n, alignCSCRings::r, alignCSCRings::s, mathSSE::sqrt(), and SiPixelPI::zero.

Referenced by fpol3dg().

187  {
188  // inversion of a positive definite symetric matrix of size n
189 
190  int i, j, k, jj;
191  double r, s;
192  double deter = 0;
193 
194  /* initialisation */
195 
196  if (n > matdim) {
197  printf("ERROR : trying to use TSFit::inverms with size %d( max allowed %d\n", n, matdim);
198  return -999.;
199  }
200 
201  int zero = 0;
202  memset((char *)al, zero, 8 * n * n);
203  memset((char *)be, zero, 8 * n * n);
204  /*
205  for(i=0;i<n;i++){
206  for(j=0;j<n;j++){
207  al[i][j] = 0.;
208  be[i][j] = 0.;
209  }
210  }
211  */
212  /* decomposition en vecteurs sur une base orthonormee */
213  al[0][0] = sqrt(g[0][0]);
214  for (i = 1; i < n; i++) {
215  al[i][0] = g[0][i] / al[0][0];
216  for (j = 1; j <= i; j++) {
217  s = 0.;
218  for (k = 0; k <= j - 1; k++) {
219  s = s + al[i][k] * al[j][k];
220  }
221  r = g[i][j] - s;
222  if (j < i)
223  al[i][j] = r / al[j][j];
224  if (j == i)
225  al[i][j] = sqrt(r);
226  }
227  }
228  /* inversion de la matrice al */
229  be[0][0] = 1. / al[0][0];
230  for (i = 1; i < n; i++) {
231  be[i][i] = 1. / al[i][i];
232  for (j = 0; j < i; j++) {
233  jj = i - j - 1;
234  s = 0.;
235  for (k = jj + 1; k <= i; k++) {
236  s = s + be[i][k] * al[k][jj];
237  }
238  be[i][jj] = -s / al[jj][jj];
239  }
240  }
241  /* calcul de la matrice ginv */
242  for (i = 0; i < n; i++) {
243  for (j = 0; j < n; j++) {
244  s = 0.;
245  for (k = 0; k < n; k++) {
246  s = s + be[k][i] * be[k][j];
247  }
248  ginv[i][j] = s;
249  }
250  }
251 
252  return deter;
253 }
double al[matdim][matdim]
Definition: TSFit.h:36
The Signals That Services Can Subscribe To This is based on ActivityRegistry and is current per Services can connect to the signals distributed by the ActivityRegistry in order to monitor the activity of the application Each possible callback has some defined which we here list in angle e g
Definition: Activities.doc:4
double be[matdim][matdim]
Definition: TSFit.h:36
static constexpr int matdim
Definition: TSFit.h:10
T sqrt(T t)
Definition: SSEVec.h:19

◆ set_params()

void TSFit::set_params ( int  n_samples,
int  niter,
int  n_presmpl,
int  sample_min,
int  sample_max,
double  time_of_max,
double  chi2_max,
int  nsbm,
int  nsam 
)

Definition at line 39 of file TSFit.cc.

References alpha_th, avtm, beta_th, MuonTCMETValueMapProducer_cff::chi2_max, iinf, isup, dqmdumpme::k, n_presamples, n_samples_aft_max, n_samples_bef_max, nbr_iter_fit, nbs, norme, sample_flag, and xki2_max.

Referenced by EcalTestPulseAnalyzer::analyze().

47  {
48  // parameters initialisation of the fit package
49  // nsbm n_samples_bef_max, nsam n_samples_aft_max
50 
51  nbs = n_samples;
52  nbr_iter_fit = niter;
53  n_presamples = n_presmpl;
54  iinf = sample_min;
55  isup = sample_max;
56  avtm = time_of_max;
58  n_samples_bef_max = nsbm;
59  n_samples_aft_max = nsam;
60 
61  norme = 0.;
62  alpha_th = 2.20;
63  beta_th = 1.11;
64 
65  int k;
66 
67  for (k = 0; k <= nbs; k++) {
68  sample_flag[k] = 0;
69  }
70 
71  for (k = sample_min; k <= sample_max; k++) {
72  sample_flag[k] = 1;
73  }
74  /*
75  int lim1 = ( iinf > n_presamples ) ? n_presamples : iinf;
76  for(int k=lim1;k<=sample_max;k++){
77  sample_flag[k] = 2;
78  }
79  */
80  // printf( "sample_fag : " );
81  // for(k=0;k<=nbs;k++){
82  // printf( "%1d ", sample_flag[k] );
83  // }
84  // printf( "\n" );
85 }
double avtm
Definition: TSFit.h:24
int n_samples_bef_max
Definition: TSFit.h:25
double alpha_th
Definition: TSFit.h:33
int nbs
Definition: TSFit.h:21
double xki2_max
Definition: TSFit.h:30
int n_presamples
Definition: TSFit.h:22
double norme
Definition: TSFit.h:30
double beta_th
Definition: TSFit.h:33
int nbr_iter_fit
Definition: TSFit.h:34
int sample_flag[SDIM]
Definition: TSFit.h:48
int n_samples_aft_max
Definition: TSFit.h:26
int isup
Definition: TSFit.h:23
int iinf
Definition: TSFit.h:23

Member Data Documentation

◆ acc

double TSFit::acc[SDIM]
private

Definition at line 52 of file TSFit.h.

◆ adcp

double TSFit::adcp[SDIM]
private

Definition at line 54 of file TSFit.h.

◆ adfmx

double TSFit::adfmx[SDIM]
private

Definition at line 53 of file TSFit.h.

Referenced by fit_third_degree_polynomial().

◆ al

double TSFit::al[matdim][matdim]
private

Definition at line 36 of file TSFit.h.

Referenced by inverms().

◆ alpha_th

double TSFit::alpha_th
private

Definition at line 33 of file TSFit.h.

Referenced by set_params().

◆ avtm

double TSFit::avtm
private

Definition at line 24 of file TSFit.h.

Referenced by set_params().

◆ be

double TSFit::be[matdim][matdim]
private

Definition at line 36 of file TSFit.h.

Referenced by inverms().

◆ beta_th

double TSFit::beta_th
private

Definition at line 33 of file TSFit.h.

Referenced by set_params().

◆ corel

double TSFit::corel[SDIM]
private

Definition at line 56 of file TSFit.h.

◆ cov

double TSFit::cov[matdim][matdim]
private

Definition at line 35 of file TSFit.h.

Referenced by fpol3dg().

◆ der

double TSFit::der[SDIM][5]
private

Definition at line 60 of file TSFit.h.

◆ diminpar

constexpr unsigned int TSFit::diminpar = 10
static

Definition at line 11 of file TSFit.h.

◆ dimoutpar

constexpr unsigned int TSFit::dimoutpar = 10
static

Definition at line 12 of file TSFit.h.

◆ errmat

double TSFit::errmat[SDIM][SDIM]
private

Definition at line 47 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), fpol3dg(), and init_errmat().

◆ f

double TSFit::f[SDIM]
private

◆ ff

double TSFit::ff[SDIM][4]
private

Definition at line 59 of file TSFit.h.

Referenced by fpol3dg().

◆ iinf

int TSFit::iinf
private

Definition at line 23 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

◆ invcov

double TSFit::invcov[matdim][matdim]
private

Definition at line 35 of file TSFit.h.

Referenced by fpol3dg().

◆ isup

int TSFit::isup
private

Definition at line 23 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

◆ maskp3

double TSFit::maskp3[SDIM]
private

Definition at line 55 of file TSFit.h.

Referenced by fit_third_degree_polynomial().

◆ matdim

constexpr int TSFit::matdim = 5
static

Definition at line 10 of file TSFit.h.

Referenced by inverms().

◆ n_presamples

int TSFit::n_presamples
private

Definition at line 22 of file TSFit.h.

Referenced by set_params().

◆ n_samples_aft_max

int TSFit::n_samples_aft_max
private

Definition at line 26 of file TSFit.h.

Referenced by set_params().

◆ n_samples_bef_max

int TSFit::n_samples_bef_max
private

Definition at line 25 of file TSFit.h.

Referenced by set_params().

◆ nbcor

double TSFit::nbcor[SDIM]
private

Definition at line 57 of file TSFit.h.

◆ nbr_iter_fit

int TSFit::nbr_iter_fit
private

Definition at line 34 of file TSFit.h.

Referenced by set_params().

◆ nbs

int TSFit::nbs
private

Definition at line 21 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

◆ nmxu_sto

int TSFit::nmxu_sto
private

Definition at line 32 of file TSFit.h.

◆ norme

double TSFit::norme
private

Definition at line 30 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

◆ npar_moni

constexpr unsigned int TSFit::npar_moni = 4
static

Definition at line 13 of file TSFit.h.

◆ parfp3

double TSFit::parfp3[dimoutpar]
private

Definition at line 40 of file TSFit.h.

Referenced by fit_third_degree_polynomial().

◆ PLSHDIM

constexpr unsigned int TSFit::PLSHDIM = 650
static

Definition at line 9 of file TSFit.h.

◆ plshdim

int TSFit::plshdim

Definition at line 64 of file TSFit.h.

Referenced by TSFit().

◆ sample_flag

int TSFit::sample_flag[SDIM]
private

Definition at line 48 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

◆ SDIM

constexpr unsigned int TSFit::SDIM = 14
static

Definition at line 8 of file TSFit.h.

◆ sdim

int TSFit::sdim

Definition at line 63 of file TSFit.h.

Referenced by init_errmat(), and TSFit().

◆ t

double TSFit::t[SDIM]
private

Definition at line 49 of file TSFit.h.

Referenced by fpol3dg().

◆ xki2_max

double TSFit::xki2_max
private

Definition at line 30 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

◆ z

double TSFit::z[SDIM]
private

Definition at line 50 of file TSFit.h.

Referenced by geometryXMLparser.Alignable::pos(), and ntupleDataFormat._HitObject::r3D().