|
|
|
#include <RecoLocalCalo/EcalRecAlgos/interface/EcalUncalibRecHitRecAnalFitAlgo.h>
Public Member Functions | |
| virtual EcalUncalibratedRecHit | makeRecHit (const C &dataFrame, const double *pedestals, const double *gainRatios, const EcalWeightSet::EcalWeightMatrix **weights, const EcalWeightSet::EcalChi2WeightMatrix **chi2Matrix) |
| Compute parameters. | |
| virtual | ~EcalUncalibRecHitRecAnalFitAlgo () |
Private Member Functions | |
| double | pedestalFunction (double *var, double *par) |
| double | pulseShapeFunction (double *var, double *par) |
Definition at line 25 of file EcalUncalibRecHitRecAnalFitAlgo.h.
| virtual EcalUncalibRecHitRecAnalFitAlgo< C >::~EcalUncalibRecHitRecAnalFitAlgo | ( | ) | [inline, virtual] |
Definition at line 46 of file EcalUncalibRecHitRecAnalFitAlgo.h.
00048 : 00049 // destructor 00050 virtual ~EcalUncalibRecHitRecAnalFitAlgo<C>() { };
| virtual EcalUncalibratedRecHit EcalUncalibRecHitRecAnalFitAlgo< C >::makeRecHit | ( | const C & | dataFrame, | |
| const double * | pedestals, | |||
| const double * | gainRatios, | |||
| const EcalWeightSet::EcalWeightMatrix ** | weights, | |||
| const EcalWeightSet::EcalChi2WeightMatrix ** | chi2Matrix | |||
| ) | [inline, virtual] |
Compute parameters.
Implements EcalUncalibRecHitRecAbsAlgo< C >.
Definition at line 54 of file EcalUncalibRecHitRecAnalFitAlgo.h.
Referenced by EcalAnalFitUncalibRecHitProducer::produce().
00058 { 00059 double amplitude_(-1.), pedestal_(-1.), jitter_(-1.), chi2_(-1.); 00060 00061 // Get time samples 00062 //HepMatrix frame(C::MAXSAMPLES, 1); 00063 double frame[C::MAXSAMPLES]; 00064 // int gainId0 = dataFrame.sample(0).gainId(); 00065 int gainId0 = 1; 00066 int iGainSwitch = 0; 00067 double maxsample(-1); 00068 int imax(-1); 00069 00070 for(int iSample = 0; iSample < C::MAXSAMPLES; iSample++) { 00071 int gainId = dataFrame.sample(iSample).gainId(); 00072 if ( gainId == 0 ) gainId = 3; 00073 if (gainId != gainId0) iGainSwitch++ ; 00074 if (!iGainSwitch) 00075 frame[iSample] = double(dataFrame.sample(iSample).adc()); 00076 else 00077 frame[iSample] = double(((double)(dataFrame.sample(iSample).adc()) - pedestals[gainId-1]) * gainRatios[gainId-1]); 00078 00079 if( frame[iSample]>maxsample ) { 00080 maxsample= frame[iSample]; 00081 imax=iSample; 00082 } 00083 } 00084 00085 // Compute parameters 00086 //std::cout << "EcalUncalibRecHitRecAnalFitAlgo::makeRecHit() not yey implemented. returning dummy rechit" << std::endl; 00087 00088 // prepare TGraph for analytic fit 00089 double xarray[10]={0.,1.,2.,3.,4.,5.,6.,7.,8.,9.}; 00090 TGraph graph(10,xarray,frame); 00091 00092 // fit functions 00093 TF1 pulseShape = TF1("pulseShape", 00094 "[0]*pow((x - [3])/[1],[2])*exp(-[2]*(x - [1] - [3])/[1])", 00095 imax-1.,imax+3.); 00096 TF1 pedestal = TF1("pedestal","[0]",0.,2.); 00097 00098 //TF1 pulseShape = TF1("pulseShape",pulseShapeFunction,imax-1.,imax+3.); 00099 //TF1 pedestal = TF1("pedestal",pedestalFunction,0.,2.); 00100 TF1 pluseAndPed = TF1("pulseAndPed","pedestal+pulseShape"); 00101 00102 //pulseShape parameters 00103 // Amplitude 00104 double FIT_A=(double)maxsample; //Amplitude 00105 pulseShape.SetParameter(0,FIT_A); 00106 pulseShape.SetParName(0,"Amplitude"); 00107 // T peak 00108 double FIT_Tp=(double)imax; //T peak 00109 pulseShape.SetParameter(1,FIT_Tp); 00110 pulseShape.SetParName(1,"t_{P}"); 00111 // Alpha 00112 double FIT_ALFA=1.5; //Alpha 00113 pulseShape.SetParameter(2,FIT_ALFA); 00114 pulseShape.SetParName(2,"\\alpha"); 00115 // T off 00116 double FIT_To=3.; //T off 00117 pulseShape.SetParameter(3,FIT_To); 00118 pulseShape.SetParName(3,"t_{0}"); 00119 00120 // pedestal 00121 pedestal.SetParameter(0,frame[0]); 00122 pedestal.SetParName(0,"Pedestal"); 00123 00124 00125 00126 graph.Fit("pulseShape","QRM"); 00127 //TF1 *pulseShape2=graph.GetFunction("pulseShape"); 00128 00129 if ( std::string(gMinuit->fCstatu.Data()) == std::string("CONVERGED ") ) { 00130 00131 double amplitude_value=pulseShape.GetParameter(0); 00132 00133 graph.Fit("pedestal","QRL"); 00134 //TF1 *pedestal2=graph.GetFunction("pedestal"); 00135 double pedestal_value=pedestal.GetParameter(0); 00136 00137 if (!iGainSwitch) 00138 amplitude_ = amplitude_value - pedestal_value; 00139 else 00140 amplitude_ = amplitude_value; 00141 00142 pedestal_ = pedestal_value; 00143 jitter_ = pulseShape.GetParameter(3); 00144 chi2_ = 1.; // successful fit 00145 00146 /* 00147 std::cout << "separate fits\nA: " << amplitude_value << ", Ped: " << pedestal_value 00148 << ", t0: " << jitter_ << ", tp: " << pulseShape.GetParameter(1) 00149 << ", alpha: " << pulseShape.GetParameter(2) 00150 << std::endl; 00151 */ 00152 00153 } 00154 00155 return EcalUncalibratedRecHit( dataFrame.id(), amplitude_, pedestal_, jitter_ - 6, chi2_); 00156 }
| double EcalUncalibRecHitRecAnalFitAlgo< C >::pedestalFunction | ( | double * | var, | |
| double * | par | |||
| ) | [inline, private] |
Definition at line 42 of file EcalUncalibRecHitRecAnalFitAlgo.h.
00042 { 00043 double x = var[0]; 00044 double ped = par[0]; 00045 return ped; 00046 };
| double EcalUncalibRecHitRecAnalFitAlgo< C >::pulseShapeFunction | ( | double * | var, | |
| double * | par | |||
| ) | [inline, private] |
Definition at line 31 of file EcalUncalibRecHitRecAnalFitAlgo.h.
00031 { 00032 double x = var[0]; 00033 double ampl = par[0]; 00034 double tp = par[1]; 00035 double alpha = par[2]; 00036 double t0 = par[3]; 00037 00038 double f = pow( (x-t0)/tp , alpha ) * exp( -alpha*(x-tp-t0)/tp ); 00039 return ampl*f; 00040 };
1.5.4