source: fact/tools/rootmacros/PulseTemplates/FPulseOverlay.C@ 13559

Last change on this file since 13559 was 13558, checked in by Jens Buss, 13 years ago
now histograms will be plotted again
File size: 35.3 KB
Line 
1/********************** FPulseOverlay ***********************
2* read FACT raw data
3* remove spikes
4* calculate baseline
5* find peaks (CFD and normal discriminator)
6+ search for Peaks in data
7+ put peaks into different histos depending on Amplitude
8+ draw histos for single, double, tripple, ....Photonpulses
9+ draw Parameterdevelopment
10 so it can be used for detection of other peaks
11*************************************************************/
12//----------------------------------------------------------------------------
13// root libraries
14//----------------------------------------------------------------------------
15
16#include <TROOT.h>
17#include <TCanvas.h>
18#include <TProfile.h>
19#include <TTimer.h>
20#include <TH1F.h>
21#include <TH2F.h>
22#include <Getline.h>
23#include <TLine.h>
24#include <TBox.h>
25#include <TMath.h>
26#include <TFile.h>
27#include <TStyle.h>
28#include <TString.h>
29#include <TObjArray.h>
30
31#include <stdio.h>
32#include <stdint.h>
33#include <cstdio>
34#include <iostream>
35#include <fstream>
36using namespace std;
37
38#define NPIX 1440
39#define NCELL 1024
40#define FAD_MAX_SAMPLES 1024
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.h"
57//#include "../factfir.C"
58
59//#include "../DrsCalibration.C"
60#include "../DrsCalibration.h"
61
62#include "../SpikeRemoval.h"
63//#include "../SpikeRemoval.C"
64
65#include "rootfilehandler.h"
66//#include "rootfilehandler.C"
67
68#include "pixel.h"
69//#include "pixel.C"
70
71
72//----------------------------------------------------------------------------
73// Decleration of global variables
74//----------------------------------------------------------------------------
75
76bool breakout=false;
77
78vector<int16_t> AllPixelDataVector;
79vector<int16_t> StartCellVector;
80unsigned int CurrentEventID;
81size_t PXLxROI;
82UInt_t gRegionOfInterest;
83UInt_t NumberOfPixels;
84TString histotitle;
85
86size_t TriggerOffsetROI, RC;
87size_t drs_n;
88vector<float> drs_basemean;
89vector<float> drs_gainmean;
90vector<float> drs_triggeroffsetmean;
91
92vector<float> Ameas(FAD_MAX_SAMPLES); // copy of the data (measured amplitude
93vector<float> Vcorr(FAD_MAX_SAMPLES); // corrected Values
94vector<float> Vslide(FAD_MAX_SAMPLES); // sliding average result
95vector<float> Vcfd(FAD_MAX_SAMPLES); // CDF result
96vector<float> Vcfd2(FAD_MAX_SAMPLES); // CDF result + 2nd sliding average
97
98int gNEvents=0;
99float gGainMean = 9;
100float gBSLMean = 0;
101
102//typedef struct{
103// double maxAmpl;
104// int countOfMax;
105// } OverlayMaximum;
106
107//typedef struct{
108// TString name;
109// TString title;
110// TString xTitle;
111// TString yTitle;
112// float yMax;
113// float yMin;
114//} histAttrib_t;
115
116//histAttrib_t* gHistoAttrib[maxPulseOrder];
117//histAttrib_t* gProfileAttrib[maxPulseOrder];
118
119
120// histograms
121const int Number_Of_Debug_Histo_Types = 7;
122
123const unsigned int
124 Ameas_ = 0,
125 N1mean_ = 1,
126 Vcorr_ = 2,
127 Vtest_ = 3,
128 Vslide_ = 4,
129 Vcfd_ = 5,
130 Vcfd2_ = 6;
131
132//const char* gHistoTypes[8] = {
133// "hMaxOverlay",
134// "hPixelMax",
135// "hEdgeOverlay",
136// "hMaxProfile",
137// "hAllPixelOverlay",
138// "hAllPixelMax",
139// "hAllPixelMaxGaus",
140// "hAllPixelProfile"
141//};
142
143//----------------------------------------------------------------------------
144// Initialisation of histograms
145//----------------------------------------------------------------------------
146
147 // Temporary Objects
148 TH1F* debugHistos = NULL;
149 //TH2F* hOverlayTemp = NULL;
150 //TH1D* hProjPeak = NULL;
151 TH1F* hTesthisto = NULL;
152 TH2F* hTesthisto2 = NULL;
153
154 //Histogram Parameters
155 Int_t gPixelOverlayXaxisLeft = 0;
156 Int_t gPixelOverlayXaxisRight = 0;
157
158 //Root-File Objects
159// TObjArray* hList[pixelSetSize] = {NULL};
160 TObjArray* hRootList = NULL;
161
162 TCanvas** cgpPixelPulses = NULL;
163 TCanvas** cgpDistributions = NULL;
164 TCanvas* cFiltered = NULL;
165 TCanvas* cgpTestHistos = NULL;
166
167//----------------------------------------------------------------------------
168// Functions
169//----------------------------------------------------------------------------
170
171void BookDebugHistos(int );
172void BookTestHistos( int );
173
174//void DeletePixelHistos(int );
175//void SaveHistograms( const char*, const char*, TObjArray*, int );
176
177void FillHistograms(Pixel*, vector<Region>*, int, int, int , int );
178void DrawTestHistograms( int);
179void ProduceDebugHistos( vector<Region> *pZXings);
180bool UseThisPulse( int, int, int, int, int, int);
181void UpdateCanvases( int, int);
182void DeletePixelCanvases( int, int );
183
184//----------------------------------------------------------------------------
185//----------------------------------------------------------------------------
186// MAIN - Funtion
187//----------------------------------------------------------------------------
188//----------------------------------------------------------------------------
189int FPulseOverlay(
190 const char* datafilename = "/fact/raw/2011/11/09/20111109_006.fits.gz",
191 const char* drsfilename = "/fact/raw/2011/11/09/20111109_003.drs.fits.gz",
192 const char* OutRootFileName = "/home_nfs/isdc/jbbuss/analysis/FPulseTemplate/test.root",
193 bool ProduceGraphic = true,
194 bool spikeDebug = false,
195 bool debugPixel = true,
196 bool testmode = false,
197 bool saveResults = false,
198 int verbosityLevel = 1, // different verbosity levels can be implemented here
199 int firstevent = 1, // Hast to be between 1 and infinity NOT 0!!!!
200 int nevents = 1000,
201 int firstpixel = 15,
202 int npixel = 3,
203 int pixelSetSize = 40,
204 int maxPulseOrder = 4,
205 int AmplWindowWidth = 14, //Width of Window for selection of pulses to histograms
206 TString histoOptions = "S",
207 float GainMean = 9,
208 float BSLMean = -1, //4 Histogramms will be drawn, decide how far pulsehights differ from eachother
209 int avg1 = 8,
210 int avg2 = 8,
211 int OverlayWindowLeft = 70,
212 int OverlayWindowRight = 230,
213 int refresh_rate = 500 //refresh rate for canvases
214 )
215{
216
217//----------------------------------------------------------------------------
218// Save-Root-File Settings
219//----------------------------------------------------------------------------
220 if (saveResults)
221 {
222 CreateRootFile( OutRootFileName, verbosityLevel );
223 }
224
225//----------------------------------------------------------------------------
226// root global Settings
227//----------------------------------------------------------------------------
228
229 gGainMean = GainMean;
230 gBSLMean = BSLMean;
231 gPixelOverlayXaxisLeft = OverlayWindowLeft;
232 gPixelOverlayXaxisRight = OverlayWindowRight;
233
234 gStyle->SetPalette(1,0);
235 gROOT->SetStyle("Plain");
236// gPad->SetGrid();
237
238//----------------------------------------------------------------------------
239// root Canvas Settings
240//----------------------------------------------------------------------------
241 //Canvas Pad numbering
242 int pulsesCanvasFrameNrs[4] = {
243 1, // Top left
244 2, // Top right
245 3, // bottom left
246 4 // bootom right
247 };
248
249 //Canvas Pad numbering
250 int distributionCanvasFrameNrs[4] = {
251 1, // Top left
252 2, // Top right
253 3, // bottom left
254 4 // bootom right
255 };
256
257 if (ProduceGraphic)
258 {
259
260 //Canvases
261 cgpPixelPulses = new TCanvas*[maxPulseOrder];
262 cgpDistributions = new TCanvas*[maxPulseOrder];
263
264 //TCanvas* gpcDevelopment = NULL;
265 TString cName = "";
266 TString cTitle = "";
267
268 //naming of pulse canvases
269 for (
270 int pulse_order = maxPulseOrder - 1;
271 pulse_order >= 0 ;
272 pulse_order--
273 )
274 {
275 cName ="cgpDistributions";
276 cName += pulse_order;
277 cTitle ="Distributions of Pulses with Order of: ";
278 cTitle += pulse_order;
279 cgpDistributions[pulse_order] = new TCanvas(cName,cTitle, 0,pulse_order*20,800,800);
280 cgpDistributions[pulse_order]->Divide(2, 2);
281 cName ="cgpPixelPulses";
282 cName += pulse_order;
283 cTitle ="Overlays of Pulses with Order of: ";
284 cTitle += pulse_order;
285 cgpPixelPulses[pulse_order] = new TCanvas(cName,cTitle, 0,pulse_order*20,800,800);
286 cgpPixelPulses[pulse_order]->Divide(2, 2);
287 }
288
289
290 // Create (pointer to) Canvases, which are used in every run,
291 // also in 'non-debug' runs
292 // Canvases only need if spike Debug, but I want to deklare
293 // the pointers anyway ...
294 if (spikeDebug)
295 {
296 cFiltered = new TCanvas("cFiltered","filtered DRS Waveforms", 1,310,400,300);
297 cFiltered->Divide(1, 3);
298 }
299
300 //additional Test histograms
301 cgpTestHistos = new TCanvas( "cgpTestHistos", "Test Histograms", 801, 0, 800, 800 );
302 cgpTestHistos->Divide(2,0);
303 }
304
305//-----------------------------------------------------------------------------
306// Filter-Settings
307//-----------------------------------------------------------------------------
308// CFD filter settings
309 int k_cfd = 10;
310 vector<double> a_cfd(k_cfd, 0);
311 double b_cfd = 1.;
312 a_cfd[0]=-0.75;
313 a_cfd[k_cfd-1]=1.;
314
315//-----------------------------------------------------------------------------
316// prepare datafile
317//-----------------------------------------------------------------------------
318
319// Open the data file
320
321 fits * datafile = NULL;
322 // Opens the raw data file and 'binds' the variables given as
323 // Parameters to the data file. So they are filled with
324 // raw data as soon as datafile->GetRow(int) is called.
325 gNEvents = openDataFits(
326 datafilename,
327 &datafile,
328 AllPixelDataVector,
329 StartCellVector,
330 CurrentEventID,
331 gRegionOfInterest,
332 NumberOfPixels,
333 PXLxROI,
334 verbosityLevel
335 );
336
337 if (gNEvents == 0)
338 {
339 cout << "return code of OpenDataFile:" << datafilename<< endl;
340 cout << "is zero -> aborting." << endl;
341 return 1;
342 }
343
344 if (verbosityLevel > 0)
345 {
346 cout << endl <<"number of events in file: "<< gNEvents << "\t";
347 }
348
349 if ( nevents == -1 || nevents > gNEvents ) nevents = gNEvents; // -1 means all!
350
351 if (verbosityLevel > 0)
352 {
353 cout <<"of, which "<< nevents << " will be processed"<< endl;
354 cout <<"Total # of Pixel: "<< NumberOfPixels << "\t";
355 }
356
357 if ( npixel == -1 || npixel > (int)NumberOfPixels ) npixel = (int)NumberOfPixels; // -1 means all!
358
359 if (verbosityLevel > 0)
360 {
361 cout <<"of, which "<< npixel << " will be processed"<< endl;
362 }
363
364 //Get the DRS calibration
365 RC = openCalibFits(
366 drsfilename,
367 drs_basemean,
368 drs_gainmean,
369 drs_triggeroffsetmean,
370 TriggerOffsetROI
371 );
372
373 if (RC == 0)
374 {
375 cout << "return code of openCalibFits:" << drsfilename << endl;
376 cout << "is zero -> aborting." << endl;
377 return 1;
378 }
379
380//-----------------------------------------------------------------------------
381// Book the histograms
382//-----------------------------------------------------------------------------
383
384 if (spikeDebug)
385 BookDebugHistos(verbosityLevel );
386
387 if (testmode)
388 BookTestHistos( verbosityLevel );
389
390
391//-----------------------------------------------------------------------------
392// Main Cycle
393//-----------------------------------------------------------------------------
394
395//----------------------------------------------------------------------------
396// Initialize Pixel
397//----------------------------------------------------------------------------
398 Pixel** pixel = new Pixel*[NPIX];
399
400 for (int i = 0 ; i < NPIX; i++)
401 {
402 pixel[i] = NULL;
403 }
404//-------------------------------------
405// Loop over Pixel Sets
406//-------------------------------------
407 for ( int firstPixelOfSet = firstpixel;
408 firstPixelOfSet < firstpixel + npixel;
409 firstPixelOfSet += pixelSetSize )
410 {
411
412 if (verbosityLevel == 0)
413 {
414 cout << "------------------------------------------------" << endl
415 << "...processing Set from Pixel "
416 << firstPixelOfSet
417 << " to Pixel "
418 << firstPixelOfSet+pixelSetSize-1 << endl;
419 }
420
421//--------------------------------------------------------------------
422// Loop over every Event of Pixel Set
423//--------------------------------------------------------------------
424 for ( int ev = firstevent ; ev < firstevent + nevents; ev++)
425 {
426 // Get an Event --> consists of 1440 Pixel ...erm....data
427 datafile->GetRow( ev - 1 );
428
429 if (verbosityLevel == 1)
430 {
431 cout << "-------------------------------------" << endl
432 << "...processing Set from Pixel "
433 << firstPixelOfSet
434 << " to Pixel "
435 << firstPixelOfSet+pixelSetSize-1
436 << "... Event: " << CurrentEventID
437 << "/" << nevents << endl;
438 }
439
440//--------------------------------------------------------------------
441// Loops over every Pixel of a Set of Pixels
442//--------------------------------------------------------------------
443 for ( int pixelID = firstPixelOfSet;
444 pixelID < firstPixelOfSet + pixelSetSize
445 && pixelID < firstpixel + npixel;
446 pixelID++ )
447 {
448 if (verbosityLevel > 1)
449 {
450 cout << "-------------------------------------" << endl
451 << "...processing Set from Pixel "
452 << firstPixelOfSet
453 << " to Pixel "
454 << firstPixelOfSet+pixelSetSize-1
455 << "... Event: " << CurrentEventID
456 << "/" << nevents << endl
457 << " Pixel: " << pixelID
458 << "/" << firstpixel + npixel -1 << endl;
459 }
460
461//-------------------------------------
462// Create individual Pixel
463//-------------------------------------
464 if (ev == firstevent)
465 {
466 pixel[pixelID] = new Pixel(
467 pixelID,
468 maxPulseOrder,
469 verbosityLevel,
470 gPixelOverlayXaxisLeft,
471 gPixelOverlayXaxisRight,
472 gBSLMean,
473 gGainMean,
474 histoOptions
475 );
476 }
477
478//-------------------------------------
479// Apply Calibration
480//-------------------------------------
481 if (verbosityLevel > 2) cout << "...applying DrsCalibration";
482 applyDrsCalibration(
483 Ameas,
484 pixelID,
485 12,
486 12,
487 drs_basemean,
488 drs_gainmean,
489 drs_triggeroffsetmean,
490 gRegionOfInterest,
491 AllPixelDataVector,
492 StartCellVector
493 );
494 if (verbosityLevel > 2) cout << "...done " << endl;
495
496//-------------------------------------
497// Apply Filters
498//-------------------------------------
499 // finds spikes in the raw data, and interpolates the value
500 // spikes are: 1 or 2 slice wide, positive non physical artifacts
501 if (verbosityLevel > 2) cout << "...removeing Spikes";
502 removeSpikes (Ameas, Vcorr);
503 if (verbosityLevel > 2) cout << "...done " << endl;
504
505 // filter Vcorr with sliding average using FIR filter function
506 if (verbosityLevel > 2) cout << "...applying sliding average filter";
507 sliding_avg(Vcorr, Vslide, avg1);
508 if (verbosityLevel > 2) cout << "...done " << endl;
509
510 // filter Vslide with CFD using FIR filter function
511 if (verbosityLevel > 2) cout << "...apllying factfir filter";
512 factfir(b_cfd , a_cfd, k_cfd, Vslide, Vcfd);
513 if (verbosityLevel > 2) cout << "...done " << endl;
514
515 // filter Vcfd with sliding average using FIR filter function
516 if (verbosityLevel > 2) cout << "...applying 2nd sliding average filter";
517 sliding_avg(Vcfd, Vcfd2, avg2);
518 if (verbosityLevel > 2) cout << "...done " << endl;
519
520//-------------------------------------
521// Search vor Zero crossings
522//-------------------------------------
523 if (verbosityLevel > 2) cout << endl << "...searching zero crossings" ;
524 // peaks in Ameas[] are found by searching for zero crossings
525 // in Vcfd2
526 // first Argument 1 means ... *rising* edge
527 // second Argument 1 means ... search with stepsize 1 ... 10 is okay as well
528 vector<Region>* pZXings = zerosearch( Vcfd2 , 1 , 8);
529 // pZXings means "zero cross ings"
530 EnlargeRegion(*pZXings, 10, 10);
531 findAbsMaxInRegions(*pZXings, Vslide);
532 removeMaximaBelow( *pZXings, 3.0);
533 removeRegionWithMaxOnEdge( *pZXings, 2);
534 removeRegionOnFallingEdge( *pZXings, 100);
535 findTimeOfHalfMaxLeft(*pZXings, Vcorr, gBSLMean, 5, 10, verbosityLevel );
536 if (verbosityLevel > 2) cout << "...done" << endl;
537
538//-----------------------------------------------------------------------------
539// Fill Overlay Histos
540//-----------------------------------------------------------------------------
541 if (verbosityLevel > 2)
542 {
543 cout << endl << "Filling Histograms of Pixel# [Set] " << pixelID
544 << " Pixel# [Got] " << pixel[pixelID]->mChid;
545 }
546 FillHistograms(
547 pixel[pixelID],
548 pZXings,
549 AmplWindowWidth,
550 ev,
551// histoOptions,
552 maxPulseOrder,
553 verbosityLevel
554 );
555//-----------------------------------------------------------------------------
556// Spike Debug
557//-----------------------------------------------------------------------------
558 if ( spikeDebug )
559 {
560 ProduceDebugHistos( pZXings );
561 }// end of if(spikeDebug)
562
563 delete pZXings;
564 if (breakout) break;
565//-------------------------------------
566// Draw 1. Set of Pixel Histograms
567//-------------------------------------
568 if ((ev % refresh_rate) == 0)
569 {
570 if (ProduceGraphic)
571 {
572 if (debugPixel)
573 {
574// for ( int order = 0;
575// order < maxPulseOrder;
576// order++
577// )
578// {
579 pixel[firstPixelOfSet]->DrawHistograms(
580 cgpPixelPulses,
581 pulsesCanvasFrameNrs
582 );
583// }
584 if (testmode)
585 {
586 DrawTestHistograms(
587 verbosityLevel
588 );
589 }
590
591 UpdateCanvases(
592 verbosityLevel,
593 maxPulseOrder
594 );
595 }
596 }
597 }
598
599 if (breakout) break;
600
601 if (verbosityLevel > 1)
602 {
603 cout << endl << "...End of Pixel"
604 << endl << "-------------------------------------"<< endl;
605 }
606 }//End of Loop over Set
607//-------------------------------------
608// end of Loops over Set
609//-------------------------------------
610
611 }//End of Loop over Events
612//-------------------------------------
613// end of Loops over Events
614//-------------------------------------
615
616
617
618//-------------------------------------
619// Draw Pixel Histogramms of Overlay Spectra
620//-------------------------------------
621
622 for ( int pixelID = firstPixelOfSet;
623 pixelID < firstPixelOfSet + pixelSetSize
624 && pixelID < firstpixel + npixel;
625 pixelID++ )
626 {
627 //here is what happends at the end of each loop over all Events
628
629 if (ProduceGraphic)
630 {
631 pixel[pixelID]->DrawHistograms(
632 cgpPixelPulses,
633 pulsesCanvasFrameNrs
634 );
635
636 UpdateCanvases(
637 verbosityLevel,
638 maxPulseOrder
639 );
640 }
641
642 //Save Histograms of Pixels into Output rootfile
643 if (verbosityLevel > 2)
644 {
645 cout << endl << "Saving Pixel# [Set] " << pixelID
646 << " Pixel# [Got] " << pixel[pixelID]->mChid;
647 }
648 pixel[pixelID]->SavePixelHistograms( OutRootFileName, saveResults );
649
650
651 if (debugPixel)
652 {
653 //Process gui events asynchronously during input
654 cout << endl;
655 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
656 timer.TurnOn();
657 TString input = Getline("Type 'q' to exit, <return> to go on: ");
658 timer.TurnOff();
659 if (input=="q\n") {
660 break;
661 }
662 }
663
664 if (verbosityLevel > 2)
665 {
666 cout << endl << "Deleting Pixel# [Set] " << pixelID
667 << " Pixel# [Got] " << pixel[pixelID]->mChid;
668 }
669 delete pixel[pixelID];
670
671 if (verbosityLevel > 2)
672 {
673 cout << endl << "...End of Pixel-Set"
674 << endl << "------------------------------------------------"
675 << endl;
676 }
677 }
678 }
679 // End of Loop over all Pixels
680
681//-------------------------------------
682// Draw Histograms
683//-------------------------------------
684
685 SaveHistograms( //save histograms of generell results into output root file
686 OutRootFileName,
687 "root",
688 hRootList,
689 saveResults,
690 verbosityLevel
691 );
692
693 if (ProduceGraphic)
694 {
695 UpdateCanvases(
696 verbosityLevel,
697 maxPulseOrder
698 );
699
700 //Process gui events asynchronously during input
701 cout << endl;
702 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
703 timer.TurnOn();
704 Getline("Press <return> to go on");
705 timer.TurnOff();
706
707 DeletePixelCanvases(
708 maxPulseOrder,
709 verbosityLevel
710 );
711 delete cgpTestHistos;
712
713 if (spikeDebug)
714 delete cFiltered;
715 }
716
717 delete datafile;
718 delete[] pixel;
719
720
721 return( 0 );
722}
723//----------------------------------------------------------------------------
724// end of main function
725//-----------------------------------------------------------------------------
726
727
728
729
730//-----------------------------------------------------------------------------
731// Funktions
732//-----------------------------------------------------------------------------
733
734/*
735ich erzeuge sowas wie 36 pixel mit ca 20 Histogrammen pro pixel
736wie viel speicher braucht das?
737Ich brauche eine möglichkeit jedes dieser histogramme zu identifizieren
738am besten ein array oder eine liste oder einen abstracten datentyp mit den histogrammen als attribute.
739ne klasse wäre super
740
741*/
742
743void BookTestHistos( int verbosityLevel )
744{
745 if (verbosityLevel > 2) cout << endl << "...book pixel histograms" << endl;
746
747 hTesthisto = new TH1F (
748 "hTesthisto",
749 "Deviation of rising edge and maximum",
750 600,
751 -10.1,
752 10.1
753 );
754 hTesthisto->GetXaxis()->SetTitle( "Timeslices [a.u.]" );
755 hTesthisto->GetYaxis()->SetTitle( "counts" );
756
757 hTesthisto2 = new TH2F (
758 "hTesthisto2",
759 "Deviation of rising edge and maximum by event #",
760 gNEvents,
761 250,
762 800,
763 600,
764 -10.1,
765 10.1
766 );
767// hTesthisto2->GetXaxis()->SetTitle( "Timeslices [a.u.]" );
768// hTesthisto2->GetYaxis()->SetTitle( "counts" );
769// hTesthisto2->SetMarkerStyle( 2 );
770
771 if (verbosityLevel > 2) cout << "...done" << endl;
772}
773//end of BookTestHistos
774//----------------------------------------------------------------------------
775
776
777void
778BookDebugHistos( int verbosityLevel )
779{
780 if (verbosityLevel > 2) cout << endl << "...book histograms" << endl;
781 hRootList = new TObjArray;
782 debugHistos = new TH1F[ Number_Of_Debug_Histo_Types ];
783 for ( int type = 0; type < Number_Of_Debug_Histo_Types; type++){
784 debugHistos[ type ].SetBins(1024, 0, 1024);
785 debugHistos[ type ].SetLineColor(1);
786 debugHistos[ type ].SetLineWidth(2);
787 debugHistos[ type ].SetStats(false);
788
789 // set X axis paras
790 debugHistos[ type ].GetXaxis()->SetLabelSize(0.08);
791 debugHistos[ type ].GetXaxis()->SetTitleSize(0.08);
792 debugHistos[ type ].GetXaxis()->SetTitleOffset(1.0);
793 debugHistos[ type ].GetXaxis()->SetTitle(Form("Time slice [%.1f ns/slice]", 1./2.));
794
795 // set Y axis paras
796 debugHistos[ type ].GetYaxis()->SetLabelSize(0.08);
797 debugHistos[ type ].GetYaxis()->SetTitleSize(0.08);
798 debugHistos[ type ].GetYaxis()->SetTitleOffset(0.3);
799 debugHistos[ type ].GetYaxis()->SetTitle("Amplitude [mV]");
800 }
801 if (verbosityLevel > 2) cout << "...done" << endl;
802}
803// end of BookDebugHistos
804//----------------------------------------------------------------------------
805
806
807void
808FillHistograms(
809 Pixel* CurrentPixel,
810 vector<Region>* pZXings,
811 int AmplWindowWidth,
812 int eventNumber,
813 int maxPulseOrder,
814 int verbosityLevel
815 )
816{
817 if (verbosityLevel > 2) cout << endl << "...filling pulse histograms" ;
818 if (verbosityLevel > 3) cout << endl << "...EventNR " << eventNumber ;
819 bool use_this_peak=false;
820 int order_of_pulse=0;
821 vector<Region>::iterator reg;
822 //Loop over all found zerocrossings in Region
823 for (reg = pZXings->begin() ; reg < pZXings->end() ; reg++)
824 {
825 //skip those who are at the Rim of the ROI
826 if (reg->maxPos < 12 || (unsigned int) reg->maxPos > gRegionOfInterest-12)
827 {
828 if (verbosityLevel > 2) cout << endl << "\t...out of range" << endl;
829 continue;
830 }
831
832 // Define axis range of Histogram
833 int Left = reg->maxPos - gPixelOverlayXaxisLeft;
834 int Right = reg->maxPos + gPixelOverlayXaxisRight;
835 int EdgeLeft = reg->halfRisingEdgePos - gPixelOverlayXaxisLeft;
836 int EdgeRight = reg->halfRisingEdgePos + gPixelOverlayXaxisRight;
837
838 if (Left < 0) Left = 0;
839 if (EdgeLeft < 0) EdgeLeft = 0;
840 if (Right > (int)Ameas.size() ) Right = Ameas.size();
841 if (EdgeRight > (int)Ameas.size() ) EdgeRight = Ameas.size();
842
843
844 if (verbosityLevel > 2) cout << endl << "\t...choosing Histogram" ;
845 //determine order of pulse and which histogram shall be filled
846 for ( int order = 0; order < maxPulseOrder; order++)
847 {
848 use_this_peak = UseThisPulse(
849 reg->maxPos,
850 order,
851 AmplWindowWidth,
852 gGainMean,
853 maxPulseOrder,
854 verbosityLevel
855 );
856 if ( use_this_peak)
857 {
858 order_of_pulse = order;
859 break;
860 }
861 }
862 //Fill histograms
863 if (use_this_peak)
864 {
865 //Histograms for Distributions
866 if (CurrentPixel->mOptions.Contains("S") )
867 {
868 if (verbosityLevel > 2)
869 cout << endl << "\t...filling Edge Histogram" ;
870 CurrentPixel->hSlopeRisingEdge[order_of_pulse]->Fill( reg->slopeOfRisingEdge ) ;
871 }
872
873 //Histograms for Maximum Overlay
874 for ( int pos = Left; pos < Right; pos++)
875 {
876 CurrentPixel->hMaxOverlay[order_of_pulse]->Fill( pos - (reg->maxPos), Ameas[pos]) ;
877 CurrentPixel->hMaxProfile[order_of_pulse]->Fill( pos - (reg->maxPos), Ameas[pos]) ;
878 }
879 //Histograms for Edge Overlay
880 for ( int pos = EdgeLeft; pos < EdgeRight; pos++)
881 {
882 CurrentPixel->hEdgeOverlay[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), Ameas[pos]) ;
883 CurrentPixel->hEdgeProfile[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), Ameas[pos]) ;
884 }
885
886 if (verbosityLevel > 2) cout << "...done" << endl;
887 }
888
889 //Breakout option
890 if (breakout)
891 {
892 break;
893 }
894 }
895 if (verbosityLevel > 2) cout << "...done" << endl;
896}
897// end of FillHistograms
898//----------------------------------------------------------------------------
899
900bool UseThisPulse(
901 int maxPos,
902 int order,
903 int AmplWindowWidth,
904 int gGainMean,
905 int maxPulseOrder,
906 int verbosityLevel
907 )
908{
909 //determine if pulse is of given order
910bool use_this_peak = false;
911 if (Ameas[maxPos] >= ((gGainMean*(order+1)) - (AmplWindowWidth/2))
912 && Ameas[maxPos] < ((gGainMean*(order+1)) + AmplWindowWidth/2))
913 {
914 if (verbosityLevel > 3) cout << "...#" << order ;
915 use_this_peak = true;
916 }
917 else
918 {
919 if (verbosityLevel > 3)
920 {
921 if (order >= (maxPulseOrder - 1))
922 {
923 cout << "...NONE" << endl ;
924 }
925 }
926 use_this_peak = false;
927 }
928
929return use_this_peak;
930}
931
932
933void
934DrawTestHistograms(
935 int verbosityLevel
936 )
937{
938 if (verbosityLevel > 2) cout << endl << "...drawing pulse histograms TODO CANVAS" ;
939
940// cgpTestHistos->cd(1);
941// hTesthisto->Draw();
942// cgpTestHistos->cd(2);
943// hTesthisto2->Draw("BOX");
944}
945// end of DrawTestHistograms
946//----------------------------------------------------------------------------
947
948void
949UpdateCanvases(
950 int verbosityLevel,
951 int max_pulse_order
952 )
953{
954 if (verbosityLevel > 3) cout << endl << "...updating canvases" ;
955 for (int pulse_order = 0; pulse_order < max_pulse_order; pulse_order++)
956 {
957 cgpPixelPulses[pulse_order]->Modified();
958 cgpPixelPulses[pulse_order]->Update();
959 cgpTestHistos->Modified();
960 cgpTestHistos->Update();
961 }
962}
963// end of UpdateCanvases
964//----------------------------------------------------------------------------
965
966void
967DeletePixelCanvases(
968 int maxPulseOrder,
969 int verbosityLevel
970 )
971{
972 if (verbosityLevel > 2) cout << endl << "...delete pixel Canvases" ;
973 for (int pulse_order = 0; pulse_order < maxPulseOrder; pulse_order++ )
974 {
975 delete cgpPixelPulses[pulse_order];
976 delete cgpDistributions[pulse_order];
977 }
978 delete[] cgpPixelPulses;
979 delete[] cgpDistributions;
980}
981
982//----------------------------------------------------------------------------
983//----------------------------------------------------------------------------
984
985void
986ProduceDebugHistos(
987 vector<Region>* pZXings
988 )
989{
990 // TODO do this correct. The vectors should be the rigt ones... this is just luck
991 debugHistos[Ameas_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
992 debugHistos[Vcorr_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
993 debugHistos[Vslide_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
994 debugHistos[Vcfd_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
995 debugHistos[Vcfd2_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
996
997 for ( unsigned int sl = 0; sl < gRegionOfInterest; sl++)
998 {
999 debugHistos[Ameas_].SetBinContent(sl, Ameas[sl]);
1000 debugHistos[Vcorr_].SetBinContent(sl, Vcorr[sl]);
1001 debugHistos[Vslide_].SetBinContent( sl, Vslide[sl] );
1002 debugHistos[Vcfd_].SetBinContent( sl, Vcfd[sl] );
1003 debugHistos[Vcfd2_].SetBinContent( sl, Vcfd2[sl] );
1004 }
1005
1006
1007 cFiltered->cd(1);
1008 gPad->SetGrid();
1009 debugHistos[Ameas_].Draw();
1010
1011 cFiltered->cd(2);
1012 gPad->SetGrid();
1013 debugHistos[Vcorr_].Draw();
1014
1015 cFiltered->cd(3);
1016 gPad->SetGrid();
1017 debugHistos[Vslide_].Draw();
1018
1019 TBox *OneBox;
1020 vector<TBox*> MyBoxes;
1021 for (unsigned int i=0; i<pZXings->size(); i++){
1022 OneBox = new TBox(
1023 pZXings->at(i).maxPos -10 ,
1024 pZXings->at(i).maxVal -0.5,
1025 pZXings->at(i).maxPos +10 ,
1026 pZXings->at(i).maxVal +0.5);
1027 OneBox->SetLineColor(kBlue);
1028 OneBox->SetLineWidth(1);
1029 OneBox->SetFillStyle(0);
1030 OneBox->SetFillColor(kRed);
1031 MyBoxes.push_back(OneBox);
1032 OneBox->Draw();
1033 }
1034
1035// cFiltered->cd(3);
1036// gPad->SetGrid();
1037// debugHistos[Vcfd2_].Draw();
1038// TLine *zeroline = new TLine(0, 0, 1024, 0);
1039// zeroline->SetLineColor(kBlue);
1040// zeroline->Draw();
1041
1042 cFiltered->Update();
1043
1044
1045 //Process gui events asynchronously during input
1046 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
1047 timer.TurnOn();
1048 TString input = Getline("Type 'q' to exit, <return> to go on: ");
1049 timer.TurnOff();
1050 if (input=="q\n") {
1051 breakout=true;
1052 }
1053
1054 //TODO!!!!!!!!!
1055 // do some Garbage collection ...
1056 // all the Objects on the heap should be deleted here.
1057
1058}// end of if(spikeDebug)
1059
1060
1061int main()
1062{
1063
1064FPulseOverlay();
1065return 0;
1066
1067}
Note: See TracBrowser for help on using the repository browser.