Changeset 3775 for trunk/MagicSoft/Mars/macros
- Timestamp:
- 04/17/04 04:45:15 (21 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/macros/pedestalstudies.C
r3316 r3775 22 22 ! 23 23 \* ======================================================================== */ 24 25 //const TString pedfile = "/remote/home/pc2/operator/Crab20040214/20040215_16743_P_CrabOn_E.root"; 26 const TString pedfile = "../20040215_16770_P_OffCrab4_E.root"; 27 28 //const TString pedfile = "/mnt/users/mdoro/Mars/Data/20040201_14418_P_OffMrk421-1_E.root"; 29 //const TString pedfile = "/mnt/Data/rootdata/CrabNebula/2004_02_10/20040210_14607_P_CrabNebula_E.root"; 30 //const TString pedfile = "/mnt/Data/rootdata/CrabNebula/2004_01_26/20040125_10412_P_Crab-On_E.root"; 31 //const TString pedfile = "/mnt/Data/rootdata/Miscellaneous/2003_12_19/20031218_03522_P_Park_E.root"; 32 33 void pedestalstudies(TString pedname=pedfile) 24 ////////////////////////////////////////////////////////////////////////////// 25 // 26 // pedestalstudies.C 27 // 28 // macro to observe the pedestals and pedestalRMS with the number of FADC 29 // slices summed up. 30 // 31 // In order to use this macro, you have to uncomment the following 32 // line in MPedCalcPedRun (line 214): 33 // 34 // fNumHiGainSamples = runheader->GetNumSamplesHiGain() & ~1; 35 // 36 // 37 ///////////////////////////////////////////////////////////////////////////////// 38 const TString pedfile = "./20040303_20123_P_NewCalBoxTestLidOpen_E.root"; 39 40 void pedestalstudies(const TString pedname=pedfile) 34 41 { 35 42 36 gStyle->SetOptStat(1111); 37 gStyle->SetOptFit(); 38 39 MStatusDisplay *display = new MStatusDisplay; 40 display->SetUpdateTime(500); 41 display->Resize(850,700); 42 43 // 44 // Create a empty Parameter List and an empty Task List 45 // The tasklist is identified in the eventloop by its name 46 // 47 MParList plist; 48 MTaskList tlist; 49 plist.AddToList(&tlist); 50 51 // 52 // Now setup the tasks and tasklist for the pedestals: 53 // --------------------------------------------------- 54 // 55 56 MReadMarsFile read("Events", pedname); 57 read.DisableAutoScheme(); 58 59 MGeomApply geomapl; 60 MExtractSignal sigcalc; 61 62 // 63 // Set the extraction range higher: 64 // 65 //sigcalc.SetRange(1,14,1,14); 66 67 MPedCalcPedRun pedcalc; 68 69 // 70 // Additionally to calculating the pedestals, 71 // you can fill histograms and look at them 72 // 73 MFillH fill("MHPedestalCam", "MExtractedSignalCam"); 74 75 tlist.AddToList(&read); 76 tlist.AddToList(&geomapl); 77 tlist.AddToList(&sigcalc); 78 tlist.AddToList(&pedcalc); 79 tlist.AddToList(&fill); 80 81 MGeomCamMagic geomcam; 82 MPedestalCam pedcam; 83 MHPedestalCam hpedcam; 84 plist.AddToList(&geomcam); 85 plist.AddToList(&pedcam); 86 plist.AddToList(&hpedcam); 87 88 // 89 // Create and setup the eventloop 90 // 91 MEvtLoop evtloop; 92 93 evtloop.SetParList(&plist); 94 evtloop.SetDisplay(display); 95 96 // 97 // Execute first analysis 98 // 99 if (!evtloop.Eventloop()) 43 Int_t loops = 13; 44 Int_t stepsize = 2; 45 46 gStyle->SetOptStat(1111); 47 gStyle->SetOptFit(); 48 49 TArrayF *hmeandiffinn = new TArrayF(loops); 50 TArrayF *hrmsdiffinn = new TArrayF(loops); 51 TArrayF *hmeandiffout = new TArrayF(loops); 52 TArrayF *hrmsdiffout = new TArrayF(loops); 53 TArrayF *hmeaninn = new TArrayF(loops); 54 TArrayF *hmeanout = new TArrayF(loops); 55 TArrayF *hrmsinn = new TArrayF(loops); 56 TArrayF *hrmsout = new TArrayF(loops); 57 TArrayF *hmuinn = new TArrayF(loops); 58 TArrayF *hmuout = new TArrayF(loops); 59 TArrayF *hsigmainn = new TArrayF(loops); 60 TArrayF *hsigmaout = new TArrayF(loops); 61 62 TArrayF *hmeandiffinnerr = new TArrayF(loops); 63 TArrayF *hrmsdiffinnerr = new TArrayF(loops); 64 TArrayF *hmeandiffouterr = new TArrayF(loops); 65 TArrayF *hrmsdiffouterr = new TArrayF(loops); 66 TArrayF *hmeaninnerr = new TArrayF(loops); 67 TArrayF *hmeanouterr = new TArrayF(loops); 68 TArrayF *hrmsinnerr = new TArrayF(loops); 69 TArrayF *hrmsouterr = new TArrayF(loops); 70 TArrayF *hmuinnerr = new TArrayF(loops); 71 TArrayF *hmuouterr = new TArrayF(loops); 72 TArrayF *hsigmainnerr = new TArrayF(loops); 73 TArrayF *hsigmaouterr = new TArrayF(loops); 74 75 76 MStatusDisplay *display = new MStatusDisplay; 77 display->SetUpdateTime(500); 78 display->Resize(850,700); 79 80 // 81 // Create a empty Parameter List and an empty Task List 82 // The tasklist is identified in the eventloop by its name 83 // 84 MParList plist; 85 MTaskList tlist; 86 plist.AddToList(&tlist); 87 88 for (Int_t samples=2;samples<stepsize*loops+1;samples=samples+stepsize) 89 { 90 91 plist.Reset(); 92 tlist.Reset(); 93 94 // 95 // Now setup the tasks and tasklist for the pedestals: 96 // --------------------------------------------------- 97 // 98 99 MReadMarsFile read("Events", pedname); 100 read.DisableAutoScheme(); 101 102 MGeomApply geomapl; 103 // 104 // Set the extraction range higher: 105 // 106 //sigcalc.SetRange(1,14,1,14); 107 MPedCalcPedRun pedcalc; 108 pedcalc.SetNumHiGainSamples((Byte_t)samples); 109 110 MExtractSignal sigcalc; 111 sigcalc.SetRange(0,samples-1,0,samples-1); 112 113 // 114 // Additionally to calculating the pedestals, 115 // you can fill histograms and look at them 116 // 117 MFillH fill("MHPedestalCam", "MExtractedSignalCam"); 118 fill.SetNameTab(Form("%s%2d","PedCam",samples)); 119 120 tlist.AddToList(&read); 121 tlist.AddToList(&geomapl); 122 tlist.AddToList(&sigcalc); 123 tlist.AddToList(&pedcalc); 124 tlist.AddToList(&fill); 125 126 MGeomCamMagic geomcam; 127 MPedestalCam pedcam; 128 MBadPixelsCam badcam; 129 badcam.AsciiRead("badpixels.dat"); 130 131 MHPedestalCam hpedcam; 132 MCalibrationPedCam cpedcam; 133 134 plist.AddToList(&geomcam); 135 plist.AddToList(&pedcam); 136 plist.AddToList(&hpedcam); 137 plist.AddToList(&cpedcam); 138 plist.AddToList(&badcam); 139 140 // 141 // Create and setup the eventloop 142 // 143 MEvtLoop evtloop; 144 145 evtloop.SetParList(&plist); 146 evtloop.SetDisplay(display); 147 148 // 149 // Execute first analysis 150 // 151 if (!evtloop.Eventloop()) 100 152 return; 101 102 tlist.PrintStatistics(); 103 104 // 105 // Look at one specific pixel, after all the histogram manipulations: 106 // 107 // hpedcam[9].DrawClone("fourierevents"); 108 109 110 MHCamera dispped0 (geomcam, "Ped;Pedestal", "Mean per Slice"); 111 MHCamera dispped1 (geomcam, "Ped;PedestalErr", "Mean Error per Slice"); 112 MHCamera dispped2 (geomcam, "Ped;PedestalRms", "RMS per Slice"); 113 MHCamera dispped3 (geomcam, "Ped;PedestalRmsErr", "RMS Error per Slice"); 114 115 MHCamera dispped4 (geomcam, "Ped;Mean", "Fitted Mean per Slice"); 116 MHCamera dispped5 (geomcam, "Ped;MeanErr", "Fitted Error of Mean per Slice"); 117 MHCamera dispped6 (geomcam, "Ped;Sigma", "Fitted Sigma per Slice"); 118 MHCamera dispped7 (geomcam, "Ped;SigmaErr", "Fitted Error of Sigma per Slice"); 119 MHCamera dispped8 (geomcam, "Ped;Prob", "Probability of Fit"); 120 MHCamera dispped9 (geomcam, "Ped;DeltaPedestalMean", "Rel. Diff. Mean per Slice (Calc.-Fitte)"); 121 MHCamera dispped10 (geomcam, "Ped;DeltaPedestalMeanError", "Rel. Diff. Mean Error per Slice (Calc.-Fitted)"); 122 MHCamera dispped11 (geomcam, "Ped;DeltaRmsSigma", "Rel. Diff. RMS per Slice (Calc.-Fitted)"); 123 MHCamera dispped12 (geomcam, "Ped;DeltaRmsSigmaError", "Rel. Diff. RMS Error per Slice (Calc.-Fitted)"); 124 MHCamera dispped13 (geomcam, "Ped;FitOK", "Gaus Fit not OK"); 125 MHCamera dispped14 (geomcam, "Ped;FourierOK", "Fourier Analysis not OK"); 126 127 dispped0.SetCamContent( pedcam, 0); 128 dispped0.SetCamError( pedcam, 1); 129 dispped1.SetCamContent( pedcam, 1); 130 dispped2.SetCamContent( pedcam, 2); 131 dispped2.SetCamError( pedcam, 3); 132 dispped3.SetCamContent( pedcam, 3); 133 134 dispped4.SetCamContent( hpedcam, 0); 135 dispped4.SetCamError( hpedcam, 1); 136 dispped5.SetCamContent( hpedcam, 1); 137 dispped6.SetCamContent( hpedcam, 2); 138 dispped6.SetCamError( hpedcam, 3); 139 dispped7.SetCamContent( hpedcam, 3); 140 dispped8.SetCamContent( hpedcam, 4); 141 dispped9.SetCamContent( hpedcam, 5); 142 dispped9.SetCamError( hpedcam, 6); 143 dispped10.SetCamContent(hpedcam, 7); 144 dispped11.SetCamContent(hpedcam, 8); 145 dispped11.SetCamError( hpedcam, 9); 146 dispped12.SetCamContent(hpedcam, 10); 147 dispped13.SetCamContent(hpedcam, 11); 148 dispped14.SetCamContent(hpedcam, 12); 149 150 dispped0.SetYTitle("Calc. Pedestal per slice [FADC counts]"); 151 dispped1.SetYTitle("Calc. Pedestal Error per slice [FADC counts]"); 152 dispped2.SetYTitle("Calc. Pedestal RMS per slice [FADC counts]"); 153 dispped3.SetYTitle("Calc. Pedestal RMS Error per slice [FADC counts]"); 154 dispped4.SetYTitle("Fitted Mean per slice [FADC counts]"); 155 dispped5.SetYTitle("Error of Fitted Mean per slice [FADC counts]"); 156 dispped6.SetYTitle("Fitted Sigma per slice [FADC counts]"); 157 dispped7.SetYTitle("Error of Fitted Sigma per slice [FADC counts]"); 158 dispped8.SetYTitle("Fit Probability [1]"); 159 dispped9.SetYTitle("Rel. Diff. Pedestal Calc.-Fitted per slice [1]"); 160 dispped10.SetYTitle("Rel. Diff. Pedestal Error Calc.-Fitted per slice [1]"); 161 dispped11.SetYTitle("Rel. Diff. Pedestal RMS Calc.-Fitted per slice [1]"); 162 dispped12.SetYTitle("Rel. Diff. Pedestal RMS Error Calc.-Fitted per slice [1]"); 163 dispped13.SetYTitle("[1]"); 164 dispped14.SetYTitle("[1]"); 165 166 // Histogram values 167 TCanvas &b1 = display->AddTab("Ped.Calc."); 168 b1.Divide(4,3); 169 170 CamDraw(b1,dispped0,pedcam,1,4,1); 171 CamDraw(b1,dispped1,pedcam,2,4,2); 172 CamDraw(b1,dispped2,pedcam,3,4,2); 173 CamDraw(b1,dispped3,pedcam,4,4,2); 174 175 // Fitted values 176 TCanvas &b2 = display->AddTab("Ped.Fit"); 177 b2.Divide(4,3); 178 179 CamDraw(b2,dispped4,hpedcam,1,4,1); 180 CamDraw(b2,dispped5,hpedcam,2,4,2); 181 CamDraw(b2,dispped6,hpedcam,3,4,2); 182 CamDraw(b2,dispped7,hpedcam,4,4,2); 183 184 185 // Fits Probability 186 TCanvas &b3 = display->AddTab("Ped.Fit Prob."); 187 b3.Divide(1,3); 188 189 CamDraw(b3,dispped8,hpedcam,1,1,3); 190 191 // Differences 192 TCanvas &c4 = display->AddTab("Rel.Diff.Calc.-Fit"); 193 c4.Divide(4,3); 194 195 CamDraw(c4,dispped9,hpedcam,1,4,1); 196 CamDraw(c4,dispped10,hpedcam,2,4,1); 197 CamDraw(c4,dispped11,hpedcam,3,4,1); 198 CamDraw(c4,dispped12,hpedcam,4,4,1); 199 200 // Defects 201 TCanvas &c5 = display->AddTab("Defects"); 202 c5.Divide(2,2); 203 204 CamDraw(c5,dispped13,hpedcam,1,2,0); 205 CamDraw(c5,dispped14,hpedcam,2,2,0); 153 154 // 155 // Look at one specific pixel, after all the histogram manipulations: 156 // 157 /* 158 MHGausEvents &hpix = hpedcam.GetAverageHiGainArea(0); 159 hpix.DrawClone("fourierevents"); 160 161 MHGausEvents &lpix = hpedcam.GetAverageHiGainArea(1); 162 lpix.DrawClone("fourierevents"); 163 164 hpedcam[170].DrawClone("fourierevents"); 165 166 */ 167 168 MHCamera dispped0 (geomcam, "Ped;Pedestal", "Mean per Slice"); 169 MHCamera dispped2 (geomcam, "Ped;PedestalRms", "RMS per Slice"); 170 MHCamera dispped4 (geomcam, "Ped;Mean", "Fitted Mean per Slice"); 171 MHCamera dispped6 (geomcam, "Ped;Sigma", "Fitted Sigma per Slice"); 172 MHCamera dispped9 (geomcam, "Ped;DeltaPedMean", "Rel. Diff. Mean per Slice (Fit-Calc.)"); 173 MHCamera dispped11 (geomcam, "Ped;DeltaRmsSigma", "Rel. Diff. RMS per Slice (Fit-Calc.)"); 174 175 dispped0.SetCamContent( pedcam, 0); 176 dispped0.SetCamError( pedcam, 1); 177 dispped2.SetCamContent( pedcam, 2); 178 dispped2.SetCamError( pedcam, 3); 179 180 dispped4.SetCamContent( hpedcam, 0); 181 dispped4.SetCamError( hpedcam, 1); 182 dispped6.SetCamContent( hpedcam, 2); 183 dispped6.SetCamError( hpedcam, 3); 184 dispped9.SetCamContent( hpedcam, 5); 185 dispped9.SetCamError( hpedcam, 6); 186 dispped11.SetCamContent(hpedcam, 8); 187 dispped11.SetCamError( hpedcam, 9); 188 189 dispped0.SetYTitle("Calc. Pedestal per slice [FADC counts]"); 190 dispped2.SetYTitle("Calc. Pedestal RMS per slice [FADC counts]"); 191 dispped4.SetYTitle("Fitted Mean per slice [FADC counts]"); 192 dispped6.SetYTitle("Fitted Sigma per slice [FADC counts]"); 193 dispped9.SetYTitle("Rel. Diff. Pedestal per slice Fit-Calc [1]"); 194 dispped11.SetYTitle("Rel. Diff. Pedestal RMS per slice Fit-Calc [1]"); 195 196 197 // Histogram values 198 TCanvas &b1 = display->AddTab(Form("%s%d","MeanRMS",samples)); 199 b1.Divide(4,3); 200 201 CamDraw(b1,dispped0,1,4,*hmeaninn,*hmeanout,*hmeaninnerr,*hmeanouterr,samples,stepsize); 202 CamDraw(b1,dispped2,2,4,*hrmsinn,*hrmsout,*hrmsinnerr,*hrmsouterr,samples,stepsize); 203 CamDraw(b1,dispped4,3,4,*hmuinn,*hmuout,*hmuinnerr,*hmuouterr,samples,stepsize); 204 CamDraw(b1,dispped6,4,4,*hsigmainn,*hsigmaout,*hsigmainnerr,*hsigmaouterr,samples,stepsize); 205 206 display->SaveAsGIF(3*((samples-1)/stepsize)+2,Form("%s%d","MeanRmsSamples",samples)); 207 208 // Differences 209 TCanvas &c4 = display->AddTab(Form("%s%d","RelDiff",samples)); 210 c4.Divide(2,3); 211 212 CamDraw(c4,dispped9,1,2,*hmeandiffinn,*hmeandiffout,*hmeandiffinnerr,*hmeandiffouterr,samples,stepsize); 213 CamDraw(c4,dispped11,2,2,*hrmsdiffinn,*hrmsdiffout,*hrmsdiffinnerr,*hrmsdiffouterr,samples,stepsize); 214 215 display->SaveAsGIF(3*((samples-1)/stepsize)+3,Form("%s%d","RelDiffSamples",samples)); 216 217 } 218 219 TF1 *logg = new TF1("logg","[1]+TMath::Log(x-[0])",1.,30.,2); 220 logg->SetParameters(1.,3.5); 221 logg->SetParLimits(0,-1.,3.); 222 logg->SetParLimits(1,-1.,7.); 223 logg->SetLineColor(kRed); 224 225 TCanvas *canvas = new TCanvas("PedstudInner","Pedestal Studies Inner Pixels",600,900); 226 canvas->Divide(2,3); 227 canvas->cd(1); 228 229 TGraphErrors *gmeaninn = new TGraphErrors(hmeaninn->GetSize(), 230 CreateXaxis(hmeaninn->GetSize(),stepsize),hmeaninn->GetArray(), 231 CreateXaxisErr(hmeaninnerr->GetSize(),stepsize),hmeaninnerr->GetArray()); 232 gmeaninn->Draw("A*"); 233 gmeaninn->SetTitle("Calculated Mean per Slice Inner Pixels"); 234 gmeaninn->GetXaxis()->SetTitle("Nr. added FADC slices"); 235 gmeaninn->GetYaxis()->SetTitle("Calculated Mean per slice"); 236 gmeaninn->Fit("pol0"); 237 gmeaninn->GetFunction("pol0")->SetLineColor(kGreen); 238 // gmeaninn->Fit(logg); 239 240 canvas->cd(2); 241 242 TGraphErrors *gmuinn = new TGraphErrors(hmuinn->GetSize(), 243 CreateXaxis(hmuinn->GetSize(),stepsize),hmuinn->GetArray(), 244 CreateXaxisErr(hmuinnerr->GetSize(),stepsize),hmuinnerr->GetArray()); 245 gmuinn->Draw("A*"); 246 gmuinn->SetTitle("Fitted Mean per Slice Inner Pixels"); 247 gmuinn->GetXaxis()->SetTitle("Nr. added FADC slices"); 248 gmuinn->GetYaxis()->SetTitle("Fitted Mean per Slice"); 249 gmuinn->Fit("pol0"); 250 gmuinn->GetFunction("pol0")->SetLineColor(kGreen); 251 //gmuinn->Fit(logg); 252 253 254 canvas->cd(3); 255 256 TGraphErrors *grmsinn = new TGraphErrors(hrmsinn->GetSize(), 257 CreateXaxis(hrmsinn->GetSize(),stepsize),hrmsinn->GetArray(), 258 CreateXaxisErr(hrmsinnerr->GetSize(),stepsize),hrmsinnerr->GetArray()); 259 grmsinn->Draw("A*"); 260 grmsinn->SetTitle("Calculated Rms per Slice Inner Pixels"); 261 grmsinn->GetXaxis()->SetTitle("Nr. added FADC slices"); 262 grmsinn->GetYaxis()->SetTitle("Calculated Rms per Slice"); 263 //grmsinn->Fit("pol2"); 264 //grmsinn->GetFunction("pol2")->SetLineColor(kRed); 265 grmsinn->Fit(logg); 266 267 canvas->cd(4); 268 269 TGraphErrors *gsigmainn = new TGraphErrors(hsigmainn->GetSize(), 270 CreateXaxis(hsigmainn->GetSize(),stepsize),hsigmainn->GetArray(), 271 CreateXaxisErr(hsigmainnerr->GetSize(),stepsize),hsigmainnerr->GetArray()); 272 gsigmainn->Draw("A*"); 273 gsigmainn->SetTitle("Fitted Sigma per Slice Inner Pixels"); 274 gsigmainn->GetXaxis()->SetTitle("Nr. added FADC slices"); 275 gsigmainn->GetYaxis()->SetTitle("Fitted Sigma per Slice"); 276 // gsigmainn->Fit("pol2"); 277 // gsigmainn->GetFunction("pol2")->SetLineColor(kRed); 278 gsigmainn->Fit(logg); 279 280 canvas->cd(5); 281 282 TGraphErrors *gmeandiffinn = new TGraphErrors(hmeandiffinn->GetSize(), 283 CreateXaxis(hmeandiffinn->GetSize(),stepsize),hmeandiffinn->GetArray(), 284 CreateXaxisErr(hmeandiffinnerr->GetSize(),stepsize),hmeandiffinnerr->GetArray()); 285 gmeandiffinn->Draw("A*"); 286 gmeandiffinn->SetTitle("Rel. Difference Mean per Slice Inner Pixels"); 287 gmeandiffinn->GetXaxis()->SetTitle("Nr. added FADC slices"); 288 gmeandiffinn->GetYaxis()->SetTitle("Rel. Difference Mean per Slice"); 289 //gmeandiffinn->Fit("pol2"); 290 //gmeandiffinn->GetFunction("pol2")->SetLineColor(kBlue); 291 gmeandiffinn->Fit(logg); 292 293 294 canvas->cd(6); 295 296 TGraphErrors *grmsdiffinn = new TGraphErrors(hrmsdiffinn->GetSize(), 297 CreateXaxis(hrmsdiffinn->GetSize(),stepsize),hrmsdiffinn->GetArray(), 298 CreateXaxisErr(hrmsdiffinnerr->GetSize(),stepsize),hrmsdiffinnerr->GetArray()); 299 grmsdiffinn->Draw("A*"); 300 grmsdiffinn->SetTitle("Rel. Difference Sigma per Slice-RMS Inner Pixels"); 301 grmsdiffinn->GetXaxis()->SetTitle("Nr. added FADC slices"); 302 grmsdiffinn->GetYaxis()->SetTitle("Rel. Difference Sigma per Slice-RMS"); 303 //grmsdiffinn->Fit("pol2"); 304 //grmsdiffinn->GetFunction("pol2")->SetLineColor(kBlue); 305 grmsdiffinn->Fit(logg); 306 307 308 TCanvas *canvas2 = new TCanvas("PedstudOut","Pedestal Studies Outer Pixels",600,900); 309 canvas2->Divide(2,3); 310 canvas2->cd(1); 311 312 313 canvas2->cd(1); 314 315 TGraphErrors *gmeanout = new TGraphErrors(hmeanout->GetSize(), 316 CreateXaxis(hmeanout->GetSize(),stepsize),hmeanout->GetArray(), 317 CreateXaxisErr(hmeanouterr->GetSize(),stepsize),hmeanouterr->GetArray()); 318 gmeanout->Draw("A*"); 319 gmeanout->SetTitle("Calculated Mean per Slice Outer Pixels"); 320 gmeanout->GetXaxis()->SetTitle("Nr. added FADC slices"); 321 gmeanout->GetYaxis()->SetTitle("Calculated Mean per Slice"); 322 gmeanout->Fit("pol0"); 323 gmeanout->GetFunction("pol0")->SetLineColor(kGreen); 324 //gmeanout->Fit(logg); 325 326 canvas2->cd(2); 327 328 TGraphErrors *gmuout = new TGraphErrors(hmuout->GetSize(), 329 CreateXaxis(hmuout->GetSize(),stepsize),hmuout->GetArray(), 330 CreateXaxisErr(hmuouterr->GetSize(),stepsize),hmuouterr->GetArray()); 331 gmuout->Draw("A*"); 332 gmuout->SetTitle("Fitted Mean per Slice Outer Pixels"); 333 gmuout->GetXaxis()->SetTitle("Nr. added FADC slices"); 334 gmuout->GetYaxis()->SetTitle("Fitted Mean per Slice"); 335 gmuout->Fit("pol0"); 336 gmuout->GetFunction("pol0")->SetLineColor(kGreen); 337 //gmuout->Fit(logg); 338 339 canvas2->cd(3); 340 341 TGraphErrors *grmsout = new TGraphErrors(hrmsout->GetSize(), 342 CreateXaxis(hrmsout->GetSize(),stepsize),hrmsout->GetArray(), 343 CreateXaxisErr(hrmsouterr->GetSize(),stepsize),hrmsouterr->GetArray()); 344 grmsout->Draw("A*"); 345 grmsout->SetTitle("Calculated Rms per Slice Outer Pixels"); 346 grmsout->GetXaxis()->SetTitle("Nr. added FADC slices"); 347 grmsout->GetYaxis()->SetTitle("Calculated Rms per Slice"); 348 //grmsout->Fit("pol2"); 349 //grmsout->GetFunction("pol2")->SetLineColor(kRed); 350 grmsout->Fit(logg); 351 352 canvas2->cd(4); 353 354 TGraphErrors *gsigmaout = new TGraphErrors(hsigmaout->GetSize(), 355 CreateXaxis(hsigmaout->GetSize(),stepsize),hsigmaout->GetArray(), 356 CreateXaxisErr(hsigmaouterr->GetSize(),stepsize),hsigmaouterr->GetArray()); 357 gsigmaout->Draw("A*"); 358 gsigmaout->SetTitle("Fitted Sigma per Slice Outer Pixels"); 359 gsigmaout->GetXaxis()->SetTitle("Nr. added FADC slices"); 360 gsigmaout->GetYaxis()->SetTitle("Fitted Sigma per Slice"); 361 //gsigmaout->Fit("pol2"); 362 //gsigmaout->GetFunction("pol2")->SetLineColor(kRed); 363 gsigmaout->Fit(logg); 364 365 366 canvas2->cd(5); 367 368 TGraphErrors *gmeandiffout = new TGraphErrors(hmeandiffout->GetSize(), 369 CreateXaxis(hmeandiffout->GetSize(),stepsize),hmeandiffout->GetArray(), 370 CreateXaxisErr(hmeandiffouterr->GetSize(),stepsize),hmeandiffouterr->GetArray()); 371 gmeandiffout->Draw("A*"); 372 gmeandiffout->SetTitle("Rel. Difference Mean per Slice Outer Pixels"); 373 gmeandiffout->GetXaxis()->SetTitle("Nr. added FADC slices"); 374 gmeandiffout->GetYaxis()->SetTitle("Rel. Difference Mean per Slice"); 375 //gmeandiffout->Fit("pol2"); 376 //gmeandiffout->GetFunction("pol2")->SetLineColor(kBlue); 377 gmeandiffout->Fit(logg); 378 379 canvas2->cd(6); 380 381 TGraphErrors *grmsdiffout = new TGraphErrors(hrmsdiffout->GetSize(), 382 CreateXaxis(hrmsdiffout->GetSize(),stepsize),hrmsdiffout->GetArray(), 383 CreateXaxisErr(hrmsdiffouterr->GetSize(),stepsize),hrmsdiffouterr->GetArray()); 384 grmsdiffout->Draw("A*"); 385 grmsdiffout->SetTitle("Rel. Difference Sigma per Slice-RMS Outer Pixels"); 386 grmsdiffout->GetXaxis()->SetTitle("Nr. added FADC slices"); 387 grmsdiffout->GetYaxis()->SetTitle("Rel. Difference Sigma per Slice-RMS"); 388 //grmsdiffout->Fit("pol2"); 389 //grmsdiffout->GetFunction("pol2")->SetLineColor(kBlue); 390 grmsdiffout->Fit(logg); 391 206 392 207 393 } 208 394 209 void CamDraw(TCanvas &c, MHCamera &cam, MCamEvent &evt, Int_t i, Int_t j, Int_t fit) 395 396 void CamDraw(TCanvas &c, MHCamera &cam, Int_t i, Int_t j, TArrayF &a1, TArrayF &a2, 397 TArrayF &a1err, TArrayF &a2err, Int_t samp, Int_t stepsize) 210 398 { 211 399 212 400 c.cd(i); 213 gPad->SetBorderMode(0);214 401 MHCamera *obj1=(MHCamera*)cam.DrawCopy("hist"); 215 // obj1->AddNotify(evt);402 obj1->SetDirectory(NULL); 216 403 217 404 c.cd(i+j); 218 gPad->SetBorderMode(0);219 405 obj1->Draw(); 220 406 ((MHCamera*)obj1)->SetPrettyPalette(); 221 407 222 if (fit != 0) 408 c.cd(i+2*j); 409 TH1D *obj2 = (TH1D*)obj1->Projection(); 410 obj2->SetDirectory(NULL); 411 412 // obj2->Sumw2(); 413 obj2->Draw(); 414 obj2->SetBit(kCanDelete); 415 416 const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst()); 417 const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast()); 418 const Double_t integ = obj2->Integral("width")/2.5066283; 419 const Double_t mean = obj2->GetMean(); 420 const Double_t rms = obj2->GetRMS(); 421 const Double_t width = max-min; 422 423 if (rms == 0. || width == 0. ) 424 return; 425 426 TArrayI s0(6); 427 s0[0] = 6; 428 s0[1] = 1; 429 s0[2] = 2; 430 s0[3] = 3; 431 s0[4] = 4; 432 s0[5] = 5; 433 434 TArrayI inner(1); 435 inner[0] = 0; 436 437 TArrayI outer(1); 438 outer[0] = 1; 439 440 // Just to get the right (maximum) binning 441 TH1D *half[2]; 442 half[0] = obj1->ProjectionS(s0, inner, "Inner"); 443 half[1] = obj1->ProjectionS(s0, outer, "Outer"); 444 445 half[0]->SetDirectory(NULL); 446 half[1]->SetDirectory(NULL); 447 448 for (int i=0; i<2; i++) 223 449 { 224 c.cd(i+2*j); 225 gPad->SetBorderMode(0); 226 TH1D *obj2 = (TH1D*)obj1->Projection(); 227 228 // obj2->Sumw2(); 229 obj2->Draw(); 230 obj2->SetBit(kCanDelete); 231 232 const Double_t min = obj2->GetBinCenter(obj2->GetXaxis()->GetFirst()); 233 const Double_t max = obj2->GetBinCenter(obj2->GetXaxis()->GetLast()); 234 const Double_t integ = obj2->Integral("width")/2.5066283; 235 const Double_t mean = obj2->GetMean(); 236 const Double_t rms = obj2->GetRMS(); 237 const Double_t width = max-min; 238 239 if (rms == 0. || width == 0. ) 240 return; 241 242 switch (fit) 450 half[i]->SetLineColor(kRed+i); 451 half[i]->SetDirectory(0); 452 half[i]->SetBit(kCanDelete); 453 half[i]->Draw("same"); 454 half[i]->Fit("gaus","Q+"); 455 456 if (i==0) 243 457 { 244 case 1: 245 TF1 *sgaus = new TF1("sgaus","gaus(0)",min,max); 246 sgaus->SetBit(kCanDelete); 247 sgaus->SetParNames("Area","#mu","#sigma"); 248 sgaus->SetParameters(integ/rms,mean,rms); 249 sgaus->SetParLimits(0,0.,integ); 250 sgaus->SetParLimits(1,min,max); 251 sgaus->SetParLimits(2,0,width/1.5); 252 obj2->Fit("sgaus","QLR"); 253 obj2->GetFunction("sgaus")->SetLineColor(kYellow); 254 break; 255 256 case 2: 257 TString dgausform = "([0]-[3])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"; 258 dgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; 259 TF1 *dgaus = new TF1("dgaus",dgausform.Data(),min,max); 260 dgaus->SetBit(kCanDelete); 261 dgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}","A_{2}","#mu_{2}","#sigma_{2}"); 262 dgaus->SetParameters(integ,(min+mean)/2.,width/4., 263 integ/width/2.,(max+mean)/2.,width/4.); 264 // The left-sided Gauss 265 dgaus->SetParLimits(0,integ-1.5,integ+1.5); 266 dgaus->SetParLimits(1,min+(width/10.),mean); 267 dgaus->SetParLimits(2,0,width/2.); 268 // The right-sided Gauss 269 dgaus->SetParLimits(3,0,integ); 270 dgaus->SetParLimits(4,mean,max-(width/10.)); 271 dgaus->SetParLimits(5,0,width/2.); 272 obj2->Fit("dgaus","QLRM"); 273 obj2->GetFunction("dgaus")->SetLineColor(kYellow); 274 break; 275 276 case 3: 277 TString tgausform = "([0]-[3]-[6])/[2]*exp(-0.5*(x-[1])*(x-[1])/[2]/[2])"; 278 tgausform += "+[3]/[5]*exp(-0.5*(x-[4])*(x-[4])/[5]/[5])"; 279 tgausform += "+[6]/[8]*exp(-0.5*(x-[7])*(x-[7])/[8]/[8])"; 280 TF1 *tgaus = new TF1("tgaus",tgausform.Data(),min,max); 281 tgaus->SetBit(kCanDelete); 282 tgaus->SetParNames("A_{tot}","#mu_{1}","#sigma_{1}", 283 "A_{2}","#mu_{2}","#sigma_{2}", 284 "A_{3}","#mu_{3}","#sigma_{3}"); 285 tgaus->SetParameters(integ,(min+mean)/2,width/4., 286 integ/width/3.,(max+mean)/2.,width/4., 287 integ/width/3.,mean,width/2.); 288 // The left-sided Gauss 289 tgaus->SetParLimits(0,integ-1.5,integ+1.5); 290 tgaus->SetParLimits(1,min+(width/10.),mean); 291 tgaus->SetParLimits(2,width/15.,width/2.); 292 // The right-sided Gauss 293 tgaus->SetParLimits(3,0.,integ); 294 tgaus->SetParLimits(4,mean,max-(width/10.)); 295 tgaus->SetParLimits(5,width/15.,width/2.); 296 // The Gauss describing the outliers 297 tgaus->SetParLimits(6,0.,integ); 298 tgaus->SetParLimits(7,min,max); 299 tgaus->SetParLimits(8,width/4.,width/1.5); 300 obj2->Fit("tgaus","QLRM"); 301 obj2->GetFunction("tgaus")->SetLineColor(kYellow); 302 break; 303 case 4: 304 obj2->Fit("pol0","Q"); 305 obj2->GetFunction("pol0")->SetLineColor(kYellow); 306 break; 307 case 9: 308 break; 309 default: 310 obj2->Fit("gaus","Q"); 311 obj2->GetFunction("gaus")->SetLineColor(kYellow); 312 break; 458 a1[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParameter(1); 459 a1err[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParError(1); 460 if (a1err[(samp-1)/stepsize] > 3.) 461 a1err[(samp-1)/stepsize] = 1.; 313 462 } 314 315 gPad->Modified(); 316 gPad->Update(); 317 463 if (i==1) 464 { 465 a2[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParameter(1); 466 a2err[(samp-1)/stepsize] = half[i]->GetFunction("gaus")->GetParError(1); 467 if (a2err[(samp-1)/stepsize] > 3.) 468 a2err[(samp-1)/stepsize] = 1.; 469 } 318 470 } 471 472 319 473 } 474 475 // ----------------------------------------------------------------------------- 476 // 477 // Create the x-axis for the event graph 478 // 479 Float_t *CreateXaxis(Int_t n, Int_t step) 480 { 481 482 Float_t *xaxis = new Float_t[n]; 483 484 for (Int_t i=0;i<n;i++) 485 xaxis[i] = 2. + step*i; 486 487 return xaxis; 488 489 } 490 491 // ----------------------------------------------------------------------------- 492 // 493 // Create the x-axis for the event graph 494 // 495 Float_t *CreateXaxisErr(Int_t n, Int_t step) 496 { 497 498 Float_t *xaxis = new Float_t[n]; 499 500 for (Int_t i=0;i<n;i++) 501 xaxis[i] = step/2.; 502 503 return xaxis; 504 505 }
Note:
See TracChangeset
for help on using the changeset viewer.