Index: /fact/tools/rootmacros/FPulseTemplate.C =================================================================== --- /fact/tools/rootmacros/FPulseTemplate.C (revision 12742) +++ /fact/tools/rootmacros/FPulseTemplate.C (revision 12742) @@ -0,0 +1,679 @@ +/********************** FPulseTemplate *********************** +* read FACT raw data +* remove spikes +* calculate baseline +* find peaks (CFD and normal discriminator) +* compute pulse height and pulse integral spektrum of the peaks ++ search for Peaks in data ++ put peaks into different histos depending und Amplitude ++ draw histos for single, double, tripple, ....above 100mV Photonpulses ++ draw Parameterdevelopment depending on photon quantity ++ form a template (shape depending on number of photons) + so it can be used for detection of other peaks +*************************************************************/ + +//----------------------------------------------------------------------------- +// root libraries +//----------------------------------------------------------------------------- + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include +#include + +#define NPIX 1440 +#define NCELL 1024 +#define FAD_MAX_SAMPLES 1024 + +#define HAVE_ZLIB + +//----------------------------------------------------------------------------- +// rootmacros +//----------------------------------------------------------------------------- + +#include "fits.h" + +#include "openFits.h" +#include "openFits.c" + +#include "discriminator.h" +#include "discriminator.C" +#include "zerosearch.h" +#include "zerosearch.C" +#include "factfir.C" + +#include "DrsCalibration.C" +#include "DrsCalibration.h" + +#include "SpikeRemoval.h" +#include "SpikeRemoval.C" + +//----------------------------------------------------------------------------- +// Decleration of global variables +//----------------------------------------------------------------------------- + +bool breakout=false; + +int NEvents; +vector AllPixelDataVector; +vector StartCellVector; +unsigned int CurrentEventID; +size_t PXLxROI; +UInt_t RegionOfInterest; +UInt_t NumberOfPixels; + +size_t TriggerOffsetROI, RC; +size_t drs_n; +vector drs_basemean; +vector drs_gainmean; +vector drs_triggeroffsetmean; + +vector Ameas(FAD_MAX_SAMPLES); // copy of the data (measured amplitude +vector Vcorr(FAD_MAX_SAMPLES); // corrected Values +vector Vslide(FAD_MAX_SAMPLES); // sliding average result +vector Vcfd(FAD_MAX_SAMPLES); // CDF result +vector Vcfd2(FAD_MAX_SAMPLES); // CDF result + 2nd sliding average + +typedef struct{ + float maxAmpl; + float countOfMax; + } OverlayMaximum; + +vector SingleMaximum; +// histograms +const int NumberOfDebugHistoTypes = 7; + +const unsigned int + Ameas_ = 0, + N1mean_ = 1, + Vcorr_ = 2, + Vtest_ = 3, + Vslide_ = 4, + Vcfd_ = 5, + Vcfd2_ = 6; + +//----------------------------------------------------------------------------- +// Initialisation of histograms +//----------------------------------------------------------------------------- + +TH1F* h; +TH1F *debugHistos; +TH2F *hOverlayTemp = NULL; +TH2F *hSinglePeakOverlay;//histogrammm for overlay of detected Peaks +TH2F *hDoublePeakOverlay; +TH2F *hTripplePeakOverlay; +TH2F *hLargePeakOverlay; + +TH1F *hMaximumTemp = NULL; +TH1F *hSinglePeakMaximum; +TH1F *hDoublePeakMaximum; +TH1F *hTripplePeakMaximum; +TH1F *hLargePeakMaximum; + +int hPeakOverlayXaxisLeft,hPeakOverlayXaxisRight; + +TObjArray hList; +TObjArray hListBaseline; + +//----------------------------------------------------------------------------- +// Functions +//----------------------------------------------------------------------------- + +void BookHistos( ); +void SaveHistograms( const char * ); + +//----------------------------------------------------------------------------- +//----------------------------------------------------------------------------- +// MAIN - Funtion +//----------------------------------------------------------------------------- +//----------------------------------------------------------------------------- +int FPulseTemplate( + char *datafilename = "../data/2011/11/09/20111109_006.fits.gz", + const char *drsfilename = "../data/2011/11/09/20111109_003.drs.fits.gz", + const char *OutRootFileName = "../analysis/20111109_006.fits.gz.peaktypes.root", + int firstevent = 0, + int nevents = -1, + int firstpixel = 0, + int npixel = -1, + bool spikeDebug = false, + int avg1 = 14, + int avg2 = 8, + int OverlayWindowLeft = 50, + int OverlayWindowRight = 250, + int verbosityLevel = 1, // different verbosity levels can be implemented here + bool ProduceGraphic = true + ) + +{ +//----------------------------------------------------------------------------- +// histogramm Settings +//----------------------------------------------------------------------------- + + hPeakOverlayXaxisLeft = OverlayWindowLeft; + hPeakOverlayXaxisRight = OverlayWindowRight; + + gStyle->SetPalette(1,0); + gROOT->SetStyle("Plain"); + +// Create (pointer to) Canvases, which are used in every run, +// also in 'non-debug' runs + // Canvases only need if spike Debug, but I want to deklare + // the pointers anyway ... + TCanvas *cFiltered = NULL; + if (spikeDebug){ + cFiltered = new TCanvas("cFiltered","filtered DRS Waveforms", 1,310,400,300); + cFiltered->Divide(1, 3); + } + // Canvases to show the peak template + TCanvas *cPulsetypes = NULL; + cPulsetypes = new TCanvas("cPulsetypes", "Overlay of detected Peaks", 1, 1, 1400, 700); + cPulsetypes->Divide(4,2); + +//----------------------------------------------------------------------------- +// Filter-Settings +//----------------------------------------------------------------------------- +// CFD filter settings + int k_cfd = 10; + vector a_cfd(k_cfd, 0); + double b_cfd = 1.; + a_cfd[0]=-0.75; + a_cfd[k_cfd-1]=1.; + +//----------------------------------------------------------------------------- +// prepare datafile +//----------------------------------------------------------------------------- + +// Open the data file + + fits * datafile; + // Opens the raw data file and 'binds' the variables given as + // Parameters to the data file. So they are filled with + // raw data as soon as datafile->GetRow(int) is called. + NEvents = openDataFits( datafilename, &datafile, + AllPixelDataVector, StartCellVector, CurrentEventID, + RegionOfInterest, NumberOfPixels, PXLxROI, verbosityLevel); + if (NEvents == 0){ + cout << "return code of OpenDataFile:" << datafilename<< endl; + cout << "is zero -> aborting." << endl; + return 1; + } + + if (verbosityLevel > 0) + cout <<"number of events in file: "<< NEvents << "\t"; + if ( nevents == -1 || nevents > NEvents ) nevents = NEvents; // -1 means all! + if (verbosityLevel > 0) + cout <<"of, which "<< nevents << "will be processed"<< endl; + + if (verbosityLevel > 0) + cout <<"Total # of Pixel: "<< NumberOfPixels << "\t"; + if ( npixel == -1 || npixel > (int)NumberOfPixels ) npixel = (int)NumberOfPixels; // -1 means all! + if (verbosityLevel > 0) + cout <<"of, which "<< npixel << "will be processed"<< endl; +//Get the DRS calibration + RC = openCalibFits( drsfilename, drs_basemean, drs_gainmean, drs_triggeroffsetmean, TriggerOffsetROI); + if (RC == 0){ + cout << "return code of openCalibFits:" << drsfilename << endl; + cout << "is zero -> aborting." << endl; + return 1; + } +// Book the histograms + BookHistos( ); + +//----------------------------------------------------------------------------- +// Loops over Every Event and Pixel +//----------------------------------------------------------------------------- + for ( int ev = firstevent; ev < firstevent + nevents; ev++) { + // Get an Event --> consists of 1440 Pixel ...erm....data + datafile->GetRow( ev ); + +//------------------------------------- +// Loop over every Pixel of Event +//------------------------------------- + + for ( int pix = firstpixel; pix < firstpixel + npixel; pix++ ){ + if (verbosityLevel > 0){ + if (pix == firstpixel){ + cout << "Processing Event: " << CurrentEventID << "/" << nevents << endl; + } + } + + applyDrsCalibration( Ameas,pix,12,12, + drs_basemean, drs_gainmean, drs_triggeroffsetmean, + RegionOfInterest, AllPixelDataVector, StartCellVector); + + // finds spikes in the raw data, and interpolates the value + // spikes are: 1 or 2 slice wide, positive non physical artifacts + removeSpikes (Ameas, Vcorr); + + // filter Vcorr with sliding average using FIR filter function + sliding_avg(Vcorr, Vslide, avg1); + // filter Vslide with CFD using FIR filter function + factfir(b_cfd , a_cfd, k_cfd, Vslide, Vcfd); + // filter Vcfd with sliding average using FIR filter function + sliding_avg(Vcfd, Vcfd2, avg2); + + // peaks in Ameas[] are found by searching for zero crossings + // in Vcfd2 + // first Argument 1 means ... *rising* edge + // second Argument 1 means ... search with stepsize 1 ... 10 is okay as well + vector * zXings = zerosearch( Vcfd2 , 1 , 8); + // zXings means "zero cross ings" + EnlargeRegion(*zXings, 10, 10); + findAbsMaxInRegions(*zXings, Vslide); + removeMaximaBelow( *zXings, 3.0); + removeRegionWithMaxOnEdge( *zXings, 2); + removeRegionOnFallingEdge( *zXings, 100); + +//----------------------------------------------------------------------------- +// Fill Overlay Histos +//----------------------------------------------------------------------------- + vector::iterator it; + for (it = zXings->begin() ; it < zXings->end() ; it++){ + if (it->maxPos < 12 || it->maxPos > RegionOfInterest-12) + continue; + //domstest->Fill(it->maxPos); + int Left = it->maxPos - OverlayWindowLeft; + int Right = it->maxPos + OverlayWindowRight; + if (Left < 0) + Left =0; + if (Right > (int)Vcorr.size() ) + Right=Vcorr.size() ; + if (Vslide[it->maxPos] >= 5 && Vslide[it->maxPos] < 13) hOverlayTemp = &*hSinglePeakOverlay; + else if (Vslide[it->maxPos] >= 13 && Vslide[it->maxPos] < 23) hOverlayTemp = hDoublePeakOverlay; + else if (Vslide[it->maxPos] >= 23 && Vslide[it->maxPos] < 33) hOverlayTemp = hTripplePeakOverlay; + else if (Vslide[it->maxPos] >= 33) hOverlayTemp = hLargePeakOverlay; + + for ( int pos = Left; pos < Right; pos++){ +// if (Vslide[it->maxPos] >= 5 && Vslide[it->maxPos] < 15) hSinglePeakOverlay->Fill( pos - (it->maxPos), Vslide[pos]); +// else if (Vslide[it->maxPos] >= 15 && Vslide[it->maxPos] < 25) hDoublePeakOverlay->Fill( pos - (it->maxPos), Vslide[pos]); +// else if (Vslide[it->maxPos] >= 25 && Vslide[it->maxPos] < 35) hTripplePeakOverlay->Fill( pos - (it->maxPos), Vslide[pos]); +// else if (Vslide[it->maxPos] >= 35) hOverlayTemp = hLargePeakOverlay; + hOverlayTemp->Fill( pos - (it->maxPos), Vslide[pos]) ; + } + } + +//----------------------------------------------------------------------------- +// Spike Debug +//----------------------------------------------------------------------------- + if ( spikeDebug ){ + // TODO do this correct. The vectors should be the rigt ones... this is just luck + debugHistos[Ameas_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5); + debugHistos[Vcorr_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5); + debugHistos[Vslide_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5); + debugHistos[Vcfd_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5); + debugHistos[Vcfd2_].GetXaxis()->Set(Ameas.size() , -0.5 , Ameas.size()-0.5); + + for ( unsigned int sl = 0; sl < RegionOfInterest; sl++){ + debugHistos[Ameas_].SetBinContent(sl, Ameas[sl]); + debugHistos[Vcorr_].SetBinContent(sl, Vcorr[sl]); + debugHistos[Vslide_].SetBinContent( sl, Vslide[sl] ); + debugHistos[Vcfd_].SetBinContent( sl, Vcfd[sl] ); + debugHistos[Vcfd2_].SetBinContent( sl, Vcfd2[sl] ); + } + + + cFiltered->cd(1); + gPad->SetGrid(); + debugHistos[Vcorr_].Draw(); + + cFiltered->cd(2); + gPad->SetGrid(); + debugHistos[Vslide_].Draw(); + + TBox *OneBox; + vector MyBoxes; + for (unsigned int i=0; isize(); i++){ + OneBox = new TBox( + zXings->at(i).maxPos -10 , + zXings->at(i).maxVal -0.5, + zXings->at(i).maxPos +10 , + zXings->at(i).maxVal +0.5); + OneBox->SetLineColor(kBlue); + OneBox->SetLineWidth(1); + OneBox->SetFillStyle(0); + OneBox->SetFillColor(kRed); + MyBoxes.push_back(OneBox); + OneBox->Draw(); + } + + cFiltered->cd(3); + gPad->SetGrid(); + debugHistos[Vcfd2_].Draw(); + TLine *zeroline = new TLine(0, 0, 1024, 0); + zeroline->SetLineColor(kBlue); + zeroline->Draw(); + + cFiltered->Update(); + + //OVERLAY PEAKS + cPulsetypes->cd(1); + gPad->SetGrid(); + hSinglePeakOverlay->Draw("COLZ"); + + cPulsetypes->cd(2); + gPad->SetGrid(); + hDoublePeakOverlay->Draw("COLZ"); + + cPulsetypes->cd(3); + gPad->SetGrid(); + hTripplePeakOverlay->Draw("COLZ"); + + cPulsetypes->cd(4); + gPad->SetGrid(); + hLargePeakOverlay->Draw("COLZ"); + + cPulsetypes->Modified(); + cPulsetypes->Update(); + + //Process gui events asynchronously during input + TTimer timer("gSystem->ProcessEvents();", 50, kFALSE); + timer.TurnOn(); + TString input = Getline("Type 'q' to exit, to go on: "); + timer.TurnOff(); + if (input=="q\n") { + breakout=true; + } + + //TODO!!!!!!!!! + // do some Garbage collection ... + // all the Objects on the heap should be deleted here. + + }// end of if(spikeDebug) + + delete zXings; + if (breakout) + break; + + } + //------------------------------------- + // end of loop over pixels + //------------------------------------- + if (ProduceGraphic){ + if (ev % 500 == 0){ + + //OVERLAY PEAKS + cPulsetypes->cd(1); + gPad->SetGrid(); + hSinglePeakOverlay->Draw("COLZ"); + + cPulsetypes->cd(2); + gPad->SetGrid(); + hDoublePeakOverlay->Draw("COLZ"); + + cPulsetypes->cd(3); + gPad->SetGrid(); + hTripplePeakOverlay->Draw("COLZ"); + + cPulsetypes->cd(4); + gPad->SetGrid(); + hLargePeakOverlay->Draw("COLZ"); + // cPulsetypes->cd(5); + // hSinglePeakMaximum->Draw(""); + cPulsetypes->Modified(); + cPulsetypes->Update(); + } + } + if (breakout) + break; + } + +//------------------------------------- +// end of loop over Events +//------------------------------------- + +//------------------------------------- +// Histogramm of Maximas in Overlay Spectra +//------------------------------------- +cout << "producing...Maximum peaks" << endl; + for (unsigned int cnt = 1 ; cnt <= 4 ; cnt ++){ + cout << "peak type: " << cnt << endl; + if (cnt == 1){ + hMaximumTemp = hSinglePeakMaximum; + hOverlayTemp = hSinglePeakOverlay; + } + else if (cnt == 2){ + hMaximumTemp = hDoublePeakMaximum; + hOverlayTemp = hDoublePeakOverlay; + } + else if (cnt == 3){ + hMaximumTemp = hTripplePeakMaximum; + hOverlayTemp = hTripplePeakOverlay; + } + else if (cnt == 4){ + hMaximumTemp = hLargePeakMaximum; + hOverlayTemp = hLargePeakOverlay; + } + // resize vector SingleMaximum to number of timeslices in Overlay Spectra + SingleMaximum.resize(hOverlayTemp->GetNbinsX()); + //put maximumvalue of every Y-projection of every timeslice into vector + for (int TimeSlice = 0; TimeSlice <= hOverlayTemp->GetNbinsX(); TimeSlice++ ){ + TString histotitle; + histotitle += "hproj_py"; + histotitle += cnt ; + TH1D* hProjPeak = hOverlayTemp->ProjectionY(histotitle, TimeSlice, TimeSlice, ""); + SingleMaximum[ TimeSlice ].maxAmpl = (hProjPeak->GetMaximumBin() * 0.5) - 3.5; + SingleMaximum[ TimeSlice ].countOfMax = hProjPeak->GetBinContent( hProjPeak->GetMaximumBin() ); + hMaximumTemp->SetBinContent(TimeSlice, SingleMaximum[ TimeSlice ].maxAmpl ); + } + hMaximumTemp->Fit("landau", "", "", -50, 250); + } +//------------------------------------- +// Draw Histograms +//------------------------------------- + if (ProduceGraphic){ + //OVERLAY PEAKS +cout << "producing...So" << endl; + cPulsetypes->cd(1); + gPad->SetGrid(); + hSinglePeakOverlay->ResetStats(); + hSinglePeakOverlay->Draw("COLZ"); +cout << "producing...DO" << endl; + cPulsetypes->cd(2); + gPad->SetGrid(); + hDoublePeakOverlay->ResetStats(); + hDoublePeakOverlay->Draw("COLZ"); +cout << "producing...TO" << endl; + cPulsetypes->cd(3); + gPad->SetGrid(); + hTripplePeakOverlay->ResetStats(); + hTripplePeakOverlay->Draw("COLZ"); +cout << "producing...LO" << endl; + cPulsetypes->cd(4); + gPad->SetGrid(); + hLargePeakOverlay->ResetStats(); + hLargePeakOverlay->Draw("COLZ"); +cout << "producing...SM" << endl; + cPulsetypes->cd(5); + hSinglePeakMaximum->ResetStats(); + hSinglePeakMaximum->Draw(""); +cout << "producing...DM" << endl; + cPulsetypes->cd(6); + hDoublePeakMaximum->ResetStats(); + hDoublePeakMaximum->Draw(""); +cout << "producing...TM" << endl; + cPulsetypes->cd(7); + hTripplePeakMaximum->ResetStats(); + hTripplePeakMaximum->Draw(""); +cout << "producing...LM" << endl; + cPulsetypes->cd(8); + hLargePeakMaximum->ResetStats(); + hLargePeakMaximum->Draw(""); +cout << "producing...Done" << endl; + + + + + cPulsetypes->Modified(); + cPulsetypes->Update(); + +// cPixelPeakOverlay->cd(); +// hPixelPeakOverlay->Draw("COLZ"); +// cPixelPeakOverlay->Modified(); +// cPixelPeakOverlay->Update(); + + } + SaveHistograms( OutRootFileName ); + + return( 0 ); +} + +// booking and parameter settings for all histos +void BookHistos( ){ + + debugHistos = new TH1F[ NumberOfDebugHistoTypes ]; + for ( int type = 0; type < NumberOfDebugHistoTypes; type++){ + debugHistos[ type ].SetBins(1024, 0, 1024); + debugHistos[ type ].SetLineColor(1); + debugHistos[ type ].SetLineWidth(2); + + // set X axis paras + debugHistos[ type ].GetXaxis()->SetLabelSize(0.1); + debugHistos[ type ].GetXaxis()->SetTitleSize(0.1); + debugHistos[ type ].GetXaxis()->SetTitleOffset(1.2); + debugHistos[ type ].GetXaxis()->SetTitle(Form("Time slice (%.1f ns/slice)", 1./2.)); + + // set Y axis paras + debugHistos[ type ].GetYaxis()->SetLabelSize(0.1); + debugHistos[ type ].GetYaxis()->SetTitleSize(0.1); + debugHistos[ type ].GetYaxis()->SetTitleOffset(0.3); + debugHistos[ type ].GetYaxis()->SetTitle("Amplitude (a.u.)"); + } +/* + ProjHistos = new TH1D[ NumberOfTimeSlices ]; + for ( int type = 0; type < NumberOfTimeSlices; type++){ + ProjHistos[ type ].SetBins(1024, 0, 1024); + ProjHistos[ type ].SetLineColor(1); + ProjHistos[ type ].SetLineWidth(2); + + // set X axis paras + ProjHistos[ type ].GetXaxis()->SetLabelSize(0.1); + ProjHistos[ type ].GetXaxis()->SetTitleSize(0.1); + ProjHistos[ type ].GetXaxis()->SetTitleOffset(1.2); + ProjHistos[ type ].GetXaxis()->SetTitle(Form("Time slice (%.1f ns/slice)", 1./2.)); + + // set Y axis paras + ProjHistos[ type ].GetYaxis()->SetLabelSize(0.1); + ProjHistos[ type ].GetYaxis()->SetTitleSize(0.1); + ProjHistos[ type ].GetYaxis()->SetTitleOffset(0.3); + ProjHistos[ type ].GetYaxis()->SetTitle("Amplitude (a.u.)"); + } +*/ +//----------------------------------------------------------------------------- + hSinglePeakOverlay = new TH2F("hSinglePeakOverlay", "Overlay of detected Single Photon Peaks", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 , + 512, -55.5, 300.5); + hSinglePeakOverlay->SetAxisRange(-5.5, 35.5, "Y"); + hSinglePeakOverlay->GetXaxis()->SetTitle( "Timeslices" ); + hSinglePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + //hSinglePeakOverlay->SetBit(TH2F::kCanRebin); + hList.Add( hSinglePeakOverlay ); + + hSinglePeakMaximum = new TH1F("hSinglePeakMaximum", "Single Peak derived with Maximum of above Spektrum", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , + (-1*hPeakOverlayXaxisLeft)-0.5, + hPeakOverlayXaxisRight-0.5 + ); + hSinglePeakMaximum->SetAxisRange(-0.5, 25.5, "Y"); + hSinglePeakMaximum->GetXaxis()->SetTitle( "Timeslices" ); + hSinglePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + hList.Add( hSinglePeakMaximum ); +//----------------------------------------------------------------------------- + + hDoublePeakOverlay = new TH2F("hDoublePeakOverlay", "Overlay of detected Double Photon Peaks", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 , + 512, -55.5, 300.5); + hDoublePeakOverlay->SetAxisRange(-5.5, 35.5, "Y"); + hDoublePeakOverlay->GetXaxis()->SetTitle( "Timeslices" ); + hDoublePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + //hSinglePeakOverlay->SetBit(TH2F::kCanRebin); + hList.Add( hDoublePeakOverlay );; + + hDoublePeakMaximum = new TH1F("hSinglePeakMaximum", "Double Peak derived with Maximum of above Spektrum", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , + (-1*hPeakOverlayXaxisLeft)-0.5, + hPeakOverlayXaxisRight-0.5 + ); + hDoublePeakMaximum->SetAxisRange(-0.5, 25.5, "Y"); + hDoublePeakMaximum->GetXaxis()->SetTitle( "Timeslices" ); + hDoublePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + hList.Add( hDoublePeakMaximum ); + //----------------------------------------------------------------------------- + + hTripplePeakOverlay= new TH2F("hTripplePeakOverlay", "Overlay of detected Tripple Photon Peaks", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 , + 512, -55.5, 300.5); + hTripplePeakOverlay->SetAxisRange(-5.5, 35.5, "Y"); + hTripplePeakOverlay->GetXaxis()->SetTitle( "Timeslices" ); + hTripplePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + //hTripplePeakOverlay->SetBit(TH2F::kCanRebin); + hList.Add( hTripplePeakOverlay );; + + hTripplePeakMaximum = new TH1F("hSinglePeakMaximum", "Triple Peak derived with Maximum of above Spektrum", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , + (-1*hPeakOverlayXaxisLeft)-0.5, + hPeakOverlayXaxisRight-0.5 + ); + hTripplePeakMaximum->SetAxisRange(-0.5, 25.5, "Y"); + hTripplePeakMaximum->GetXaxis()->SetTitle( "Timeslices" ); + hTripplePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + hList.Add( hTripplePeakMaximum ); + //----------------------------------------------------------------------------- + + hLargePeakOverlay= new TH2F("hLargePeakOverlay", "Overlay of detected Multi Photon Peaks", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 , + 512, -5.5, 300.5 ); + hLargePeakOverlay->SetAxisRange(-5.5, 200.5, "Y"); + hLargePeakOverlay->GetXaxis()->SetTitle( "Timeslices" ); + hLargePeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + //hLargePeakOverlay->SetBit(TH2F::kCanRebin); + hList.Add( hLargePeakOverlay );; + + hLargePeakMaximum = new TH1F("hLargePeakMaximum", "Peak derived with Maximum of above Spektrum", + hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , + (-1*hPeakOverlayXaxisLeft)-0.5, + hPeakOverlayXaxisRight-0.5 + ); + hLargePeakMaximum->SetAxisRange(-0.5, 50.5, "Y"); + hLargePeakMaximum->GetXaxis()->SetTitle( "Timeslices" ); + hLargePeakMaximum->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); + hList.Add( hLargePeakMaximum ); + //----------------------------------------------------------------------------- + +// hPixelPeakOverlay = new TH2F("hPixelPeakOverlay", "Maximum of Statistic of overlayed Peaks", +// hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 , +// 512, -55.5, 200.5 ); +// hPixelPeakOverlay->GetXaxis()->SetTitle( "Timeslices" ); +// hPixelPeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); +// hList.Add( hPixelPeakOverlay ); + +// hEventPeakOverlay = new TH2F("hEventPeakOverlay", "Overlay of detected Peaks of all Pixel of one Event", +// hPeakOverlayXaxisLeft + hPeakOverlayXaxisRight , (-1*hPeakOverlayXaxisLeft)-0.5, hPeakOverlayXaxisRight-0.5 , +// 4096, -48.5, 200.5 ); +// hEventPeakOverlay->GetXaxis()->SetTitle( "Timeslices" ); +// hEventPeakOverlay->GetYaxis()->SetTitle( "Amplitude [a.u.]" ); +// hList.Add( hEventPeakOverlay ); + +} + +void SaveHistograms(const char * loc_fname ){ + // create the histogram file (replace if already existing) + TFile tf( loc_fname, "RECREATE"); + + hList.Write(); // write the major histograms into the top level directory + + tf.Close(); // close the file +} // end of SaveHistograms( char * loc_fname ) +