Changeset 2920 for trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h
- Timestamp:
- 01/26/04 19:32:58 (21 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/mcalib/MHCalibrationBlindPixel.h
r2913 r2920 71 71 72 72 //Getters 73 const Double_t GetLambda() const { return fLambda; }74 const Double_t GetLambdaCheck() const { return fLambdaCheck;}75 const Double_t GetMu0() const { return fMu0; }76 const Double_t GetMu1() const { return fMu1; }77 const Double_t GetSigma0() const { return fSigma0; }78 const Double_t GetSigma1() const { return fSigma1; }79 80 const Double_t GetLambdaErr() const { return fLambdaErr; }81 const Double_t GetLambdaCheckErr() 82 const Double_t GetMu0Err() const { return fMu0Err; }83 const Double_t GetMu1Err() const { return fMu1Err; }84 const Double_t GetSigma0Err() const { return fSigma0Err; }85 const Double_t GetSigma1Err() const { return fSigma1Err; }86 87 const Double_t GetChiSquare() const { return fChisquare; }88 const Double_t GetProb() const { return fProb;}89 const Int_t GetNdf() const { return fNdf;}90 91 const Double_t GetMeanTime() const { return fMeanTime;}92 const Double_t GetMeanTimeErr() const { return fMeanTimeErr; }93 const Double_t GetSigmaTime() const { return fSigmaTime; }94 const Double_t GetSigmaTimeErr() const { return fSigmaTimeErr;}95 96 const Bool_t IsFitOK() { return fFitOK;}97 98 const TH1I *GetHBlindPixelChargevsN() const { return fHBlindPixelChargevsN;}99 const TH1F *GetHBlindPixelPSD() const { return fHBlindPixelPSD;}73 const Double_t GetLambda() const { return fLambda; } 74 const Double_t GetLambdaCheck() const { return fLambdaCheck; } 75 const Double_t GetMu0() const { return fMu0; } 76 const Double_t GetMu1() const { return fMu1; } 77 const Double_t GetSigma0() const { return fSigma0; } 78 const Double_t GetSigma1() const { return fSigma1; } 79 80 const Double_t GetLambdaErr() const { return fLambdaErr; } 81 const Double_t GetLambdaCheckErr() const { return fLambdaCheckErr; } 82 const Double_t GetMu0Err() const { return fMu0Err; } 83 const Double_t GetMu1Err() const { return fMu1Err; } 84 const Double_t GetSigma0Err() const { return fSigma0Err; } 85 const Double_t GetSigma1Err() const { return fSigma1Err; } 86 87 const Double_t GetChiSquare() const { return fChisquare; } 88 const Double_t GetProb() const { return fProb; } 89 const Int_t GetNdf() const { return fNdf; } 90 91 const Double_t GetMeanTime() const { return fMeanTime; } 92 const Double_t GetMeanTimeErr() const { return fMeanTimeErr; } 93 const Double_t GetSigmaTime() const { return fSigmaTime; } 94 const Double_t GetSigmaTimeErr() const { return fSigmaTimeErr; } 95 96 const Bool_t IsFitOK() const { return fFitOK; } 97 98 const TH1I *GetHBlindPixelChargevsN() const { return fHBlindPixelChargevsN; } 99 const TH1F *GetHBlindPixelPSD() const { return fHBlindPixelPSD; } 100 100 101 101 // Draws … … 131 131 inline static Double_t fFitFuncMichele(Double_t *x, Double_t *par) 132 132 { 133 134 const Double_t qe = 0.247; 135 136 Double_t lambda = par[0]; 137 138 Double_t b1 = par[1]; 139 Double_t b2 = par[2]; 140 141 Double_t sum = 0.; 133 134 Double_t lambda1cat = par[0]; 135 Double_t lambda1dyn = par[1]; 136 Double_t mu0 = par[2]; 137 Double_t mu1cat = par[3]; 138 Double_t mu1dyn = par[4]; 139 Double_t sigma0 = par[5]; 140 Double_t sigma1cat = par[6]; 141 Double_t sigma1dyn = par[7]; 142 143 Double_t sumcat = 0.; 144 Double_t sumdyn = 0.; 142 145 Double_t arg = 0.; 143 146 144 return TMath::Exp(-1.*lambda)*par[5]*sum; 145 147 if (mu1cat < mu0) 148 return fNoWay; 149 150 if (sigma1cat < sigma0) 151 return fNoWay; 152 153 // if (sigma1cat < sigma1dyn) 154 // return NoWay; 155 156 //if (mu1cat < mu1dyn) 157 // return NoWay; 158 159 // if (lambda1cat < lambda1dyn) 160 // return NoWay; 161 162 Double_t mu2cat = (2.*mu1cat)-mu0; 163 Double_t mu2dyn = (2.*mu1dyn)-mu0; 164 Double_t mu3cat = (3.*mu1cat)-(2.*mu0); 165 Double_t mu3dyn = (3.*mu1dyn)-(2.*mu0); 166 167 Double_t sigma2cat = TMath::Sqrt((2.*sigma1cat*sigma1cat) - (sigma0*sigma0)); 168 Double_t sigma2dyn = TMath::Sqrt((2.*sigma1dyn*sigma1dyn) - (sigma0*sigma0)); 169 Double_t sigma3cat = TMath::Sqrt((3.*sigma1cat*sigma1cat) - (2.*sigma0*sigma0)); 170 Double_t sigma3dyn = TMath::Sqrt((3.*sigma1dyn*sigma1dyn) - (2.*sigma0*sigma0)); 171 172 Double_t lambda2cat = lambda1cat*lambda1cat; 173 Double_t lambda2dyn = lambda1dyn*lambda1dyn; 174 Double_t lambda3cat = lambda2cat*lambda1cat; 175 Double_t lambda3dyn = lambda2dyn*lambda1dyn; 176 177 // k=0: 178 arg = (x[0] - mu0)/sigma0; 179 sumcat = TMath::Exp(-0.5*arg*arg)/sigma0; 180 sumdyn =sumcat; 181 182 // k=1cat: 183 arg = (x[0] - mu1cat)/sigma1cat; 184 sumcat += lambda1cat*TMath::Exp(-0.5*arg*arg)/sigma1cat; 185 // k=1dyn: 186 arg = (x[0] - mu1dyn)/sigma1dyn; 187 sumdyn += lambda1dyn*TMath::Exp(-0.5*arg*arg)/sigma1dyn; 188 189 // k=2cat: 190 arg = (x[0] - mu2cat)/sigma2cat; 191 sumcat += 0.5*lambda2cat*TMath::Exp(-0.5*arg*arg)/sigma2cat; 192 // k=2dyn: 193 arg = (x[0] - mu2dyn)/sigma2dyn; 194 sumdyn += 0.5*lambda2dyn*TMath::Exp(-0.5*arg*arg)/sigma2dyn; 195 196 197 // k=3cat: 198 arg = (x[0] - mu3cat)/sigma3cat; 199 sumcat += 0.1666666667*lambda3cat*TMath::Exp(-0.5*arg*arg)/sigma3cat; 200 // k=3dyn: 201 arg = (x[0] - mu3dyn)/sigma3dyn; 202 sumdyn += 0.1666666667*lambda3dyn*TMath::Exp(-0.5*arg*arg)/sigma3dyn; 203 204 sumcat = TMath::Exp(-1.*lambda1cat)*sumcat; 205 sumdyn = TMath::Exp(-1.*lambda1dyn)*sumdyn; 206 207 return par[8]*(sumcat+sumdyn)/2.; 208 146 209 } 147 210 … … 339 402 340 403 } 404 405 inline static Double_t fPolya(Double_t *x, Double_t *par) 406 { 407 408 const Double_t QEcat = 0.247; // mean quantum efficiency 409 const Double_t sqrt2 = 1.4142135623731; 410 const Double_t sqrt3 = 1.7320508075689; 411 const Double_t sqrt4 = 2.; 412 413 const Double_t lambda = par[0]; // mean number of photons 414 415 const Double_t excessPoisson = par[1]; // non-Poissonic noise contribution 416 const Double_t delta1 = par[2]; // amplification first dynode 417 const Double_t delta2 = par[3]; // amplification subsequent dynodes 418 419 const Double_t electronicAmpl = par[4]; // electronic amplification and conversion to FADC charges 420 421 const Double_t pmtAmpl = delta1*delta2*delta2*delta2*delta2*delta2; // total PMT gain 422 const Double_t A = 1. + excessPoisson - QEcat 423 + 1./delta1 424 + 1./delta1/delta2 425 + 1./delta1/delta2/delta2; // variance contributions from PMT and QE 426 427 const Double_t totAmpl = QEcat*pmtAmpl*electronicAmpl; // Total gain and conversion 428 429 const Double_t mu0 = par[7]; // pedestal 430 const Double_t mu1 = totAmpl; // single phe position 431 const Double_t mu2 = 2*totAmpl; // double phe position 432 const Double_t mu3 = 3*totAmpl; // triple phe position 433 const Double_t mu4 = 4*totAmpl; // quadruple phe position 434 435 const Double_t sigma0 = par[5]; 436 const Double_t sigma1 = electronicAmpl*pmtAmpl*TMath::Sqrt(QEcat*A); 437 const Double_t sigma2 = sqrt2*sigma1; 438 const Double_t sigma3 = sqrt3*sigma1; 439 const Double_t sigma4 = sqrt4*sigma1; 440 441 const Double_t lambda2 = lambda*lambda; 442 const Double_t lambda3 = lambda2*lambda; 443 const Double_t lambda4 = lambda3*lambda; 444 445 //-- calculate the area---- 446 Double_t arg = (x[0] - mu0)/sigma0; 447 Double_t sum = TMath::Exp(-0.5*arg*arg)/sigma0; 448 449 // k=1: 450 arg = (x[0] - mu1)/sigma1; 451 sum += lambda*TMath::Exp(-0.5*arg*arg)/sigma1; 452 453 // k=2: 454 arg = (x[0] - mu2)/sigma2; 455 sum += 0.5*lambda2*TMath::Exp(-0.5*arg*arg)/sigma2; 456 457 // k=3: 458 arg = (x[0] - mu3)/sigma3; 459 sum += 0.1666666667*lambda3*TMath::Exp(-0.5*arg*arg)/sigma3; 460 461 // k=4: 462 arg = (x[0] - mu4)/sigma4; 463 sum += 0.041666666666667*lambda4*TMath::Exp(-0.5*arg*arg)/sigma4; 464 465 return TMath::Exp(-1.*lambda)*par[6]*sum; 466 } 467 468 341 469 342 470 ClassDef(MHCalibrationBlindPixel, 1) // Histograms from the Calibration Blind Pixel
Note:
See TracChangeset
for help on using the changeset viewer.
