CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
List of all members | Public Member Functions | 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)
 
virtual ~TSFit ()
 

Public Attributes

int plshdim
 
int sdim
 

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 16 of file TSFit.h.

Constructor & Destructor Documentation

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

Definition at line 20 of file TSFit.cc.

References plshdim, sdim, and findQualityFiles::size.

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

Definition at line 72 of file TSFit.h.

72 {}

Member Function Documentation

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

Definition at line 258 of file TSFit.cc.

References adfmx, errmat, fpol3dg(), i, iinf, isup, gen::k, maskp3, nbs, norme, parfp3, sample_flag, tmp, and xki2_max.

260  {
261  // third degree polynomial fit of the pulse summit.
262  // samples are contained in array bdc and must be pedestal
263  // substracted.
264  // only samples having sample_flag >= 1 are used.
265  // the unit of time is one clock unit, that is to say 25 ns.
266  // output: ret_dat[0] = pulse height
267  // ret_dat[1] position of maximum in the sample frame in clock units
268  // ret_dat[2] adc value of the highest sample
269  // ret_dat[3] number of the highest sample
270  // ret_dat[4] lower sample number used for fitting
271  // ret_dat[5] upper sample number used for fitting
272  // errmat inverse of the error matrix
273 
274  int i;
275  int nus;
276  double xki2;
277  double tm, tmp, amp;
278 
279  static double nevt;
280 
281  ret_dat[0] = -999.;
282  ret_dat[1] = -999.;
283 
284  // search the maximum
285  double val_max = 0.;
286  int imax = 0;
287  for(i=0;i<nbs;i++){
288  if( sample_flag[i] == 0 )continue;
289  if( bdc[i] > val_max ){
290  val_max = bdc[i];
291  imax = i;
292  }
293  }
294 
295  if( (val_max*val_max) * errmat[imax][imax] < 16. )return -118;
296 
297  // if( imax != 9 )printf( "imax : %d !!!!!!!!!!!!!!!!!!!!!!!!!!!\n", imax );
298 
299  if( norme == 0. )norme = val_max;
300 
301  // look for samples above 1/3 of maximum before and 1/2 after
302  double val2 = val_max / 2.;
303  double val3 = val_max / 2.;
304  int ilow = iinf;
305  int ihig = 0;
306 
307  for(i=iinf;i<=isup;i++){
308  if( sample_flag[i] >= 1 ){
309  if( ( bdc[i] < val3 ) && ( i < imax ) )ilow = i;
310  if( bdc[i] > val2 )ihig = i;
311  }
312  }
313 
314  ilow++;
315 
316  //ilow = imax - 1;
317 
318  /* le test suivant, apparemment idiot, est mis a cause des sequences 0. 2048. qui apparaissent dans certains mauvais evts JPP 11/09/00 */
319 
320  if( ihig == ilow)return -105;
321  if( ilow == imax )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. )return -101.;
353  tmp = parfp3[6] + (float)ilow;
354 
355  /*
356  validation of fit quality. Most of the time the fit is done with
357  four samples, therefore there is no possible ki2 check. When more than
358  4 samples are used the ki2 is often bad. So, in order to suppress some
359  events with bad samples, a consistency check on the position of the
360  maximum and minimum of the 3rd degree polynomial is used.
361  */
362 
363  if( xki2 > xki2_max ){
364  return -102.;
365  }
366  if( (tm < (double)ilow ) || (tm > (double)ihig)){
367  return -103.;
368  }
369 
370  if( (tmp > (double)ilow ) && (tmp < (double)ihig - 1.) ){
371  return -104.;
372  }
373 
374  nevt += 1.;
375 
376  ret_dat[0] = amp;
377  ret_dat[1] = tm;
378  ret_dat[2] = val_max;
379  ret_dat[3] = (double)imax;
380  ret_dat[4] = (double)ilow;
381  ret_dat[5] = (double)ihig;
382  ret_dat[6] = (double)tmp;
383  int k;
384  for(i=0;i<4;i++){
385  k=i+7;
386  ret_dat[k] = parfp3[i];
387  }
388 
389  return xki2;
390 }
int i
Definition: DBlmapReader.cc:9
double parfp3[dimoutpar]
Definition: TSFit.h:42
int nbs
Definition: TSFit.h:23
double xki2_max
Definition: TSFit.h:32
double fpol3dg(int, double *, double *, double *)
Definition: TSFit.cc:104
double norme
Definition: TSFit.h:32
int k[5][pyjets_maxn]
double errmat[SDIM][SDIM]
Definition: TSFit.h:50
double maskp3[SDIM]
Definition: TSFit.h:58
double adfmx[SDIM]
Definition: TSFit.h:56
int sample_flag[SDIM]
Definition: TSFit.h:51
std::vector< std::vector< double > > tmp
Definition: MVATrainer.cc:100
int isup
Definition: TSFit.h:25
int iinf
Definition: TSFit.h:25
double TSFit::fpol3dg ( int  nmxul,
double *  parom,
double *  mask,
double *  adc 
)

Definition at line 104 of file TSFit.cc.

References cov, delta, errmat, ff, h, i, invcov, inverms(), gen::k, prof2calltree::l, asciidump::s, mathSSE::sqrt(), t, and tmp.

Referenced by fit_third_degree_polynomial().

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

Definition at line 86 of file TSFit.cc.

References errmat, i, j, and sdim.

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

Definition at line 190 of file TSFit.cc.

References al, be, g, i, j, gen::k, n, csvReporter::r, asciidump::s, mathSSE::sqrt(), and zero.

Referenced by fpol3dg().

190  {
191  // inversion of a positive definite symetric matrix of size n
192 
193  int i, j, k, jj;
194  double r, s;
195  double deter=0;
196 
197  /* initialisation */
198 
199  if( n > matdim ){
200  printf(
201  "ERROR : trying to use TSFit::inverms with size %d( max allowed %d\n",
202  n, matdim );
203  return -999.;
204  }
205 
206  int zero = 0;
207  memset( (char *)al, zero, 8*n*n );
208  memset( (char *)be, zero, 8*n*n );
209  /*
210  for(i=0;i<n;i++){
211  for(j=0;j<n;j++){
212  al[i][j] = 0.;
213  be[i][j] = 0.;
214  }
215  }
216  */
217  /* decomposition en vecteurs sur une base orthonormee */
218  al[0][0] = sqrt( g[0][0] );
219  for(i=1;i<n;i++){
220  al[i][0] = g[0][i] / al[0][0];
221  for(j=1;j<=i;j++){
222  s=0.;
223  for(k=0;k<=j-1;k++){
224  s = s + al[i][k] * al[j][k];
225  }
226  r= g[i][j] - s;
227  if( j < i ) al[i][j] = r / al[j][j];
228  if( j == i ) al[i][j] = sqrt( r );
229  }
230  }
231  /* inversion de la matrice al */
232  be[0][0] = 1./al[0][0];
233  for(i=1;i<n;i++){
234  be[i][i] = 1. / al[i][i];
235  for(j=0;j<i;j++){
236  jj=i-j-1;
237  s=0.;
238  for(k=jj+1;k<=i;k++){
239  s=s+ be[i][k] * al[k][jj];
240  }
241  be[i][jj]=-s/al[jj][jj];
242  }
243  }
244  /* calcul de la matrice ginv */
245  for(i=0;i<n;i++){
246  for(j=0;j<n;j++){
247  s=0.;
248  for(k=0;k<n;k++){
249  s=s+ be[k][i] * be[k][j];
250  }
251  ginv[i][j]=s;
252  }
253  }
254 
255  return deter;
256 }
int i
Definition: DBlmapReader.cc:9
double al[matdim][matdim]
Definition: TSFit.h:38
#define matdim
Definition: TSFit.h:9
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:38
T sqrt(T t)
Definition: SSEVec.h:28
int j
Definition: DBlmapReader.cc:9
int k[5][pyjets_maxn]
string s
Definition: asciidump.py:422
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 42 of file TSFit.cc.

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

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

Member Data Documentation

double TSFit::acc[SDIM]
private

Definition at line 55 of file TSFit.h.

double TSFit::adcp[SDIM]
private

Definition at line 57 of file TSFit.h.

double TSFit::adfmx[SDIM]
private

Definition at line 56 of file TSFit.h.

Referenced by fit_third_degree_polynomial().

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

Definition at line 38 of file TSFit.h.

Referenced by inverms().

double TSFit::alpha_th
private

Definition at line 35 of file TSFit.h.

Referenced by set_params().

double TSFit::avtm
private

Definition at line 26 of file TSFit.h.

Referenced by set_params().

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

Definition at line 38 of file TSFit.h.

Referenced by inverms().

double TSFit::beta_th
private

Definition at line 35 of file TSFit.h.

Referenced by set_params().

double TSFit::corel[SDIM]
private

Definition at line 59 of file TSFit.h.

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

Definition at line 37 of file TSFit.h.

Referenced by fpol3dg().

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

Definition at line 63 of file TSFit.h.

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

Definition at line 50 of file TSFit.h.

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

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

Definition at line 62 of file TSFit.h.

Referenced by fpol3dg().

int TSFit::iinf
private

Definition at line 25 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

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

Definition at line 37 of file TSFit.h.

Referenced by fpol3dg().

int TSFit::isup
private

Definition at line 25 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

double TSFit::maskp3[SDIM]
private

Definition at line 58 of file TSFit.h.

Referenced by fit_third_degree_polynomial().

int TSFit::n_presamples
private

Definition at line 24 of file TSFit.h.

Referenced by set_params().

int TSFit::n_samples_aft_max
private

Definition at line 28 of file TSFit.h.

Referenced by set_params().

int TSFit::n_samples_bef_max
private

Definition at line 27 of file TSFit.h.

Referenced by set_params().

double TSFit::nbcor[SDIM]
private

Definition at line 60 of file TSFit.h.

int TSFit::nbr_iter_fit
private

Definition at line 36 of file TSFit.h.

Referenced by set_params().

int TSFit::nbs
private

Definition at line 23 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

int TSFit::nmxu_sto
private

Definition at line 34 of file TSFit.h.

double TSFit::norme
private

Definition at line 32 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

double TSFit::parfp3[dimoutpar]
private

Definition at line 42 of file TSFit.h.

Referenced by fit_third_degree_polynomial().

int TSFit::plshdim

Definition at line 68 of file TSFit.h.

Referenced by TSFit().

int TSFit::sample_flag[SDIM]
private

Definition at line 51 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

int TSFit::sdim

Definition at line 67 of file TSFit.h.

Referenced by init_errmat(), and TSFit().

double TSFit::t[SDIM]
private

Definition at line 52 of file TSFit.h.

Referenced by fpol3dg().

double TSFit::xki2_max
private

Definition at line 32 of file TSFit.h.

Referenced by fit_third_degree_polynomial(), and set_params().

double TSFit::z[SDIM]
private

Definition at line 53 of file TSFit.h.

Referenced by geometryXMLparser.Alignable::pos().