Changeset 2923
- Timestamp:
- 01/26/04 22:42:42 (21 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/macros/calibration.C
r2906 r2923 23 23 \* ======================================================================== */ 24 24 25 const TString pedfile = "/mnt/Data/rootdata/Miscellaneous/2003_12_19/20031218_03522_P_Park_E.root"; 26 const TString calfile = "/mnt/Data/rootdata/Miscellaneous/2003_12_19/20031218_03527_C_Park_E.root"; 25 const TString pedfile = "/mnt/Data/rootdata/Mrk421/2004_01_26/20040125_12094_P_Mrk421_E.root"; 26 const TString calfile = "/mnt/Data/rootdata/Mrk421/2004_01_26/20040125_1211*_C_Mrk421_E.root"; 27 28 //const TString pedfile = "/mnt/Data/rootdata/Miscellaneous/2003_12_19/20031218_03522_P_Park_E.root"; 29 //const TString calfile = "/mnt/Data/rootdata/Miscellaneous/2003_12_19/20031218_03527_C_Park_E.root"; 27 30 28 31 void calibration(TString pedname=pedfile, … … 114 117 115 118 // 119 // Making the step size a bit bigger, gives us 120 // faster results 121 // 122 timecalc.SetStepSize(0.5); 123 124 // 116 125 // As long, as we don't have digital modules, 117 126 // we have to set the color of the pulser LED by hand … … 148 157 // 149 158 MCalibrationBlindPix *bp = calcam.GetBlindPixel(); 159 // bp->ChangeFitFunc(MHCalibrationBlindPixel::kEPolya); 150 160 bp->ChangeFitFunc(MHCalibrationBlindPixel::kEPoisson5); 151 161 … … 156 166 // of the arrival times using a spline, uncomment the next line 157 167 // 158 //tlist2.AddToList(&timecalc);168 tlist2.AddToList(&timecalc); 159 169 tlist2.AddToList(&calcalc); 160 170 … … 164 174 MEvtLoop evtloop2; 165 175 evtloop2.SetParList(&plist2); 166 176 167 177 // 168 178 // Execute second analysis … … 183 193 calcam[17].DrawClone(); 184 194 185 MHCamEvent camevt;186 195 MHCamera disp1 (geomcam, "MCalibrationPix;Charge", "Fitted Mean Charges"); 187 196 MHCamera disp3 (geomcam, "MCalibrationPix;SigmaCharge", "Sigma of Fitted Charges"); … … 258 267 c1.Divide(2,3); 259 268 260 CamDraw(c1,disp1, &calcam,1,2,1);261 CamDraw(c1,disp3, &calcam,2,2,2);269 CamDraw(c1,disp1,calcam,1,2,1); 270 CamDraw(c1,disp3,calcam,2,2,1); 262 271 263 272 // Fit Probability … … 265 274 c2.Divide(1,3); 266 275 267 CamDraw(c2,disp5, &calcam,1,1,3);276 CamDraw(c2,disp5,calcam,1,1,3); 268 277 269 278 // Times … … 271 280 c3.Divide(3,3); 272 281 273 CamDraw(c3,disp6, &calcam,1,3,1);274 CamDraw(c3,disp7, &calcam,2,3,0);275 CamDraw(c3,disp8, &calcam,3,3,0);282 CamDraw(c3,disp6,calcam,1,3,1); 283 CamDraw(c3,disp7,calcam,2,3,0); 284 CamDraw(c3,disp8,calcam,3,3,0); 276 285 277 286 // Pedestals … … 279 288 c4.Divide(2,3); 280 289 281 CamDraw(c4,disp9, &calcam,1,2,0);282 CamDraw(c4,disp10, &calcam,2,2,1);290 CamDraw(c4,disp9,calcam,1,2,0); 291 CamDraw(c4,disp10,calcam,2,2,1); 283 292 284 293 // Reduced Sigmas … … 286 295 c5.Divide(2,3); 287 296 288 CamDraw(c5,disp11,&calcam,1,2,2); 289 CamDraw(c5,disp16,&calcam,2,2,2); 297 // CamDraw(c5,disp11,calcam,1,2,1); 298 CamDraw(c5,disp11,calcam,1,2,2); 299 CamDraw(c5,disp16,calcam,2,2,1); 290 300 291 301 // F-Factor Method … … 293 303 c6.Divide(2,3); 294 304 295 CamDraw(c6,disp12, &calcam,1,2,1);296 CamDraw(c6,disp13, &calcam,2,2,2);305 CamDraw(c6,disp12,calcam,1,2,1); 306 CamDraw(c6,disp13,calcam,2,2,1); 297 307 298 308 // Blind Pixel Method … … 300 310 c7.Divide(2, 3); 301 311 302 CamDraw(c7,disp14, &calcam,1,2,9);303 CamDraw(c7,disp15, &calcam,2,2,2);312 CamDraw(c7,disp14,calcam,1,2,9); 313 CamDraw(c7,disp15,calcam,2,2,1); 304 314 305 315 } 306 316 307 void CamDraw(TCanvas &c, MHCamera &cam, MCamEvent *evt, Int_t i, Int_t j, Int_t fit)317 void CamDraw(TCanvas &c, MHCamera &cam, MCamEvent &evt, Int_t i, Int_t j, Int_t fit) 308 318 { 309 319 … … 311 321 gPad->SetBorderMode(0); 312 322 MHCamera *obj1=(MHCamera*)cam.DrawCopy("hist"); 313 obj1->AddNotify( *evt);323 obj1->AddNotify(evt); 314 324 315 325 c.cd(i+j); … … 326 336 const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst()); 327 337 const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast()); 328 const Double_t integ = obj2->Integral("width")/2.5 ;338 const Double_t integ = obj2->Integral("width")/2.5066283; 329 339 const Double_t mean = obj2->GetMean(); 330 340 const Double_t rms = obj2->GetRMS(); 331 341 const Double_t width = max-min; 332 342 333 343 if (rms == 0. || width == 0. ) 334 344 return; … … 349 359 350 360 case 1: 351 TF1 *dgaus = new TF1("dgaus","gaus(0)+gaus(3)",min,max); 361 TString dgausform = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"; 362 dgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; 363 TF1 *dgaus = new TF1("dgaus",dgausform.Data(),min,max); 352 364 dgaus->SetBit(kCanDelete); 353 dgaus->SetParNames("A 1","#mu1","#sigma1","A2","#mu2","#sigma2");354 dgaus->SetParameters(integ /width,max-width/6.,width/4.,355 integ/width ,min+width/6.,width/4.);356 dgaus->SetParLimits(0, 0,integ);357 dgaus->SetParLimits(1,min ,max);358 dgaus->SetParLimits(2,0,width/ 3.);365 dgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}","A_{2}","#mu_{2}","#sigma_{2}"); 366 dgaus->SetParameters(integ,(min+mean)/2.,width/4., 367 integ/width/2.,(max+mean)/2.,width/4.); 368 dgaus->SetParLimits(0,integ-1.5,integ+1.5); 369 dgaus->SetParLimits(1,min+(width/10.),mean); 370 dgaus->SetParLimits(2,0,width/2.); 359 371 dgaus->SetParLimits(3,0,integ); 360 dgaus->SetParLimits(4,min,max); 361 dgaus->SetParLimits(5,0,width/3.); 362 obj2->Fit("dgaus","QLR"); 372 dgaus->SetParLimits(4,mean,max-(width/10.)); 373 cout << "mean: " << mean << " maxL " << max-(width/10.) << endl; 374 375 dgaus->SetParLimits(5,0,width/2.); 376 obj2->Fit("dgaus","QLRM"); 363 377 obj2->GetFunction("dgaus")->SetLineColor(kYellow); 364 378 break; 365 379 366 380 case 2: 381 TString tgausform = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"; 382 tgausform += "+[3]/[4]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; 367 383 TF1 *tgaus = new TF1("tgaus","gaus(0)+gaus(3)+gaus(6)",min,max); 368 384 tgaus->SetBit(kCanDelete);
Note:
See TracChangeset
for help on using the changeset viewer.