Changeset 12721
- Timestamp:
- 12/13/11 16:18:15 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
fact/tools/rootmacros/gainfit2.C
r12713 r12721 14 14 // - draw gain vs. pixel 15 15 // - find Pixels with deviating gainvalue (e.g. dead , crazy Pixels) 16 // and write their number to log-file (e.g. date_run_gainfit.txt)16 // and write their number to txt-file (e.g. date_run_gainfit.txt) 17 17 // - save histograms into root-file (e.g. date_run_gainfit.root) 18 // - Works also with ROI 300 18 19 // 19 20 //+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ … … 65 66 66 67 //Verbosity for Histograms 67 TString verbos[ 3]={"Q","", ""};68 TString verbos[4]={"0Q","Q", "", "V"}; 68 69 //----------------------------------------------------------------------------- 69 70 // Decleration of Histograms … … 96 97 97 98 int gainfit2( 98 const char * InputRootFileName = "../analysis/fpeak/201111 09_006/20111109_006_fpeak.root",99 const char * InputBaselineFileName = "../analysis/fpeak/201111 09_006/20111109_006_fbsl.root",100 const char * OutTextFileName = "../analysis/fpeak/201111 09_006/20111109_006_gainfit.txt",101 const char * RootOutFileName = "../analysis/fpeak/201111 09_006/20111109_006_gainfit.root",102 bool showHistos = false,103 bool debug = false,104 int verbosityLevel = 3)99 const char * InputRootFileName = "../analysis/fpeak/20111110_005/20111110_005_fpeak.root", 100 const char * InputBaselineFileName = "../analysis/fpeak/20111110_005/20111110_005_fbsl.root", 101 const char * OutTextFileName = "../analysis/fpeak/20111110_005/20111110_005_gainfit.txt", 102 const char * RootOutFileName = "../analysis/fpeak/20111110_005/20111110_005_gainfit.root", 103 bool showHistos = true, 104 bool debug = true, 105 int verbosityLevel = 2) 105 106 { 107 108 //----------------------------------------------------------------------------- 109 // Fit Verbosity 110 //----------------------------------------------------------------------------- 111 112 UInt_t fitverb = verbosityLevel; 113 if (verbosityLevel >= 3) fitverb =3; 106 114 107 115 //----------------------------------------------------------------------------- … … 127 135 128 136 //Define Fit Functions 129 Double_t par[9]; 137 Double_t par[6]; 138 Double_t par2[6]; 130 139 TF1 *g1 = new TF1("g1","gaus",5,15); 131 140 TF1 *g2 = new TF1("g2","gaus",15,25); 132 TF1 *g3 = new TF1("g3","gaus",25,35);141 //TF1 *g3 = new TF1("g3","gaus",-5,5); 133 142 TF1 *g4 = new TF1("g4","gaus"); 134 143 //TF1 *total = new TF1("total","gaus(0)+gaus(3)+gaus(6)",7,60); 135 TF1 *fspek_fit = new TF1("fspek_fit", spek_fit, 5, 60, 5); 136 TF1 *fspek_fit2 = new TF1("fspek_fit2", spek_fit2, 5, 60, 4); 137 144 TF1 *fspek_fit = new TF1("fspek_fit", spek_fit, 5, 60, 6); 145 fspek_fit->SetParNames("1.Max", "Gain", "Sigma", "2.Max", "Baseline", "Xtalk"); 146 TF1 *fspek_fit2 = new TF1("fspek_fit2g3", spek_fit2, 5, 60, 4); 147 fspek_fit2->SetParNames("1.Max", "Gain", "Sigma", "2.Max"); 138 148 //Create Canvases 139 149 TCanvas *cGainFit = new TCanvas("cGainFit", "Fit of Amplitude Spektra", 0,0, 1200,400); … … 211 221 //Amplitude Spectrum of a Pixel 212 222 hPixelAmplSpek[ pixel ] = new TH1D; 213 hPixelAmplSpek[ pixel ] = hAmplSpek->ProjectionY( histotitle , pixel+ 0, pixel+0); //Projectipon of each Pixel out of Ampl.Spectrum223 hPixelAmplSpek[ pixel ] = hAmplSpek->ProjectionY( histotitle , pixel+1, pixel+1); //Projectipon of each Pixel out of Ampl.Spectrum 214 224 hPixelAmplSpek[ pixel ]->SetAxisRange(0, 60, "X"); 215 225 … … 234 244 gPad->SetLogy(1); 235 245 236 hPixelAmplSpek[ pixel ]->Fit(g1,verbos[ verbosityLevel-1]+"R");237 hPixelAmplSpek[ pixel ]->Fit(g2,verbos[ verbosityLevel-1]+"0R+");238 hPixelAmplSpek[ pixel ]->Fit(g3,verbos[verbosityLevel-1]+"0R+");246 hPixelAmplSpek[ pixel ]->Fit(g1,verbos[fitverb]+"R"); 247 hPixelAmplSpek[ pixel ]->Fit(g2,verbos[fitverb]+"0R+"); 248 //hPixelAmplSpek[ pixel ]->Fit(g3,verbos[fitverb]+"0R+"); 239 249 g1->GetParameters(&par[0]); 250 g1->GetParameters(&par2[0]); 240 251 g2->GetParameters(&par[3]); 241 g3->GetParameters(&par[6]); 252 g2->GetParameters(&par2[3]); 253 254 //g3->GetParameters(&par[6]); 242 255 243 256 //total->SetParameters(par); 244 257 245 //hPixelAmplSpek[ pixel ]->Fit(total,verbos[ verbosityLevel-1]+"R+");258 //hPixelAmplSpek[ pixel ]->Fit(total,verbos[fitverb]+"R+"); 246 259 //channel[ pixel ].amplSD = total->GetParameter(4) - total->GetParameter(1); 247 260 //channel[ pixel ].amplDT = total->GetParameter(7) - total->GetParameter(4); 248 261 //channel[ pixel ].amplMean = fspek_fit->GetParameter(1); 249 262 250 par[4] = -1; 263 par2[4]=par[4] = hbaseline->GetMean(); 264 cout << "starvalue of Baseline:" << par[4] << endl; 251 265 fspek_fit->SetParameters(par); 252 cout << "spek_fit" << endl;253 hPixelAmplSpek[ pixel ]->Fit(fspek_fit,verbos[ verbosityLevel-1]+"R+");266 if (verbosityLevel > 2) cout << "spek_fit" << endl; 267 hPixelAmplSpek[ pixel ]->Fit(fspek_fit,verbos[fitverb]+"R+"); 254 268 channel[ pixel ].ampl_fspek_fit = fspek_fit->GetParameter(1); 255 269 channel[ pixel ].bsl = fspek_fit->GetParameter(4); 256 cout << "Baseline:" << channel[ pixel ].bsl << endl;257 258 fspek_fit2->SetParameters(par );259 cout << "fspek_fit2" << endl;260 hPixelAmplSpek[ pixel ]->Fit(fspek_fit2,verbos[ verbosityLevel-1]+"R+");270 if (verbosityLevel > 2) cout << "Baseline:" << channel[ pixel ].bsl << endl; 271 272 fspek_fit2->SetParameters(par2); 273 if (verbosityLevel > 2) cout << "fspek_fit2" << endl; 274 hPixelAmplSpek[ pixel ]->Fit(fspek_fit2,verbos[fitverb]+"R+"); 261 275 channel[ pixel ].ampl_fspek_fit2 = fspek_fit2->GetParameter(1); 262 276 … … 278 292 g4->SetLineColor(2); 279 293 hAmplDistribution->Fill( channel[ pixel ].ampl_fspek_fit ) ; 280 hAmplDistribution->Fit(g4, verbos[ verbosityLevel-1]);294 hAmplDistribution->Fit(g4, verbos[fitverb]); 281 295 distribution.sigma = g4->GetParameter(2); 282 296 distribution.mean = g4->GetParameter(1); … … 287 301 g4->SetLineColor(kBlue); 288 302 hAmplDistribution2->Fill(channel[ pixel ].ampl_fspek_fit2) ; 289 hAmplDistribution2->Fit(g4, verbos[ verbosityLevel-1]);303 hAmplDistribution2->Fit(g4, verbos[fitverb]); 290 304 291 305 cGainFit->cd(3); … … 493 507 Double_t peak3 = cross*cross*par[0]*TMath::Exp(-0.5 * arg3 * arg3 /3); 494 508 Double_t peak4 = cross*cross*cross*par[0]*TMath::Exp(-0.5 * arg4 * arg4 /4); 509 par[5]=cross; 495 510 return peak1 + peak2 + peak3 + peak4; 496 511 } … … 500 515 //----------------------------------------------------------------------------- 501 516 502 Double_t spek_fit2(Double_t *x, Double_t *par ){517 Double_t spek_fit2(Double_t *x, Double_t *par2){ 503 518 504 519 Double_t arg = 0; … … 506 521 Double_t arg3 = 0; 507 522 Double_t arg4 = 0; 508 Double_t cross = 1;509 if (par [0] != 0) cross = par[3]/par[0];510 if (par [2] != 0) arg = (x[0] - par[1])/par[2];511 if (par [2] != 0) arg2 = (x[0] - 2*par[1])/par[2];512 if (par [2] != 0) arg3 = (x[0] - 3*par[1])/par[2];513 if (par [2] != 0) arg4 = (x[0] - 4*par[1])/par[2];514 Double_t peak1 = par [0]*TMath::Exp(-0.5 * arg * arg);515 Double_t peak2 = cross*par [0]*TMath::Exp(-0.5 * arg2 * arg2 /2);516 Double_t peak3 = cross*cross*par [0]*TMath::Exp(-0.5 * arg3 * arg3 /3);517 Double_t peak4 = cross*cross*cross*par [0]*TMath::Exp(-0.5 * arg4 * arg4 /4);523 Double_t cross = 0; 524 if (par2[0] != 0) cross = par2[3]/par2[0]; 525 if (par2[2] != 0) arg = (x[0] - par2[1])/par2[2]; 526 if (par2[2] != 0) arg2 = (x[0] - 2*par2[1])/par2[2]; 527 if (par2[2] != 0) arg3 = (x[0] - 3*par2[1])/par2[2]; 528 if (par2[2] != 0) arg4 = (x[0] - 4*par2[1])/par2[2]; 529 Double_t peak1 = par2[0]*TMath::Exp(-0.5 * arg * arg); 530 Double_t peak2 = cross*par2[0]*TMath::Exp(-0.5 * arg2 * arg2 /2); 531 Double_t peak3 = cross*cross*par2[0]*TMath::Exp(-0.5 * arg3 * arg3 /3); 532 Double_t peak4 = cross*cross*cross*par2[0]*TMath::Exp(-0.5 * arg4 * arg4 /4); 518 533 return peak1 + peak2 + peak3 + peak4; 519 534 } … … 530 545 } 531 546 else{ 532 if (verbosityLevel > 2) cout << "Not normal behavio truer in Pixel: " << pixel +1 << endl;547 if (verbosityLevel > 2) cout << "Not normal behaviour in Pixel: " << pixel +1 << endl; 533 548 return true; 534 549 }
Note:
See TracChangeset
for help on using the changeset viewer.