source: fact/tools/rootmacros/FPulseTemplate.C@ 12918

Last change on this file since 12918 was 12742, checked in by Jens Buss, 14 years ago
Macro to produce pulse templates ala tpeak for a given pulsehight
File size: 26.2 KB
Line 
1/********************** FPulseTemplate ***********************
2* read FACT raw data
3* remove spikes
4* calculate baseline
5* find peaks (CFD and normal discriminator)
6* compute pulse height and pulse integral spektrum of the peaks
7+ search for Peaks in data
8+ put peaks into different histos depending und Amplitude
9+ draw histos for single, double, tripple, ....above 100mV Photonpulses
10+ draw Parameterdevelopment depending on photon quantity
11+ form a template (shape depending on number of photons)
12 so it can be used for detection of other peaks
13*************************************************************/
14
15//-----------------------------------------------------------------------------
16// root libraries
17//-----------------------------------------------------------------------------
18
19#include <TROOT.h>
20#include <TCanvas.h>
21#include <TProfile.h>
22#include <TTimer.h>
23#include <TH1F.h>
24#include <TH2F.h>
25#include <Getline.h>
26#include <TLine.h>
27#include <TBox.h>
28#include <TMath.h>
29#include <TFile.h>
30#include <TStyle.h>
31#include <TString.h>
32
33#include <stdio.h>
34#include <stdint.h>
35#include <cstdio>
36
37#define NPIX 1440
38#define NCELL 1024
39#define FAD_MAX_SAMPLES 1024
40
41#define HAVE_ZLIB
42
43//-----------------------------------------------------------------------------
44// rootmacros
45//-----------------------------------------------------------------------------
46
47#include "fits.h"
48
49#include "openFits.h"
50#include "openFits.c"
51
52#include "discriminator.h"
53#include "discriminator.C"
54#include "zerosearch.h"
55#include "zerosearch.C"
56#include "factfir.C"
57
58#include "DrsCalibration.C"
59#include "DrsCalibration.h"
60
61#include "SpikeRemoval.h"
62#include "SpikeRemoval.C"
63
64//-----------------------------------------------------------------------------
65// Decleration of global variables
66//-----------------------------------------------------------------------------
67
68bool breakout=false;
69
70int NEvents;
71vector<int16_t> AllPixelDataVector;
72vector<int16_t> StartCellVector;
73unsigned int CurrentEventID;
74size_t PXLxROI;
75UInt_t RegionOfInterest;
76UInt_t NumberOfPixels;
77
78size_t TriggerOffsetROI, RC;
79size_t drs_n;
80vector<float> drs_basemean;
81vector<float> drs_gainmean;
82vector<float> drs_triggeroffsetmean;
83
84vector<float> Ameas(FAD_MAX_SAMPLES); // copy of the data (measured amplitude
85vector<float> Vcorr(FAD_MAX_SAMPLES); // corrected Values
86vector<float> Vslide(FAD_MAX_SAMPLES); // sliding average result
87vector<float> Vcfd(FAD_MAX_SAMPLES); // CDF result
88vector<float> Vcfd2(FAD_MAX_SAMPLES); // CDF result + 2nd sliding average
89
90typedef struct{
91 float maxAmpl;
92 float countOfMax;
93 } OverlayMaximum;
94
95vector<OverlayMaximum> SingleMaximum;
96// histograms
97const int NumberOfDebugHistoTypes = 7;
98
99const unsigned int
100 Ameas_ = 0,
101 N1mean_ = 1,
102 Vcorr_ = 2,
103 Vtest_ = 3,
104 Vslide_ = 4,
105 Vcfd_ = 5,
106 Vcfd2_ = 6;
107
108//-----------------------------------------------------------------------------
109// Initialisation of histograms
110//-----------------------------------------------------------------------------
111
112TH1F* h;
113TH1F *debugHistos;
114TH2F *hOverlayTemp = NULL;
115TH2F *hSinglePeakOverlay;//histogrammm for overlay of detected Peaks
116TH2F *hDoublePeakOverlay;
117TH2F *hTripplePeakOverlay;
118TH2F *hLargePeakOverlay;
119
120TH1F *hMaximumTemp = NULL;
121TH1F *hSinglePeakMaximum;
122TH1F *hDoublePeakMaximum;
123TH1F *hTripplePeakMaximum;
124TH1F *hLargePeakMaximum;
125
126int hPeakOverlayXaxisLeft,hPeakOverlayXaxisRight;
127
128TObjArray hList;
129TObjArray hListBaseline;
130
131//-----------------------------------------------------------------------------
132// Functions
133//-----------------------------------------------------------------------------
134
135void BookHistos( );
136void SaveHistograms( const char * );
137
138//-----------------------------------------------------------------------------
139//-----------------------------------------------------------------------------
140// MAIN - Funtion
141//-----------------------------------------------------------------------------
142//-----------------------------------------------------------------------------
143int FPulseTemplate(
144 char *datafilename = "../data/2011/11/09/20111109_006.fits.gz",
145 const char *drsfilename = "../data/2011/11/09/20111109_003.drs.fits.gz",
146 const char *OutRootFileName = "../analysis/20111109_006.fits.gz.peaktypes.root",
147 int firstevent = 0,
148 int nevents = -1,
149 int firstpixel = 0,
150 int npixel = -1,
151 bool spikeDebug = false,
152 int avg1 = 14,
153 int avg2 = 8,
154 int OverlayWindowLeft = 50,
155 int OverlayWindowRight = 250,
156 int verbosityLevel = 1, // different verbosity levels can be implemented here
157 bool ProduceGraphic = true
158 )
159
160{
161//-----------------------------------------------------------------------------
162// histogramm Settings
163//-----------------------------------------------------------------------------
164
165 hPeakOverlayXaxisLeft = OverlayWindowLeft;
166 hPeakOverlayXaxisRight = OverlayWindowRight;
167
168 gStyle->SetPalette(1,0);
169 gROOT->SetStyle("Plain");
170
171// Create (pointer to) Canvases, which are used in every run,
172// also in 'non-debug' runs
173 // Canvases only need if spike Debug, but I want to deklare
174 // the pointers anyway ...
175 TCanvas *cFiltered = NULL;
176 if (spikeDebug){
177 cFiltered = new TCanvas("cFiltered","filtered DRS Waveforms", 1,310,400,300);
178 cFiltered->Divide(1, 3);
179 }
180 // Canvases to show the peak template
181 TCanvas *cPulsetypes = NULL;
182 cPulsetypes = new TCanvas("cPulsetypes", "Overlay of detected Peaks", 1, 1, 1400, 700);
183 cPulsetypes->Divide(4,2);
184
185//-----------------------------------------------------------------------------
186// Filter-Settings
187//-----------------------------------------------------------------------------
188// CFD filter settings
189 int k_cfd = 10;
190 vector<double> a_cfd(k_cfd, 0);
191 double b_cfd = 1.;
192 a_cfd[0]=-0.75;
193 a_cfd[k_cfd-1]=1.;
194
195//-----------------------------------------------------------------------------
196// prepare datafile
197//-----------------------------------------------------------------------------
198
199// Open the data file
200
201 fits * datafile;
202 // Opens the raw data file and 'binds' the variables given as
203 // Parameters to the data file. So they are filled with
204 // raw data as soon as datafile->GetRow(int) is called.
205 NEvents = openDataFits( datafilename, &datafile,
206 AllPixelDataVector, StartCellVector, CurrentEventID,
207 RegionOfInterest, NumberOfPixels, PXLxROI, verbosityLevel);
208 if (NEvents == 0){
209 cout << "return code of OpenDataFile:" << datafilename<< endl;
210 cout << "is zero -> aborting." << endl;
211 return 1;
212 }
213
214 if (verbosityLevel > 0)
215 cout <<"number of events in file: "<< NEvents << "\t";
216 if ( nevents == -1 || nevents > NEvents ) nevents = NEvents; // -1 means all!
217 if (verbosityLevel > 0)
218 cout <<"of, which "<< nevents << "will be processed"<< endl;
219
220 if (verbosityLevel > 0)
221 cout <<"Total # of Pixel: "<< NumberOfPixels << "\t";
222 if ( npixel == -1 || npixel > (int)NumberOfPixels ) npixel = (int)NumberOfPixels; // -1 means all!
223 if (verbosityLevel > 0)
224 cout <<"of, which "<< npixel << "will be processed"<< endl;
225//Get the DRS calibration
226 RC = openCalibFits( drsfilename, drs_basemean, drs_gainmean, drs_triggeroffsetmean, TriggerOffsetROI);
227 if (RC == 0){
228 cout << "return code of openCalibFits:" << drsfilename << endl;
229 cout << "is zero -> aborting." << endl;
230 return 1;
231 }
232// Book the histograms
233 BookHistos( );
234
235//-----------------------------------------------------------------------------
236// Loops over Every Event and Pixel
237//-----------------------------------------------------------------------------
238 for ( int ev = firstevent; ev < firstevent + nevents; ev++) {
239 // Get an Event --> consists of 1440 Pixel ...erm....data
240 datafile->GetRow( ev );
241
242//-------------------------------------
243// Loop over every Pixel of Event
244//-------------------------------------
245
246 for ( int pix = firstpixel; pix < firstpixel + npixel; pix++ ){
247 if (verbosityLevel > 0){
248 if (pix == firstpixel){
249 cout << "Processing Event: " << CurrentEventID << "/" << nevents << endl;
250 }
251 }
252
253 applyDrsCalibration( Ameas,pix,12,12,
254 drs_basemean, drs_gainmean, drs_triggeroffsetmean,
255 RegionOfInterest, AllPixelDataVector, StartCellVector);
256
257 // finds spikes in the raw data, and interpolates the value
258 // spikes are: 1 or 2 slice wide, positive non physical artifacts
259 removeSpikes (Ameas, Vcorr);
260
261 // filter Vcorr with sliding average using FIR filter function
262 sliding_avg(Vcorr, Vslide, avg1);
263 // filter Vslide with CFD using FIR filter function
264 factfir(b_cfd , a_cfd, k_cfd, Vslide, Vcfd);
265 // filter Vcfd with sliding average using FIR filter function
266 sliding_avg(Vcfd, Vcfd2, avg2);
267
268 // peaks in Ameas[] are found by searching for zero crossings
269 // in Vcfd2
270 // first Argument 1 means ... *rising* edge
271 // second Argument 1 means ... search with stepsize 1 ... 10 is okay as well
272 vector<Region> * zXings = zerosearch( Vcfd2 , 1 , 8);
273 // zXings means "zero cross ings"
274 EnlargeRegion(*zXings, 10, 10);
275 findAbsMaxInRegions(*zXings, Vslide);
276 removeMaximaBelow( *zXings, 3.0);
277 removeRegionWithMaxOnEdge( *zXings, 2);
278 removeRegionOnFallingEdge( *zXings, 100);
279
280//-----------------------------------------------------------------------------
281// Fill Overlay Histos
282//-----------------------------------------------------------------------------
283 vector<Region>::iterator it;
284 for (it = zXings->begin() ; it < zXings->end() ; it++){
285 if (it->maxPos < 12 || it->maxPos > RegionOfInterest-12)
286 continue;
287 //domstest->Fill(it->maxPos);
288 int Left = it->maxPos - OverlayWindowLeft;
289 int Right = it->maxPos + OverlayWindowRight;
290 if (Left < 0)
291 Left =0;
292 if (Right > (int)Vcorr.size() )
293 Right=Vcorr.size() ;
294 if (Vslide[it->maxPos] >= 5 && Vslide[it->maxPos] < 13) hOverlayTemp = &*hSinglePeakOverlay;
295 else if (Vslide[it->maxPos] >= 13 && Vslide[it->maxPos] < 23) hOverlayTemp = hDoublePeakOverlay;
296 else if (Vslide[it->maxPos] >= 23 && Vslide[it->maxPos] < 33) hOverlayTemp = hTripplePeakOverlay;
297 else if (Vslide[it->maxPos] >= 33) hOverlayTemp = hLargePeakOverlay;
298
299 for ( int pos = Left; pos < Right; pos++){
300// if (Vslide[it->maxPos] >= 5 && Vslide[it->maxPos] < 15) hSinglePeakOverlay->Fill( pos - (it->maxPos), Vslide[pos]);
301// else if (Vslide[it->maxPos] >= 15 && Vslide[it->maxPos] < 25) hDoublePeakOverlay->Fill( pos - (it->maxPos), Vslide[pos]);
302// else if (Vslide[it->maxPos] >= 25 && Vslide[it->maxPos] < 35) hTripplePeakOverlay->Fill( pos - (it->maxPos), Vslide[pos]);
303// else if (Vslide[it->maxPos] >= 35) hOverlayTemp = hLargePeakOverlay;
304 hOverlayTemp->Fill( pos - (it->maxPos), Vslide[pos]) ;
305 }
306 }
307
308//-----------------------------------------------------------------------------
309// Spike Debug
310//-----------------------------------------------------------------------------
311 if ( spikeDebug ){
312 // TODO do this correct. The vectors should be the rigt ones... this is just luck
313 debugHistos[Ameas_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
314 debugHistos[Vcorr_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
315 debugHistos[Vslide_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
316 debugHistos[Vcfd_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
317 debugHistos[Vcfd2_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
318
319 for ( unsigned int sl = 0; sl < RegionOfInterest; sl++){
320 debugHistos[Ameas_].SetBinContent(sl, Ameas[sl]);
321 debugHistos[Vcorr_].SetBinContent(sl, Vcorr[sl]);
322 debugHistos[Vslide_].SetBinContent( sl, Vslide[sl] );
323 debugHistos[Vcfd_].SetBinContent( sl, Vcfd[sl] );
324 debugHistos[Vcfd2_].SetBinContent( sl, Vcfd2[sl] );
325 }
326
327
328 cFiltered->cd(1);
329 gPad->SetGrid();
330 debugHistos[Vcorr_].Draw();
331
332 cFiltered->cd(2);
333 gPad->SetGrid();
334 debugHistos[Vslide_].Draw();
335
336 TBox *OneBox;
337 vector<TBox*> MyBoxes;
338 for (unsigned int i=0; i<zXings->size(); i++){
339 OneBox = new TBox(
340 zXings->at(i).maxPos -10 ,
341 zXings->at(i).maxVal -0.5,
342 zXings->at(i).maxPos +10 ,
343 zXings->at(i).maxVal +0.5);
344 OneBox->SetLineColor(kBlue);
345 OneBox->SetLineWidth(1);
346 OneBox->SetFillStyle(0);
347 OneBox->SetFillColor(kRed);
348 MyBoxes.push_back(OneBox);
349 OneBox->Draw();
350 }
351
352 cFiltered->cd(3);
353 gPad->SetGrid();
354 debugHistos[Vcfd2_].Draw();
355 TLine *zeroline = new TLine(0, 0, 1024, 0);
356 zeroline->SetLineColor(kBlue);
357 zeroline->Draw();
358
359 cFiltered->Update();
360
361 //OVERLAY PEAKS
362 cPulsetypes->cd(1);
363 gPad->SetGrid();
364 hSinglePeakOverlay->Draw("COLZ");
365
366 cPulsetypes->cd(2);
367 gPad->SetGrid();
368 hDoublePeakOverlay->Draw("COLZ");
369
370 cPulsetypes->cd(3);
371 gPad->SetGrid();
372 hTripplePeakOverlay->Draw("COLZ");
373
374 cPulsetypes->cd(4);
375 gPad->SetGrid();
376 hLargePeakOverlay->Draw("COLZ");
377
378 cPulsetypes->Modified();
379 cPulsetypes->Update();
380
381 //Process gui events asynchronously during input
382 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
383 timer.TurnOn();
384 TString input = Getline("Type 'q' to exit, <return> to go on: ");
385 timer.TurnOff();
386 if (input=="q\n") {
387 breakout=true;
388 }
389
390 //TODO!!!!!!!!!
391 // do some Garbage collection ...
392 // all the Objects on the heap should be deleted here.
393
394 }// end of if(spikeDebug)
395
396 delete zXings;
397 if (breakout)
398 break;
399
400 }
401 //-------------------------------------
402 // end of loop over pixels
403 //-------------------------------------
404 if (ProduceGraphic){
405 if (ev % 500 == 0){
406
407 //OVERLAY PEAKS
408 cPulsetypes->cd(1);
409 gPad->SetGrid();
410 hSinglePeakOverlay->Draw("COLZ");
411
412 cPulsetypes->cd(2);
413 gPad->SetGrid();
414 hDoublePeakOverlay->Draw("COLZ");
415
416 cPulsetypes->cd(3);
417 gPad->SetGrid();
418 hTripplePeakOverlay->Draw("COLZ");
419
420 cPulsetypes->cd(4);
421 gPad->SetGrid();
422 hLargePeakOverlay->Draw("COLZ");
423 // cPulsetypes->cd(5);
424 // hSinglePeakMaximum->Draw("");
425 cPulsetypes->Modified();
426 cPulsetypes->Update();
427 }
428 }
429 if (breakout)
430 break;
431 }
432
433//-------------------------------------
434// end of loop over Events
435//-------------------------------------
436
437//-------------------------------------
438// Histogramm of Maximas in Overlay Spectra
439//-------------------------------------
440cout << "producing...Maximum peaks" << endl;
441 for (unsigned int cnt = 1 ; cnt <= 4 ; cnt ++){
442 cout << "peak type: " << cnt << endl;
443 if (cnt == 1){
444 hMaximumTemp = hSinglePeakMaximum;
445 hOverlayTemp = hSinglePeakOverlay;
446 }
447 else if (cnt == 2){
448 hMaximumTemp = hDoublePeakMaximum;
449 hOverlayTemp = hDoublePeakOverlay;
450 }
451 else if (cnt == 3){
452 hMaximumTemp = hTripplePeakMaximum;
453 hOverlayTemp = hTripplePeakOverlay;
454 }
455 else if (cnt == 4){
456 hMaximumTemp = hLargePeakMaximum;
457 hOverlayTemp = hLargePeakOverlay;
458 }
459 // resize vector SingleMaximum to number of timeslices in Overlay Spectra
460 SingleMaximum.resize(hOverlayTemp->GetNbinsX());
461 //put maximumvalue of every Y-projection of every timeslice into vector
462 for (int TimeSlice = 0; TimeSlice <= hOverlayTemp->GetNbinsX(); TimeSlice++ ){
463 TString histotitle;
464 histotitle += "hproj_py";
465 histotitle += cnt ;
466 TH1D* hProjPeak = hOverlayTemp->ProjectionY(histotitle, TimeSlice, TimeSlice, "");
467 SingleMaximum[ TimeSlice ].maxAmpl = (hProjPeak->GetMaximumBin() * 0.5) - 3.5;
468 SingleMaximum[ TimeSlice ].countOfMax = hProjPeak->GetBinContent( hProjPeak->GetMaximumBin() );
469 hMaximumTemp->SetBinContent(TimeSlice, SingleMaximum[ TimeSlice ].maxAmpl );
470 }
471 hMaximumTemp->Fit("landau", "", "", -50, 250);
472 }
473//-------------------------------------
474// Draw Histograms
475//-------------------------------------
476 if (ProduceGraphic){
477 //OVERLAY PEAKS
478cout << "producing...So" << endl;
479 cPulsetypes->cd(1);
480 gPad->SetGrid();
481 hSinglePeakOverlay->ResetStats();
482 hSinglePeakOverlay->Draw("COLZ");
483cout << "producing...DO" << endl;
484 cPulsetypes->cd(2);
485 gPad->SetGrid();
486 hDoublePeakOverlay->ResetStats();
487 hDoublePeakOverlay->Draw("COLZ");
488cout << "producing...TO" << endl;
489 cPulsetypes->cd(3);
490 gPad->SetGrid();
491 hTripplePeakOverlay->ResetStats();
492 hTripplePeakOverlay->Draw("COLZ");
493cout << "producing...LO" << endl;
494 cPulsetypes->cd(4);
495 gPad->SetGrid();
496 hLargePeakOverlay->ResetStats();
497 hLargePeakOverlay->Draw("COLZ");
498cout << "producing...SM" << endl;
499 cPulsetypes->cd(5);
500 hSinglePeakMaximum->ResetStats();
501 hSinglePeakMaximum->Draw("");
502cout << "producing...DM" << endl;
503 cPulsetypes->cd(6);
504 hDoublePeakMaximum->ResetStats();
505 hDoublePeakMaximum->Draw("");
506cout << "producing...TM" << endl;
507 cPulsetypes->cd(7);
508 hTripplePeakMaximum->ResetStats();
509 hTripplePeakMaximum->Draw("");
510cout << "producing...LM" << endl;
511 cPulsetypes->cd(8);
512 hLargePeakMaximum->ResetStats();
513 hLargePeakMaximum->Draw("");
514cout << "producing...Done" << endl;
515
516
517
518
519 cPulsetypes->Modified();
520 cPulsetypes->Update();
521
522// cPixelPeakOverlay->cd();
523// hPixelPeakOverlay->Draw("COLZ");
524// cPixelPeakOverlay->Modified();
525// cPixelPeakOverlay->Update();
526
527 }
528 SaveHistograms( OutRootFileName );
529
530 return( 0 );
531}
532
533// booking and parameter settings for all histos
534void BookHistos( ){
535
536 debugHistos = new TH1F[ NumberOfDebugHistoTypes ];
537 for ( int type = 0; type < NumberOfDebugHistoTypes; type++){
538 debugHistos[ type ].SetBins(1024, 0, 1024);
539 debugHistos[ type ].SetLineColor(1);
540 debugHistos[ type ].SetLineWidth(2);
541
542 // set X axis paras
543 debugHistos[ type ].GetXaxis()->SetLabelSize(0.1);
544 debugHistos[ type ].GetXaxis()->SetTitleSize(0.1);
545 debugHistos[ type ].GetXaxis()->SetTitleOffset(1.2);
546 debugHistos[ type ].GetXaxis()->SetTitle(Form("Time slice (%.1f ns/slice)", 1./2.));
547
548 // set Y axis paras
549 debugHistos[ type ].GetYaxis()->SetLabelSize(0.1);
550 debugHistos[ type ].GetYaxis()->SetTitleSize(0.1);
551 debugHistos[ type ].GetYaxis()->SetTitleOffset(0.3);
552 debugHistos[ type ].GetYaxis()->SetTitle("Amplitude (a.u.)");
553 }
554/*
555 ProjHistos = new TH1D[ NumberOfTimeSlices ];
556 for ( int type = 0; type < NumberOfTimeSlices; type++){
557 ProjHistos[ type ].SetBins(1024, 0, 1024);
558 ProjHistos[ type ].SetLineColor(1);
559 ProjHistos[ type ].SetLineWidth(2);
560
561 // set X axis paras
562 ProjHistos[ type ].GetXaxis()->SetLabelSize(0.1);
563 ProjHistos[ type ].GetXaxis()->SetTitleSize(0.1);
564 ProjHistos[ type ].GetXaxis()->SetTitleOffset(1.2);
565 ProjHistos[ type ].GetXaxis()->SetTitle(Form("Time slice (%.1f ns/slice)", 1./2.));
566
567 // set Y axis paras
568 ProjHistos[ type ].GetYaxis()->SetLabelSize(0.1);
569 ProjHistos[ type ].GetYaxis()->SetTitleSize(0.1);
570 ProjHistos[ type ].GetYaxis()->SetTitleOffset(0.3);
571 ProjHistos[ type ].GetYaxis()->SetTitle("Amplitude (a.u.)");
572 }
573*/
574//-----------------------------------------------------------------------------
575 hSinglePeakOverlay = new TH2F("hSinglePeakOverlay", "Overlay of detected Single Photon Peaks",
576 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 ,
577 512, -55.5, 300.5);
578 hSinglePeakOverlay->SetAxisRange(-5.5, 35.5, "Y");
579 hSinglePeakOverlay->GetXaxis()->SetTitle( "Timeslices" );
580 hSinglePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
581 //hSinglePeakOverlay->SetBit(TH2F::kCanRebin);
582 hList.Add( hSinglePeakOverlay );
583
584 hSinglePeakMaximum = new TH1F("hSinglePeakMaximum", "Single Peak derived with Maximum of above Spektrum",
585 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight ,
586 (-1*hPeakOverlayXaxisLeft)-0.5,
587 hPeakOverlayXaxisRight-0.5
588 );
589 hSinglePeakMaximum->SetAxisRange(-0.5, 25.5, "Y");
590 hSinglePeakMaximum->GetXaxis()->SetTitle( "Timeslices" );
591 hSinglePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
592 hList.Add( hSinglePeakMaximum );
593//-----------------------------------------------------------------------------
594
595 hDoublePeakOverlay = new TH2F("hDoublePeakOverlay", "Overlay of detected Double Photon Peaks",
596 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 ,
597 512, -55.5, 300.5);
598 hDoublePeakOverlay->SetAxisRange(-5.5, 35.5, "Y");
599 hDoublePeakOverlay->GetXaxis()->SetTitle( "Timeslices" );
600 hDoublePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
601 //hSinglePeakOverlay->SetBit(TH2F::kCanRebin);
602 hList.Add( hDoublePeakOverlay );;
603
604 hDoublePeakMaximum = new TH1F("hSinglePeakMaximum", "Double Peak derived with Maximum of above Spektrum",
605 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight ,
606 (-1*hPeakOverlayXaxisLeft)-0.5,
607 hPeakOverlayXaxisRight-0.5
608 );
609 hDoublePeakMaximum->SetAxisRange(-0.5, 25.5, "Y");
610 hDoublePeakMaximum->GetXaxis()->SetTitle( "Timeslices" );
611 hDoublePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
612 hList.Add( hDoublePeakMaximum );
613 //-----------------------------------------------------------------------------
614
615 hTripplePeakOverlay= new TH2F("hTripplePeakOverlay", "Overlay of detected Tripple Photon Peaks",
616 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 ,
617 512, -55.5, 300.5);
618 hTripplePeakOverlay->SetAxisRange(-5.5, 35.5, "Y");
619 hTripplePeakOverlay->GetXaxis()->SetTitle( "Timeslices" );
620 hTripplePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
621 //hTripplePeakOverlay->SetBit(TH2F::kCanRebin);
622 hList.Add( hTripplePeakOverlay );;
623
624 hTripplePeakMaximum = new TH1F("hSinglePeakMaximum", "Triple Peak derived with Maximum of above Spektrum",
625 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight ,
626 (-1*hPeakOverlayXaxisLeft)-0.5,
627 hPeakOverlayXaxisRight-0.5
628 );
629 hTripplePeakMaximum->SetAxisRange(-0.5, 25.5, "Y");
630 hTripplePeakMaximum->GetXaxis()->SetTitle( "Timeslices" );
631 hTripplePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
632 hList.Add( hTripplePeakMaximum );
633 //-----------------------------------------------------------------------------
634
635 hLargePeakOverlay= new TH2F("hLargePeakOverlay", "Overlay of detected Multi Photon Peaks",
636 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 ,
637 512, -5.5, 300.5 );
638 hLargePeakOverlay->SetAxisRange(-5.5, 200.5, "Y");
639 hLargePeakOverlay->GetXaxis()->SetTitle( "Timeslices" );
640 hLargePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
641 //hLargePeakOverlay->SetBit(TH2F::kCanRebin);
642 hList.Add( hLargePeakOverlay );;
643
644 hLargePeakMaximum = new TH1F("hLargePeakMaximum", "Peak derived with Maximum of above Spektrum",
645 hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight ,
646 (-1*hPeakOverlayXaxisLeft)-0.5,
647 hPeakOverlayXaxisRight-0.5
648 );
649 hLargePeakMaximum->SetAxisRange(-0.5, 50.5, "Y");
650 hLargePeakMaximum->GetXaxis()->SetTitle( "Timeslices" );
651 hLargePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
652 hList.Add( hLargePeakMaximum );
653 //-----------------------------------------------------------------------------
654
655// hPixelPeakOverlay = new TH2F("hPixelPeakOverlay", "Maximum of Statistic of overlayed Peaks",
656// hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 ,
657// 512, -55.5, 200.5 );
658// hPixelPeakOverlay->GetXaxis()->SetTitle( "Timeslices" );
659// hPixelPeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
660// hList.Add( hPixelPeakOverlay );
661
662// hEventPeakOverlay = new TH2F("hEventPeakOverlay", "Overlay of detected Peaks of all Pixel of one Event",
663// hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 ,
664// 4096, -48.5, 200.5 );
665// hEventPeakOverlay->GetXaxis()->SetTitle( "Timeslices" );
666// hEventPeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" );
667// hList.Add( hEventPeakOverlay );
668
669}
670
671void SaveHistograms(const char * loc_fname ){
672 // create the histogram file (replace if already existing)
673 TFile tf( loc_fname, "RECREATE");
674
675 hList.Write(); // write the major histograms into the top level directory
676
677 tf.Close(); // close the file
678} // end of SaveHistograms( char * loc_fname )
679
Note: See TracBrowser for help on using the repository browser.