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
PulseFitWithShape Class Reference

#include <PulseFitWithShape.h>

Inheritance diagram for PulseFitWithShape:

Public Member Functions

virtual double doFit (double *, double *cova=0)
 
double getAmpl ()
 
double getTime ()
 
virtual void init (int, int, int, int, int, std::vector< double >, double)
 
 PulseFitWithShape ()
 
virtual ~PulseFitWithShape ()
 

Public Attributes

double fAmp_fitted_max
 
double fTim_fitted_max
 

Private Attributes

std::vector< double > dshape
 
int fNb_iter
 
double fNoise
 
int fNsamples
 
int fNsamplesShape
 
int fNum_samp_after_max
 
int fNum_samp_bef_max
 
std::vector< double > pshape
 

Detailed Description

Definition at line 10 of file PulseFitWithShape.h.

Constructor & Destructor Documentation

PulseFitWithShape::PulseFitWithShape ( )

Definition at line 21 of file PulseFitWithShape.cc.

References fNoise, fNsamples, fNsamplesShape, fNum_samp_after_max, and fNum_samp_bef_max.

22 {
23 
24  fNsamples = 0;
25  fNsamplesShape = 0;
26  fNum_samp_bef_max = 0;
27  fNum_samp_after_max = 0;
28  fNoise = 0.0;
29 }
PulseFitWithShape::~PulseFitWithShape ( )
virtual

Definition at line 32 of file PulseFitWithShape.cc.

33 {
34 }

Member Function Documentation

double PulseFitWithShape::doFit ( double *  adc,
double *  cova = 0 
)
virtual

Definition at line 63 of file PulseFitWithShape.cc.

References trackerHits::c, dshape, fAmp_fitted_max, fNb_iter, fNoise, fNsamples, fNsamplesShape, fNum_samp_after_max, fNum_samp_bef_max, fTim_fitted_max, i, j, pshape, q1, q2, and indexGen::s2.

64 {
65 
66  // xpar = fit paramaters
67  // [0] = signal amplitude
68  // [1] = residual pedestal
69  // [2] = clock phase
70 
71  bool useCova=true;
72  if(cova==0) useCova=false;
73 
74  double xpar[3];
75  double chi2;
76 
77  fAmp_fitted_max = 0. ;
78  fTim_fitted_max = 0. ;
79 
80  // for now don't fit pedestal
81 
82  xpar[1]=0.0;
83 
84  // Sample noise. If the cova matrix is defined, use it :
85 
86  double noise=fNoise;
87  //if(cova[0] > 0.) noise=1./sqrt(cova[0]);
88 
89  xpar[0]=0.;
90  xpar[2]=0.;
91 
92 
93  // first locate max:
94 
95  int imax=0;
96  double amax=0.0;
97  for(int i=0; i<fNsamples; i++){
98  if (adc[i]>amax){
99  amax=adc[i];
100  imax=i;
101  }
102  }
103 
104  // Shift pulse shape and calculate its derivative:
105 
106  double qms=0.;
107  int ims=0;
108 
109  // 1) search for maximum
110 
111  for(int is=0; is<fNsamplesShape; is++){
112  if(pshape[is] > qms){
113  qms=pshape[is];
114  ims=is;
115  }
116 
117  // 2) compute shape derivative :
118 
119  if(is < fNsamplesShape-2)
120  dshape[is]= (pshape[is+2]-pshape[is])*12.5;
121  else
122  dshape[is]=dshape[is-1];
123  }
124 
125  // 3) compute pol2 max
126 
127  double sq1=pshape[ims-1];
128  double sq2=pshape[ims];
129  double sq3=pshape[ims+1];
130 
131  double a2=(sq3+sq1)/2.0-sq2;
132  double a1=sq2-sq1+a2*(1-2*ims);
133 
134 
135  double t_ims=0;
136  if(a2!=0) t_ims=-a1/(2.0*a2);
137 
138 
139  // Starting point of the fit : qmax and tmax given by a
140  // P2 fit on 3 max samples.
141 
142  double qm=0.;
143  int im=0;
144 
145  int nsamplef=fNum_samp_bef_max + fNum_samp_after_max +1 ; // number of samples used in the fit
146  int nsampleo=imax-fNum_samp_bef_max; // first sample number = sample max-fNum_samp_bef_max
147 
148  for(int is=0; is<nsamplef; is++){
149 
150  if(adc[is+nsampleo] > qm){
151  qm=adc[is+nsampleo];
152  im=nsampleo+is;
153  }
154  }
155 
156  double tm;
157  if(qm > 5.*noise){
158  if(im >= nsamplef+nsampleo) im=nsampleo+nsamplef-1;
159  double q1=adc[im-1];
160  double q2=adc[im];
161  double q3=adc[im+1];
162  double y2=(q1+q3)/2.-q2;
163  double y1=q2-q1+y2*(1-2*im);
164  double y0=q2-y1*(double)im-y2*(double)(im*im);
165  tm=-y1/2./y2;
166  xpar[0]=y0+y1*tm+y2*tm*tm;
167  xpar[2]=(double)ims/25.-tm;
168  }
169 
170  double chi2old=999999.;
171  chi2=99999.;
172  int nsfit=nsamplef;
173  int iloop=0;
174  int nloop=fNb_iter;
175  if(qm <= 5*noise)nloop=1; // Just one iteration for very low signal
176 
177  double *resi;
178  resi= new double[fNsamples];
179  for (int j=0;j<fNsamples;j++) resi[j]=0;
180 
181  while(std::fabs(chi2old-chi2) > 0.1 && iloop < nloop)
182  {
183  iloop++;
184  chi2old=chi2;
185 
186  double c=0.;
187  double y1=0.;
188  double s1=0.;
189  double s2=0.;
190  double ys1=0.;
191  double sp1=0.;
192  double sp2=0.;
193  double ssp=0.;
194  double ysp=0.;
195 
196  for(int is=0; is<nsfit; is++)
197  {
198  int iis=is;
199  iis=is+nsampleo;
200 
201  double t1=(double)iis+xpar[2];
202  double xbin=t1*25.;
203  int ibin1=(int)xbin;
204 
205  if(ibin1 < 0) ibin1=0;
206 
207  double amp1,amp11,amp12,der1,der11,der12;
208 
209  if(ibin1 <= fNsamplesShape-2){ // Interpolate shape values to get the right number :
210 
211  int ibin2=ibin1+1;
212  double xfrac=xbin-ibin1;
213  amp11=pshape[ibin1];
214  amp12=pshape[ibin2];
215  amp1=(1.-xfrac)*amp11+xfrac*amp12;
216  der11=dshape[ibin1];
217  der12=dshape[ibin2];
218  der1=(1.-xfrac)*der11+xfrac*der12;
219 
220  }else{ // Or extraoplate if we are outside the array :
221 
222  amp1=pshape[fNsamplesShape-1]+dshape[fNsamplesShape-1]*
223  (xbin-double(fNsamplesShape-1))/25.;
224  der1=dshape[fNsamplesShape-1];
225  }
226 
227  if( useCova ){ // Use covariance matrix:
228  for(int js=0; js<nsfit; js++)
229  {
230  int jjs=js;
231  jjs+=nsampleo;
232 
233  t1=(double)jjs+xpar[2];
234  xbin=t1*25.;
235  ibin1=(int)xbin;
236  if(ibin1 < 0) ibin1=0;
237  if(ibin1 > fNsamplesShape-2)ibin1=fNsamplesShape-2;
238  int ibin2=ibin1+1;
239  double xfrac=xbin-ibin1;
240  amp11=pshape[ibin1];
241  amp12=pshape[ibin2];
242  double amp2=(1.-xfrac)*amp11+xfrac*amp12;
243  der11=dshape[ibin1];
244  der12=dshape[ibin2];
245  double der2=(1.-xfrac)*der11+xfrac*der12;
246  c=c+cova[js*fNsamples+is];
247  y1=y1+adc[iis]*cova[js*fNsamples+is];
248  s1=s1+amp1*cova[js*fNsamples+is];
249  s2=s2+amp1*amp2*cova[js*fNsamples+is];
250  ys1=ys1+adc[iis]*amp2*cova[js*fNsamples+is];
251  sp1=sp1+der1*cova[is*fNsamples+js];
252  sp2=sp2+der1*der2*cova[js*fNsamples+is];
253  ssp=ssp+amp1*der2*cova[js*fNsamples+is];
254  ysp=ysp+adc[iis]*der2*cova[js*fNsamples+is];
255  }
256  }else { // Don't use covariance matrix:
257  c++;
258  y1=y1+adc[iis];
259  s1=s1+amp1;
260  s2=s2+amp1*amp1;
261  ys1=ys1+amp1*adc[iis];
262  sp1=sp1+der1;
263  sp2=sp2+der1*der1;
264  ssp=ssp+der1*amp1;
265  ysp=ysp+der1*adc[iis];
266  }
267  }
268  xpar[0]=(ysp*ssp-ys1*sp2)/(ssp*ssp-s2*sp2);
269  xpar[2]+=(ysp/xpar[0]/sp2-ssp/sp2);
270 
271  for(int is=0; is<nsfit; is++){
272  int iis=is;
273  iis=is+nsampleo;
274 
275  double t1=(double)iis+xpar[2];
276  double xbin=t1*25.;
277  int ibin1=(int)xbin;
278  if(ibin1 < 0) ibin1=0;
279 
280  if(ibin1 < 0) ibin1=0;
281  if(ibin1 > fNsamplesShape-2)ibin1=fNsamplesShape-2;
282  int ibin2=ibin1+1;
283  double xfrac=xbin-ibin1;
284  double amp11=xpar[0]*pshape[ibin1];
285  double amp12=xpar[0]*pshape[ibin2];
286  double amp1=xpar[1]+(1.-xfrac)*amp11+xfrac*amp12;
287  resi[iis]=adc[iis]-amp1;
288  }
289 
290  chi2=0.;
291  for(int is=0; is<nsfit; is++){
292  int iis=is;
293  iis+=nsampleo;
294 
295  if( useCova ){
296  for(int js=0; js<nsfit; js++){
297  int jjs=js;
298  jjs+=nsampleo;
299  chi2+=resi[iis]*resi[jjs]*cova[js*fNsamples+is];
300  }
301 
302  }else chi2+=resi[iis]*resi[iis];
303  }
304  }
305 
306  fAmp_fitted_max = xpar[0];
307  fTim_fitted_max = (double)t_ims/25.-xpar[2];
308 
309  return chi2 ;
310 
311 }
ecalMGPA::adc
int adc(sample_type sample)
get the ADC sample (12 bits)
Definition: EcalMGPASample.h:11
PulseFitWithShape::pshape
std::vector< double > pshape
Definition: PulseFitWithShape.h:39
i
int i
Definition: DBlmapReader.cc:9
q2
double q2[4]
Definition: TauolaWrapper.h:88
indexGen.s2
tuple s2
Definition: indexGen.py:106
j
int j
Definition: DBlmapReader.cc:9
q1
double q1[4]
Definition: TauolaWrapper.h:87
PulseFitWithShape::dshape
std::vector< double > dshape
Definition: PulseFitWithShape.h:40
double PulseFitWithShape::getAmpl ( )
inline

Definition at line 29 of file PulseFitWithShape.h.

References fAmp_fitted_max.

29 { return fAmp_fitted_max; }
double PulseFitWithShape::getTime ( )
inline

Definition at line 30 of file PulseFitWithShape.h.

References fTim_fitted_max.

30 { return fTim_fitted_max; }
void PulseFitWithShape::init ( int  n_samples,
int  samplb,
int  sampla,
int  niter,
int  n_samplesShape,
std::vector< double >  shape,
double  nois 
)
virtual

Definition at line 38 of file PulseFitWithShape.cc.

References gather_cfg::cout, dshape, fNb_iter, fNoise, fNsamples, fNsamplesShape, fNum_samp_after_max, fNum_samp_bef_max, i, pshape, and hitfit::return.

39 {
40 
41  fNsamples = n_samples ;
42  fNsamplesShape = n_samplesShape ;
43  fNb_iter = niter ;
44  fNum_samp_bef_max = samplb ;
45  fNum_samp_after_max = sampla ;
46 
47 
48  if( fNsamplesShape < fNum_samp_bef_max+fNum_samp_after_max+1){
49  std::cout<<"PulseFitWithShape::init: Error! Configuration samples in fit greater than total number of samples!" << std::endl;
50  }
51 
52  for(int i=0;i<fNsamplesShape;i++){
53  pshape.push_back(shape[i]);
54  dshape.push_back(0.0);
55  }
56 
57  fNoise=nois;
58  return ;
59  }
PulseFitWithShape::pshape
std::vector< double > pshape
Definition: PulseFitWithShape.h:39
i
int i
Definition: DBlmapReader.cc:9
PulseFitWithShape::dshape
std::vector< double > dshape
Definition: PulseFitWithShape.h:40
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121

Member Data Documentation

std::vector< double > PulseFitWithShape::dshape
private

Definition at line 40 of file PulseFitWithShape.h.

Referenced by doFit(), and init().

double PulseFitWithShape::fAmp_fitted_max

Definition at line 26 of file PulseFitWithShape.h.

Referenced by doFit(), and getAmpl().

int PulseFitWithShape::fNb_iter
private

Definition at line 42 of file PulseFitWithShape.h.

Referenced by doFit(), and init().

double PulseFitWithShape::fNoise
private

Definition at line 37 of file PulseFitWithShape.h.

Referenced by doFit(), init(), and PulseFitWithShape().

int PulseFitWithShape::fNsamples
private

Definition at line 35 of file PulseFitWithShape.h.

Referenced by doFit(), init(), and PulseFitWithShape().

int PulseFitWithShape::fNsamplesShape
private

Definition at line 36 of file PulseFitWithShape.h.

Referenced by doFit(), init(), and PulseFitWithShape().

int PulseFitWithShape::fNum_samp_after_max
private

Definition at line 44 of file PulseFitWithShape.h.

Referenced by doFit(), init(), and PulseFitWithShape().

int PulseFitWithShape::fNum_samp_bef_max
private

Definition at line 43 of file PulseFitWithShape.h.

Referenced by doFit(), init(), and PulseFitWithShape().

double PulseFitWithShape::fTim_fitted_max

Definition at line 27 of file PulseFitWithShape.h.

Referenced by doFit(), and getTime().

std::vector< double > PulseFitWithShape::pshape
private

Definition at line 39 of file PulseFitWithShape.h.

Referenced by doFit(), and init().

Generated for CMSSW Reference Manual by doxygen 1.8.5