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

Last change on this file since 13547 was 13547, checked in by neise, 13 years ago
deleted Pixel** pixel
File size: 34.8 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[sampleSetSize] = {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, TString, int, int);
178void DrawTestHistograms( int);
179void ProduceDebugHistos( vector<Region> *pZXings);
180
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 = "test.root",
193 bool ProduceGraphic = false,
194 bool spikeDebug = false,
195 bool debugPixel = false,
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 = 100,
201 int firstpixel = 0,
202 int npixel = -1,
203 int sampleSetSize = 40,
204 int maxPulseOrder = 3,
205 int AmplWindowWidth = 14, //Width of Window for selection of pulses to histograms
206 TString histoOptions = "S",
207 float GainMean = 8,
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 ="Pulses of Order: ";
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 ="Pulses of Order: ";
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 firstPixelOfSample = firstpixel;
408 firstPixelOfSample < firstpixel + npixel;
409 firstPixelOfSample += sampleSetSize )
410 {
411
412 if (verbosityLevel == 0)
413 {
414 cout << "------------------------------------------------" << endl
415 << "...processing Sample from Pixel "
416 << firstPixelOfSample
417 << " to Pixel "
418 << firstPixelOfSample+sampleSetSize-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 cout << endl << "processing event " << ev << endl;
428 datafile->GetRow( ev - 1 );
429
430 if (verbosityLevel == 1)
431 {
432 cout << "-------------------------------------" << endl
433 << "...processing Sample from Pixel "
434 << firstPixelOfSample
435 << " to Pixel "
436 << firstPixelOfSample+sampleSetSize-1
437 << "... Event: " << CurrentEventID
438 << "/" << nevents << endl;
439 }
440
441//--------------------------------------------------------------------
442// Loops over every Pixel of a Set of Pixels
443//--------------------------------------------------------------------
444 for ( int pixelID = firstPixelOfSample;
445 pixelID < firstPixelOfSample + sampleSetSize
446 && pixelID < firstpixel + npixel;
447 pixelID++ )
448 {
449 if (verbosityLevel > 1)
450 {
451 cout << "-------------------------------------" << endl
452 << "...processing Sample from Pixel "
453 << firstPixelOfSample
454 << " to Pixel "
455 << firstPixelOfSample+sampleSetSize-1
456 << "... Event: " << CurrentEventID
457 << "/" << nevents << endl
458 << " Pixel: " << pixelID
459 << "/" << firstpixel + npixel -1 << endl;
460 }
461
462//-------------------------------------
463// Create individual Pixel
464//-------------------------------------
465 if (ev == firstevent)
466 {
467 pixel[pixelID] = new Pixel(
468 pixelID,
469 maxPulseOrder,
470 verbosityLevel,
471 gPixelOverlayXaxisLeft,
472 gPixelOverlayXaxisRight,
473 gBSLMean,
474 gGainMean,
475 histoOptions
476 );
477 }
478
479//-------------------------------------
480// Apply Calibration
481//-------------------------------------
482 if (verbosityLevel > 2) cout << "...applying DrsCalibration";
483 applyDrsCalibration(
484 Ameas,
485 pixelID,
486 12,
487 12,
488 drs_basemean,
489 drs_gainmean,
490 drs_triggeroffsetmean,
491 gRegionOfInterest,
492 AllPixelDataVector,
493 StartCellVector
494 );
495 if (verbosityLevel > 2) cout << "...done " << endl;
496
497//-------------------------------------
498// Apply Filters
499//-------------------------------------
500 // finds spikes in the raw data, and interpolates the value
501 // spikes are: 1 or 2 slice wide, positive non physical artifacts
502 if (verbosityLevel > 2) cout << "...removeing Spikes";
503 removeSpikes (Ameas, Vcorr);
504 if (verbosityLevel > 2) cout << "...done " << endl;
505
506 // filter Vcorr with sliding average using FIR filter function
507 if (verbosityLevel > 2) cout << "...applying sliding average filter";
508 sliding_avg(Vcorr, Vslide, avg1);
509 if (verbosityLevel > 2) cout << "...done " << endl;
510
511 // filter Vslide with CFD using FIR filter function
512 if (verbosityLevel > 2) cout << "...apllying factfir filter";
513 factfir(b_cfd , a_cfd, k_cfd, Vslide, Vcfd);
514 if (verbosityLevel > 2) cout << "...done " << endl;
515
516 // filter Vcfd with sliding average using FIR filter function
517 if (verbosityLevel > 2) cout << "...applying 2nd sliding average filter";
518 sliding_avg(Vcfd, Vcfd2, avg2);
519 if (verbosityLevel > 2) cout << "...done " << endl;
520
521//-------------------------------------
522// Search vor Zero crossings
523//-------------------------------------
524 if (verbosityLevel > 2) cout << endl << "...searching zero crossings" ;
525 // peaks in Ameas[] are found by searching for zero crossings
526 // in Vcfd2
527 // first Argument 1 means ... *rising* edge
528 // second Argument 1 means ... search with stepsize 1 ... 10 is okay as well
529 vector<Region>* pZXings = zerosearch( Vcfd2 , 1 , 8);
530 // pZXings means "zero cross ings"
531 EnlargeRegion(*pZXings, 10, 10);
532 findAbsMaxInRegions(*pZXings, Vslide);
533 removeMaximaBelow( *pZXings, 3.0);
534 removeRegionWithMaxOnEdge( *pZXings, 2);
535 removeRegionOnFallingEdge( *pZXings, 100);
536 findTimeOfHalfMaxLeft(*pZXings, Vcorr, gBSLMean, 5, 10, verbosityLevel );
537 if (verbosityLevel > 2) cout << "...done" << endl;
538
539//-----------------------------------------------------------------------------
540// Fill Overlay Histos
541//-----------------------------------------------------------------------------
542 if (verbosityLevel > 2)
543 {
544 cout << endl << "Filling Histograms of Pixel# [Set] " << pixelID
545 << " Pixel# [Got] " << pixel[pixelID]->mChid;
546 }
547 FillHistograms(
548 pixel[pixelID],
549 pZXings,
550 AmplWindowWidth,
551 ev,
552 histoOptions,
553 maxPulseOrder,
554 verbosityLevel
555 );
556//-----------------------------------------------------------------------------
557// Spike Debug
558//-----------------------------------------------------------------------------
559 if ( spikeDebug )
560 {
561 ProduceDebugHistos( pZXings );
562 }// end of if(spikeDebug)
563
564 delete pZXings;
565 if (breakout) break;
566//-------------------------------------
567// Draw 1. Set of Pixel Histograms
568//-------------------------------------
569 if ((ev % refresh_rate) == 0)
570 {
571 if (ProduceGraphic)
572 {
573 if (debugPixel)
574 {
575 for ( int order = 0;
576 order < maxPulseOrder;
577 order++
578 )
579 {
580 pixel[pixelID]->DrawHistograms(
581 cgpPixelPulses[order],
582 pulsesCanvasFrameNrs
583 );
584 }
585 if (testmode)
586 {
587 DrawTestHistograms(
588 verbosityLevel
589 );
590 }
591
592 UpdateCanvases(
593 verbosityLevel,
594 maxPulseOrder
595 );
596 }
597 }
598 }
599
600 if (breakout) break;
601
602 if (verbosityLevel > 1)
603 {
604 cout << endl << "...End of Pixel"
605 << endl << "-------------------------------------"<< endl;
606 }
607 }//End of Loop over Sample
608//-------------------------------------
609// end of Loops over Sample
610//-------------------------------------
611
612 }//End of Loop over Events
613//-------------------------------------
614// end of Loops over Events
615//-------------------------------------
616
617
618
619//-------------------------------------
620// Draw Pixel Histogramms of Overlay Spectra
621//-------------------------------------
622
623 for ( int pixelID = firstPixelOfSample;
624 pixelID < firstPixelOfSample + sampleSetSize
625 && pixelID < firstpixel + npixel;
626 pixelID++ )
627 {
628 //here is what happends at the end of each loop over all Events
629// pixel[pixelID]->DrawHistograms(
630// cgpPixelPulses,
631// pulsesCanvasFrameNrs
632// );
633
634// SaveHistograms( //save histograms of generell results into output root file
635// OutRootFileName,
636// CreateSubDirName(pixel[pixelID]->mChid),
637// pixel[pixelID]->hList,
638// saveResults,
639// verbosityLevel
640// );
641
642 if (ProduceGraphic)
643 {
644 UpdateCanvases(
645 verbosityLevel,
646 maxPulseOrder
647 );
648 }
649
650 //Save Histograms of Pixels into Output rootfile
651 if (verbosityLevel > 2)
652 {
653 cout << endl << "Saving Pixel# [Set] " << pixelID
654 << " Pixel# [Got] " << pixel[pixelID]->mChid;
655 }
656 pixel[pixelID]->SavePixelHistograms( OutRootFileName, saveResults );
657
658
659 if (debugPixel)
660 {
661 //Process gui events asynchronously during input
662 cout << endl;
663 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
664 timer.TurnOn();
665 TString input = Getline("Type 'q' to exit, <return> to go on: ");
666 timer.TurnOff();
667 if (input=="q\n") {
668 break;
669 }
670 }
671
672 if (verbosityLevel > 2)
673 {
674 cout << endl << "Deleting Pixel# [Set] " << pixelID
675 << " Pixel# [Got] " << pixel[pixelID]->mChid;
676 }
677 delete pixel[pixelID];
678
679 if (verbosityLevel > 2)
680 {
681 cout << endl << "...End of Pixel-Set"
682 << endl << "------------------------------------------------"
683 << endl;
684 }
685 }
686 }
687 // End of Loop over all Pixels
688
689//-------------------------------------
690// Draw Histograms
691//-------------------------------------
692
693 SaveHistograms( //save histograms of generell results into output root file
694 OutRootFileName,
695 "root",
696 hRootList,
697 saveResults,
698 verbosityLevel
699 );
700
701 if (ProduceGraphic)
702 {
703 UpdateCanvases(
704 verbosityLevel,
705 maxPulseOrder
706 );
707
708 //Process gui events asynchronously during input
709 cout << endl;
710 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
711 timer.TurnOn();
712 Getline("Press <return> to go on");
713 timer.TurnOff();
714
715 DeletePixelCanvases(
716 maxPulseOrder,
717 verbosityLevel
718 );
719 }
720
721 delete datafile;
722 delete[] pixel;
723 return( 0 );
724}
725//----------------------------------------------------------------------------
726// end of main function
727//-----------------------------------------------------------------------------
728
729
730
731
732//-----------------------------------------------------------------------------
733// Funktions
734//-----------------------------------------------------------------------------
735
736/*
737ich erzeuge sowas wie 36 pixel mit ca 20 Histogrammen pro pixel
738wie viel speicher braucht das?
739Ich brauche eine möglichkeit jedes dieser histogramme zu identifizieren
740am besten ein array oder eine liste oder einen abstracten datentyp mit den histogrammen als attribute.
741ne klasse wäre super
742
743*/
744
745void BookTestHistos( int verbosityLevel )
746{
747 if (verbosityLevel > 2) cout << endl << "...book pixel histograms" << endl;
748
749 hTesthisto = new TH1F (
750 "hTesthisto",
751 "Deviation of rising edge and maximum",
752 600,
753 -10.1,
754 10.1
755 );
756 hTesthisto->GetXaxis()->SetTitle( "Timeslices [a.u.]" );
757 hTesthisto->GetYaxis()->SetTitle( "counts" );
758
759 hTesthisto2 = new TH2F (
760 "hTesthisto2",
761 "Deviation of rising edge and maximum by event #",
762 gNEvents,
763 250,
764 800,
765 600,
766 -10.1,
767 10.1
768 );
769// hTesthisto2->GetXaxis()->SetTitle( "Timeslices [a.u.]" );
770// hTesthisto2->GetYaxis()->SetTitle( "counts" );
771// hTesthisto2->SetMarkerStyle( 2 );
772
773 if (verbosityLevel > 2) cout << "...done" << endl;
774}
775//end of BookTestHistos
776//----------------------------------------------------------------------------
777
778
779void
780BookDebugHistos( int verbosityLevel )
781{
782 if (verbosityLevel > 2) cout << endl << "...book histograms" << endl;
783 hRootList = new TObjArray;
784 debugHistos = new TH1F[ Number_Of_Debug_Histo_Types ];
785 for ( int type = 0; type < Number_Of_Debug_Histo_Types; type++){
786 debugHistos[ type ].SetBins(1024, 0, 1024);
787 debugHistos[ type ].SetLineColor(1);
788 debugHistos[ type ].SetLineWidth(2);
789 debugHistos[ type ].SetStats(false);
790
791 // set X axis paras
792 debugHistos[ type ].GetXaxis()->SetLabelSize(0.08);
793 debugHistos[ type ].GetXaxis()->SetTitleSize(0.08);
794 debugHistos[ type ].GetXaxis()->SetTitleOffset(1.0);
795 debugHistos[ type ].GetXaxis()->SetTitle(Form("Time slice [%.1f ns/slice]", 1./2.));
796
797 // set Y axis paras
798 debugHistos[ type ].GetYaxis()->SetLabelSize(0.08);
799 debugHistos[ type ].GetYaxis()->SetTitleSize(0.08);
800 debugHistos[ type ].GetYaxis()->SetTitleOffset(0.3);
801 debugHistos[ type ].GetYaxis()->SetTitle("Amplitude [mV]");
802 }
803 if (verbosityLevel > 2) cout << "...done" << endl;
804}
805// end of BookDebugHistos
806//----------------------------------------------------------------------------
807
808
809void
810FillHistograms(
811 Pixel* CurrentPixel,
812 vector<Region>* pZXings,
813 int AmplWindowWidth,
814 int eventNumber,
815 TString histoOptions,
816 int maxPulseOrder,
817 int verbosityLevel
818 )
819{
820 if (verbosityLevel > 2) cout << endl << "...filling pulse histograms" ;
821 bool use_this_peak=false;
822 int order_of_pulse=0;
823 vector<Region>::iterator reg;
824 //Loop over all found zerocrossings in Region
825 for (reg = pZXings->begin() ; reg < pZXings->end() ; reg++)
826 {
827 //skip those who are at the Rim of the ROI
828 if (reg->maxPos < 12 || (unsigned int) reg->maxPos > gRegionOfInterest-12)
829 {
830 if (verbosityLevel > 2) cout << endl << "\t...out of range" << endl;
831 continue;
832 }
833 // Define axis range of Histogram
834 int Left = reg->maxPos - gPixelOverlayXaxisLeft;
835 int Right = reg->maxPos + gPixelOverlayXaxisRight;
836 int EdgeLeft = reg->halfRisingEdgePos - gPixelOverlayXaxisLeft;
837 int EdgeRight = reg->halfRisingEdgePos + gPixelOverlayXaxisRight;
838
839 if (Left < 0) Left = 0;
840 if (EdgeLeft < 0) EdgeLeft = 0;
841 if (Right > (int)Ameas.size() ) Right = Ameas.size();
842 if (EdgeRight > (int)Ameas.size() ) EdgeRight = Ameas.size();
843
844
845// hTesthisto2->Fill( eventNumber, reg->slopeOfRisingEdge ) ;
846
847 if (verbosityLevel > 2) cout << endl << "\t...choosing Histogram" ;
848
849 //determine order of pulse and which histogram shall be filled
850 for(int order = 0; order < maxPulseOrder; order++)
851 {
852 if (Ameas[reg->maxPos] >= ((gGainMean*(order+1)) - (AmplWindowWidth/2))
853 && Ameas[reg->maxPos] < ((gGainMean*(order+1)) + AmplWindowWidth/2))
854 {
855 if (verbosityLevel > 2) cout << "...#" << order ;
856 order_of_pulse = order;
857 use_this_peak = true;
858 break;
859 }
860 else if (order >= (maxPulseOrder - 1))
861 {
862 use_this_peak = false;
863 if (verbosityLevel > 2) cout << "...NONE" << endl ;
864 }
865 if (use_this_peak)
866 {
867 if ( histoOptions.Contains("S") )
868 {
869// if (verbosityLevel > 2)
870 cout << endl << "\t...filling Edge Histogram" ;
871 CurrentPixel->hSlopeRisingEdge[order]->Fill( reg->slopeOfRisingEdge ) ;
872 }
873 }
874 }
875 //Fill overlay und profile histograms
876 if (use_this_peak){
877 for ( int pos = Left; pos < Right; pos++){
878// if ();
879 CurrentPixel->hMaxOverlay[order_of_pulse]->Fill( pos - (reg->maxPos), Ameas[pos]) ;
880 CurrentPixel->hMaxProfile[order_of_pulse]->Fill( pos - (reg->maxPos), Ameas[pos]) ;
881 }
882 for ( int pos = EdgeLeft; pos < EdgeRight; pos++){
883 CurrentPixel->hEdgeOverlay[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), Ameas[pos]) ;
884 CurrentPixel->hEdgeProfile[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), Ameas[pos]) ;
885
886 }
887 if (verbosityLevel > 2) cout << "...done" << endl;
888 }
889 //Breakout option
890 if (breakout)
891 break;
892 }
893 if (verbosityLevel > 2) cout << "...done" << endl;
894}
895// end of FillHistograms
896//----------------------------------------------------------------------------
897
898void
899DrawTestHistograms(
900 int verbosityLevel
901 )
902{
903 if (verbosityLevel > 2) cout << endl << "...drawing pulse histograms TODO CANVAS" ;
904
905// cgpTestHistos->cd(1);
906// hTesthisto->Draw();
907// cgpTestHistos->cd(2);
908// hTesthisto2->Draw("BOX");
909}
910// end of DrawTestHistograms
911//----------------------------------------------------------------------------
912
913void
914UpdateCanvases(
915 int verbosityLevel,
916 int max_pulse_order
917 )
918{
919 if (verbosityLevel > 3) cout << endl << "...updating canvases" ;
920 for (int pulse_order = 0; pulse_order < max_pulse_order; pulse_order++)
921 {
922 cgpPixelPulses[pulse_order]->Modified();
923 cgpPixelPulses[pulse_order]->Update();
924 cgpTestHistos->Modified();
925 cgpTestHistos->Update();
926 }
927}
928// end of UpdateCanvases
929//----------------------------------------------------------------------------
930
931void
932DeletePixelCanvases(
933 int maxPulseOrder,
934 int verbosityLevel
935 )
936{
937 if (verbosityLevel > 2) cout << endl << "...delete pixel Canvases" ;
938 for (int pulse_order = 0; pulse_order < maxPulseOrder; pulse_order++ )
939 {
940 delete cgpPixelPulses[pulse_order];
941 delete cgpDistributions[pulse_order];
942 }
943 delete[] cgpPixelPulses;
944 delete[] cgpDistributions;
945}
946
947//----------------------------------------------------------------------------
948//----------------------------------------------------------------------------
949
950void
951ProduceDebugHistos(
952 vector<Region>* pZXings
953 )
954{
955 // TODO do this correct. The vectors should be the rigt ones... this is just luck
956 debugHistos[Ameas_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
957 debugHistos[Vcorr_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
958 debugHistos[Vslide_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
959 debugHistos[Vcfd_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
960 debugHistos[Vcfd2_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5);
961
962 for ( unsigned int sl = 0; sl < gRegionOfInterest; sl++)
963 {
964 debugHistos[Ameas_].SetBinContent(sl, Ameas[sl]);
965 debugHistos[Vcorr_].SetBinContent(sl, Vcorr[sl]);
966 debugHistos[Vslide_].SetBinContent( sl, Vslide[sl] );
967 debugHistos[Vcfd_].SetBinContent( sl, Vcfd[sl] );
968 debugHistos[Vcfd2_].SetBinContent( sl, Vcfd2[sl] );
969 }
970
971
972 cFiltered->cd(1);
973 gPad->SetGrid();
974 debugHistos[Ameas_].Draw();
975
976 cFiltered->cd(2);
977 gPad->SetGrid();
978 debugHistos[Vcorr_].Draw();
979
980 cFiltered->cd(3);
981 gPad->SetGrid();
982 debugHistos[Vslide_].Draw();
983
984 TBox *OneBox;
985 vector<TBox*> MyBoxes;
986 for (unsigned int i=0; i<pZXings->size(); i++){
987 OneBox = new TBox(
988 pZXings->at(i).maxPos -10 ,
989 pZXings->at(i).maxVal -0.5,
990 pZXings->at(i).maxPos +10 ,
991 pZXings->at(i).maxVal +0.5);
992 OneBox->SetLineColor(kBlue);
993 OneBox->SetLineWidth(1);
994 OneBox->SetFillStyle(0);
995 OneBox->SetFillColor(kRed);
996 MyBoxes.push_back(OneBox);
997 OneBox->Draw();
998 }
999
1000// cFiltered->cd(3);
1001// gPad->SetGrid();
1002// debugHistos[Vcfd2_].Draw();
1003// TLine *zeroline = new TLine(0, 0, 1024, 0);
1004// zeroline->SetLineColor(kBlue);
1005// zeroline->Draw();
1006
1007 cFiltered->Update();
1008
1009
1010 //Process gui events asynchronously during input
1011 TTimer timer("gSystem->ProcessEvents();", 50, kFALSE);
1012 timer.TurnOn();
1013 TString input = Getline("Type 'q' to exit, <return> to go on: ");
1014 timer.TurnOff();
1015 if (input=="q\n") {
1016 breakout=true;
1017 }
1018
1019 //TODO!!!!!!!!!
1020 // do some Garbage collection ...
1021 // all the Objects on the heap should be deleted here.
1022
1023}// end of if(spikeDebug)
1024
1025
1026int main()
1027{
1028
1029FPulseOverlay();
1030return 0;
1031
1032}
Note: See TracBrowser for help on using the repository browser.