/********************** FPulseOverlay *********************** * read FACT raw data * remove spikes * calculate baseline * find peaks (CFD and normal discriminator) + search for Peaks in data + put peaks into different histos depending on Amplitude + draw histos for single, double, tripple, ....Photonpulses + draw Parameterdevelopment 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 #include #include #include #include #include using namespace std; #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.h" //#include "../factfir.C" //#include "../DrsCalibration.C" #include "../DrsCalibration.h" #include "../SpikeRemoval.h" //#include "../SpikeRemoval.C" #include "rootfilehandler.h" //#include "rootfilehandler.C" #include "pixel.h" //#include "pixel.C" //---------------------------------------------------------------------------- // Decleration of global variables //---------------------------------------------------------------------------- bool breakout=false; vector AllPixelDataVector; vector StartCellVector; unsigned int CurrentEventID; size_t PXLxROI; UInt_t gRegionOfInterest; UInt_t NumberOfPixels; TString histotitle; 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 int gNEvents=0; float gGainMean = 9; float gBSLMean = 0; //typedef struct{ // double maxAmpl; // int countOfMax; // } OverlayMaximum; //typedef struct{ // TString name; // TString title; // TString xTitle; // TString yTitle; // float yMax; // float yMin; //} histAttrib_t; //histAttrib_t* gHistoAttrib[maxPulseOrder]; //histAttrib_t* gProfileAttrib[maxPulseOrder]; // histograms const int Number_Of_Debug_Histo_Types = 7; const unsigned int Ameas_ = 0, N1mean_ = 1, Vcorr_ = 2, Vtest_ = 3, Vslide_ = 4, Vcfd_ = 5, Vcfd2_ = 6; //const char* gHistoTypes[8] = { // "hMaxOverlay", // "hPixelMax", // "hEdgeOverlay", // "hMaxProfile", // "hAllPixelOverlay", // "hAllPixelMax", // "hAllPixelMaxGaus", // "hAllPixelProfile" //}; //---------------------------------------------------------------------------- // Initialisation of histograms //---------------------------------------------------------------------------- // Temporary Objects TH1F* debugHistos = NULL; //TH2F* hOverlayTemp = NULL; //TH1D* hProjPeak = NULL; TH1F* hTesthisto = NULL; TH2F* hTesthisto2 = NULL; //Histogram Parameters Int_t gPixelOverlayXaxisLeft = 0; Int_t gPixelOverlayXaxisRight = 0; //Root-File Objects // TObjArray* hList[pixelSetSize] = {NULL}; TObjArray* hRootList = NULL; TCanvas** cgpPixelPulses = NULL; TCanvas** cgpDistributions = NULL; TCanvas* cFiltered = NULL; TCanvas* cgpTestHistos = NULL; //---------------------------------------------------------------------------- // Functions //---------------------------------------------------------------------------- void BookDebugHistos(int ); void BookTestHistos( int ); //void DeletePixelHistos(int ); //void SaveHistograms( const char*, const char*, TObjArray*, int ); void FillHistograms(Pixel*, vector*, int, int, int , int ); void DrawTestHistograms( int); bool ProduceDebugHistos( vector *pZXings); bool UseThisPulse( int, int, int, float, int, int); void UpdateCanvases( int, int, bool); void DeletePixelCanvases( int, int ); TString SetHostsPaths(TString filename, bool Input); TString gDataDirectory; TString gHomeDirectory; TString ghostname; //---------------------------------------------------------------------------- //---------------------------------------------------------------------------- // MAIN - Funtion //---------------------------------------------------------------------------- //---------------------------------------------------------------------------- int FPulseOverlay( TString datafilename = "raw/2011/11/09/20111109_006.fits.gz", TString drsfilename = "raw/2011/11/09/20111109_003.drs.fits.gz", TString OutRootFileName = "analysis/FPulseTemplate/20111109_006/Overlay/20111109_006.root", bool ProduceGraphic = false, bool spikeDebug = false, bool debugPixel = false, bool testmode = false, bool saveResults = true, int verbosityLevel = 0, // different verbosity levels can be implemented here int firstevent = 1, // Hast to be between 1 and infinity NOT 0!!!! int nevents = 1, int firstpixel = 0, int npixel = 1, int pixelSetSize = 200, int maxPulseOrder = 4, int AmplWindowWidth = 14, //Width of Window for selection of pulses to histograms TString histoOptions = "SRM", float GainMean = 9, float BSLMean = -1, //4 Histogramms will be drawn, decide how far pulsehights differ from eachother int avg1 = 8, int avg2 = 8, int OverlayWindowLeft = 70, int OverlayWindowRight = 230, int refresh_rate = 500 //refresh rate for canvases ) { //---------------------------------------------------------------------------- // Save-Root-File Settings //---------------------------------------------------------------------------- datafilename = SetHostsPaths(datafilename , true ); drsfilename = SetHostsPaths(drsfilename , true ); OutRootFileName = SetHostsPaths(OutRootFileName , false ); if (saveResults) { CreateRootFile( OutRootFileName, false, verbosityLevel ); } //---------------------------------------------------------------------------- // root global Settings //---------------------------------------------------------------------------- gGainMean = GainMean; gBSLMean = BSLMean; gPixelOverlayXaxisLeft = OverlayWindowLeft; gPixelOverlayXaxisRight = OverlayWindowRight; gStyle->SetPalette(1,0); gROOT->SetStyle("Plain"); // gPad->SetGrid(); //---------------------------------------------------------------------------- // root Canvas Settings //---------------------------------------------------------------------------- //Canvas Pad numbering int pulsesCanvasFrameNrs[4] = { 1, // Top left 2, // Top right 3, // bottom left 4 // bootom right }; //Canvas Pad numbering int distributionCanvasFrameNrs[4] = { 1, // Top left 2, // Top right 3, // bottom left 4 // bootom right }; if (ProduceGraphic) { //Canvases cgpPixelPulses = new TCanvas*[maxPulseOrder]; cgpDistributions = new TCanvas*[maxPulseOrder]; //TCanvas* gpcDevelopment = NULL; TString cName = ""; TString cTitle = ""; //naming of pulse canvases for ( int pulse_order = maxPulseOrder - 1; pulse_order >= 0 ; pulse_order-- ) { cName ="cgpDistributions"; cName += pulse_order; cTitle ="Distributions of Pulses with Order of: "; cTitle += pulse_order; cgpDistributions[pulse_order] = new TCanvas(cName,cTitle, 720,pulse_order*20,720,720); cgpDistributions[pulse_order]->Divide(2, 2); cName ="cgpPixelPulses"; cName += pulse_order; cTitle ="Overlays of Pulses with Order of: "; cTitle += pulse_order; cgpPixelPulses[pulse_order] = new TCanvas(cName,cTitle, 0,pulse_order*20,720,720); cgpPixelPulses[pulse_order]->Divide(2, 2); } // 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 ... if (spikeDebug) { cFiltered = new TCanvas("cFiltered","filtered DRS Waveforms", 1080,420,360,360); cFiltered->Divide(1, 3); } if (testmode) { //additional Test histograms cgpTestHistos = new TCanvas( "cgpTestHistos", "Test Histograms", 360, 420, 360, 360 ); cgpTestHistos->Divide(2,0); } } //----------------------------------------------------------------------------- // 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 = NULL; // 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. gNEvents = openDataFits( datafilename, &datafile, AllPixelDataVector, StartCellVector, CurrentEventID, gRegionOfInterest, NumberOfPixels, PXLxROI, verbosityLevel ); if (gNEvents == 0) { cout << "return code of OpenDataFile:" << datafilename<< endl; cout << "is zero -> aborting." << endl; return 1; } if (verbosityLevel > 0) { cout << endl <<"number of events in file: "<< gNEvents << "\t"; } if ( nevents == -1 || nevents > gNEvents ) nevents = gNEvents; // -1 means all! if (verbosityLevel > 0) { cout <<"of, which "<< nevents << " will be processed"<< endl; 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 //----------------------------------------------------------------------------- if (testmode) BookTestHistos( verbosityLevel ); //----------------------------------------------------------------------------- // Main Cycle //----------------------------------------------------------------------------- //---------------------------------------------------------------------------- // Initialize Pixel //---------------------------------------------------------------------------- Pixel** pixel = new Pixel*[NPIX]; for (int i = 0 ; i < NPIX; i++) { pixel[i] = NULL; } int lastPixelOfSet = 0; //------------------------------------- // Loop over Pixel Sets //------------------------------------- for ( int firstPixelOfSet = firstpixel; firstPixelOfSet < firstpixel + npixel; firstPixelOfSet += pixelSetSize ) { lastPixelOfSet = firstPixelOfSet + pixelSetSize-1; if (lastPixelOfSet > 1439) { lastPixelOfSet = 1439; } if (verbosityLevel == 0) { cout << "------------------------------------------------" << endl << "...processing Set from Pixel " << firstPixelOfSet << " to Pixel " << lastPixelOfSet << endl; } //-------------------------------------------------------------------- // Loop over every Event of Pixel Set //-------------------------------------------------------------------- for ( int ev = firstevent ; ev < firstevent + nevents; ev++) { // Get an Event --> consists of 1440 Pixel ...erm....data datafile->GetRow( ev - 1 ); if (verbosityLevel == 1) { cout << "-------------------------------------" << endl << "...processing Set from Pixel " << firstPixelOfSet << " to Pixel " << lastPixelOfSet << "... Event: " << CurrentEventID << "/" << nevents << endl; } //-------------------------------------------------------------------- // Loops over every Pixel of a Set of Pixels //-------------------------------------------------------------------- for ( int pixelID = firstPixelOfSet; pixelID < lastPixelOfSet + 1 && pixelID < firstpixel + npixel; pixelID++ ) { if (verbosityLevel > 1) { cout << "-------------------------------------" << endl << "...processing Set from Pixel " << firstPixelOfSet << " to Pixel " << lastPixelOfSet << "... Event: " << CurrentEventID << "/" << nevents << endl << " Pixel: " << pixelID << "/" << firstpixel + npixel -1 << endl; } //------------------------------------- // Create individual Pixel //------------------------------------- if (ev == firstevent) { pixel[pixelID] = new Pixel( pixelID, maxPulseOrder, verbosityLevel, false, NULL, gPixelOverlayXaxisLeft, gPixelOverlayXaxisRight, gBSLMean, gGainMean, histoOptions ); } //------------------------------------- // Apply Calibration //------------------------------------- if (verbosityLevel > 2) cout << "...applying DrsCalibration"; applyDrsCalibration( Ameas, pixelID, 12, 12, drs_basemean, drs_gainmean, drs_triggeroffsetmean, gRegionOfInterest, AllPixelDataVector, StartCellVector ); if (verbosityLevel > 2) cout << "...done " << endl; //------------------------------------- // Apply Filters //------------------------------------- // finds spikes in the raw data, and interpolates the value // spikes are: 1 or 2 slice wide, positive non physical artifacts if (verbosityLevel > 2) cout << "...removeing Spikes"; removeSpikes (Ameas, Vcorr); if (verbosityLevel > 2) cout << "...done " << endl; // filter Vcorr with sliding average using FIR filter function if (verbosityLevel > 2) cout << "...applying sliding average filter"; sliding_avg(Vcorr, Vslide, avg1); if (verbosityLevel > 2) cout << "...done " << endl; // filter Vslide with CFD using FIR filter function if (verbosityLevel > 2) cout << "...apllying factfir filter"; factfir(b_cfd , a_cfd, k_cfd, Vslide, Vcfd); if (verbosityLevel > 2) cout << "...done " << endl; // filter Vcfd with sliding average using FIR filter function if (verbosityLevel > 2) cout << "...applying 2nd sliding average filter"; sliding_avg(Vcfd, Vcfd2, avg2); if (verbosityLevel > 2) cout << "...done " << endl; //------------------------------------- // Search vor Zero crossings //------------------------------------- if (verbosityLevel > 2) cout << endl << "...searching zero crossings" ; // 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* pZXings = zerosearch( Vcfd2 , 1 , 8); // pZXings means "zero cross ings" EnlargeRegion(*pZXings, 10, 10); findAbsMaxInRegions(*pZXings, Vslide); removeMaximaBelow( *pZXings, 3.0); removeRegionWithMaxOnEdge( *pZXings, 2); removeRegionOnFallingEdge( *pZXings, 100); findTimeOfHalfMaxLeft(*pZXings, Vslide, gBSLMean, 5, 10, verbosityLevel ); removeRegionWithToFlatSlope(*pZXings, 0.5); if (verbosityLevel > 2) cout << "...done" << endl; //----------------------------------------------------------------------------- // Fill Overlay Histos //----------------------------------------------------------------------------- if (verbosityLevel > 2) { cout << endl << "Filling Histograms of Pixel# [Set] " << pixelID << " Pixel# [Got] " << pixel[pixelID]->mChid; } FillHistograms( pixel[pixelID], pZXings, AmplWindowWidth, ev, // histoOptions, maxPulseOrder, verbosityLevel ); //----------------------------------------------------------------------------- // Spike Debug //----------------------------------------------------------------------------- if ( spikeDebug ) { BookDebugHistos(verbosityLevel ); breakout = ProduceDebugHistos( pZXings ); }// end of if(spikeDebug) delete pZXings; if (breakout) break; //------------------------------------- // Draw 1. Set of Pixel Histograms //------------------------------------- if ((ev % refresh_rate) == 0) { if (ProduceGraphic) { if (debugPixel) { // for ( int order = 0; // order < maxPulseOrder; // order++ // ) // { pixel[firstPixelOfSet]->DrawOverlayHistograms( cgpPixelPulses, pulsesCanvasFrameNrs ); pixel[firstPixelOfSet]->DrawDistributionHistograms( cgpDistributions, distributionCanvasFrameNrs ); // } if (testmode) { DrawTestHistograms( verbosityLevel ); } UpdateCanvases( verbosityLevel, maxPulseOrder, testmode ); } } } if (breakout) break; if (verbosityLevel > 1) { cout << endl << "...End of Pixel" << endl << "-------------------------------------"<< endl; } }//End of Loop over Set //------------------------------------- // end of Loops over Set //------------------------------------- }//End of Loop over Events //------------------------------------- // end of Loops over Events //------------------------------------- //------------------------------------- // Draw Pixel Histogramms of Overlay Spectra //------------------------------------- for ( int pixelID = firstPixelOfSet; pixelID < firstPixelOfSet + pixelSetSize && pixelID < firstpixel + npixel; pixelID++ ) { //here is what happends at the end of each loop over all Events if (ProduceGraphic) { pixel[pixelID]->DrawOverlayHistograms( cgpPixelPulses, pulsesCanvasFrameNrs ); pixel[pixelID]->DrawDistributionHistograms( cgpDistributions, distributionCanvasFrameNrs ); UpdateCanvases( verbosityLevel, maxPulseOrder, testmode ); } //Save Histograms of Pixels into Output rootfile if (verbosityLevel > 2) { cout << endl << "Saving Pixel# [Set] " << pixelID << " Pixel# [Got] " << pixel[pixelID]->mChid; } pixel[pixelID]->SavePixelHistograms( OutRootFileName, saveResults ); if (debugPixel) { //Process gui events asynchronously during input cout << endl; TTimer timer("gSystem->ProcessEvents();", 50, kFALSE); timer.TurnOn(); TString input = Getline("Type 'q' to exit, to go on: "); timer.TurnOff(); if (input=="q\n") { break; } } if (verbosityLevel > 2) { cout << endl << "Deleting Pixel# [Set] " << pixelID << " Pixel# [Got] " << pixel[pixelID]->mChid; } delete pixel[pixelID]; if (verbosityLevel > 2) { cout << endl << "...End of Pixel-Set" << endl << "------------------------------------------------" << endl; } } } // End of Loop over all Pixels //------------------------------------- // Draw Histograms //------------------------------------- SaveHistograms( //save histograms of generell results into output root file OutRootFileName, "root", hRootList, saveResults, verbosityLevel ); if (ProduceGraphic) { UpdateCanvases( verbosityLevel, maxPulseOrder, testmode ); //Process gui events asynchronously during input cout << endl; TTimer timer("gSystem->ProcessEvents();", 50, kFALSE); timer.TurnOn(); Getline("Press to go on"); timer.TurnOff(); DeletePixelCanvases( maxPulseOrder, verbosityLevel ); if (testmode) delete cgpTestHistos; if (spikeDebug) delete cFiltered; } delete datafile; delete[] pixel; return( 0 ); } //---------------------------------------------------------------------------- // end of main function //----------------------------------------------------------------------------- //----------------------------------------------------------------------------- // Funktions //----------------------------------------------------------------------------- /* ich erzeuge sowas wie 36 pixel mit ca 20 Histogrammen pro pixel wie viel speicher braucht das? Ich brauche eine möglichkeit jedes dieser histogramme zu identifizieren am besten ein array oder eine liste oder einen abstracten datentyp mit den histogrammen als attribute. ne klasse wäre super */ void BookTestHistos( int verbosityLevel ) { if (verbosityLevel > 2) cout << endl << "...book pixel histograms" << endl; hTesthisto = new TH1F ( "hTesthisto", "Deviation of rising edge and maximum", 600, -10.1, 10.1 ); hTesthisto->GetXaxis()->SetTitle( "Timeslices [a.u.]" ); hTesthisto->GetYaxis()->SetTitle( "counts" ); hTesthisto2 = new TH2F ( "hTesthisto2", "Deviation of rising edge and maximum by event #", gNEvents, 250, 800, 600, -10.1, 10.1 ); // hTesthisto2->GetXaxis()->SetTitle( "Timeslices [a.u.]" ); // hTesthisto2->GetYaxis()->SetTitle( "counts" ); // hTesthisto2->SetMarkerStyle( 2 ); if (verbosityLevel > 2) cout << "...done" << endl; } //end of BookTestHistos //---------------------------------------------------------------------------- void BookDebugHistos( int verbosityLevel ) { if (verbosityLevel > 2) cout << endl << "...book histograms" << endl; hRootList = new TObjArray; debugHistos = new TH1F[ Number_Of_Debug_Histo_Types ]; for ( int type = 0; type < Number_Of_Debug_Histo_Types; type++){ debugHistos[ type ].SetBins(1024, 0, 1024); debugHistos[ type ].SetLineColor(1); debugHistos[ type ].SetLineWidth(2); debugHistos[ type ].SetStats(false); // set X axis paras debugHistos[ type ].GetXaxis()->SetLabelSize(0.08); debugHistos[ type ].GetXaxis()->SetTitleSize(0.08); debugHistos[ type ].GetXaxis()->SetTitleOffset(1.0); debugHistos[ type ].GetXaxis()->SetTitle(Form("Time slice [%.1f ns/slice]", 1./2.)); // set Y axis paras debugHistos[ type ].GetYaxis()->SetLabelSize(0.08); debugHistos[ type ].GetYaxis()->SetTitleSize(0.08); debugHistos[ type ].GetYaxis()->SetTitleOffset(0.3); debugHistos[ type ].GetYaxis()->SetTitle("Amplitude [mV]"); // hRootList->Add( &debugHistos[ type ] ); } if (verbosityLevel > 2) cout << "...done" << endl; } // end of BookDebugHistos //---------------------------------------------------------------------------- void FillHistograms( Pixel* CurrentPixel, vector* pZXings, int AmplWindowWidth, int eventNumber, int maxPulseOrder, int verbosityLevel ) { if (verbosityLevel > 2) cout << endl << "...filling pulse histograms" ; if (verbosityLevel > 3) cout << endl << "...EventNR " << eventNumber ; bool use_this_peak=false; int order_of_pulse=0; vector::iterator reg; //Loop over all found zerocrossings in Region for (reg = pZXings->begin() ; reg < pZXings->end() ; reg++) { //skip those who are at the Rim of the ROI if (reg->maxPos < 12 || (unsigned int) reg->maxPos > gRegionOfInterest-12) { if (verbosityLevel > 2) cout << endl << "\t...out of range" << endl; continue; } // Define axis range of Histogram int Left = reg->maxPos - gPixelOverlayXaxisLeft; int Right = reg->maxPos + gPixelOverlayXaxisRight; int EdgeLeft = reg->halfRisingEdgePos - gPixelOverlayXaxisLeft; int EdgeRight = reg->halfRisingEdgePos + gPixelOverlayXaxisRight; if (Left < 0) Left = 0; if (EdgeLeft < 0) EdgeLeft = 0; if (Right > (int)Ameas.size() ) Right = Ameas.size(); if (EdgeRight > (int)Ameas.size() ) EdgeRight = Ameas.size(); if (verbosityLevel > 2) cout << endl << "\t...choosing Histogram" ; //determine order of pulse and which histogram shall be filled for ( int order = 0; order < maxPulseOrder; order++) { use_this_peak = UseThisPulse( reg->maxPos, order, AmplWindowWidth, gGainMean, maxPulseOrder, verbosityLevel ); if ( use_this_peak) { order_of_pulse = order; break; } } //Fill histograms if (use_this_peak) { //Histograms for Distributions if (CurrentPixel->mOptions.Contains("S") ) { if (verbosityLevel > 2) cout << endl << "\t...filling Edge Histogram" ; CurrentPixel->hSlopeRisingEdge[order_of_pulse]->Fill( reg->slopeOfRisingEdge ) ; } if (CurrentPixel->mOptions.Contains("R") ) { if (verbosityLevel > 2) cout << endl << "\t...filling RisingEdgeToMax Histogram" ; CurrentPixel->hRisingEdgeToMax[order_of_pulse]->Fill( reg->distanceEdgeToMax ) ; } if (CurrentPixel->mOptions.Contains("M") ) { if (verbosityLevel > 2) cout << endl << "\t...filling PosOfMax Histogram" ; CurrentPixel->hPosOfMax[order_of_pulse]->Fill( reg->maxPos ) ; } //Histograms for Maximum Overlay for ( int pos = Left; pos < Right; pos++) { CurrentPixel->hMaxOverlay[order_of_pulse]->Fill( pos - (reg->maxPos), Ameas[pos]) ; CurrentPixel->hMaxProfile[order_of_pulse]->Fill( pos - (reg->maxPos), Ameas[pos]) ; } //Histograms for Edge Overlay for ( int pos = EdgeLeft; pos < EdgeRight; pos++) { CurrentPixel->hEdgeOverlay[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), Ameas[pos]) ; CurrentPixel->hEdgeProfile[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), Ameas[pos]) ; } if (verbosityLevel > 2) cout << "...done" << endl; } //Breakout option if (breakout) { break; } } if (verbosityLevel > 2) cout << "...done" << endl; } // end of FillHistograms //---------------------------------------------------------------------------- bool UseThisPulse( int maxPos, int order, int AmplWindowWidth, float gGainMean, int maxPulseOrder, int verbosityLevel ) { //determine if pulse is of given order bool use_this_peak = false; if (Ameas[maxPos] >= ((gGainMean*(order+1)) - (AmplWindowWidth/2)) && Ameas[maxPos] < ((gGainMean*(order+1)) + AmplWindowWidth/2)) { if (verbosityLevel > 3) cout << "...#" << order ; use_this_peak = true; } else { if (verbosityLevel > 3) { if (order >= (maxPulseOrder - 1)) { cout << "...NONE" << endl ; } } use_this_peak = false; } return use_this_peak; } void DrawTestHistograms( int verbosityLevel ) { if (verbosityLevel > 2) cout << endl << "...drawing pulse histograms TODO CANVAS" ; // cgpTestHistos->cd(1); // hTesthisto->Draw(); // cgpTestHistos->cd(2); // hTesthisto2->Draw("BOX"); } // end of DrawTestHistograms //---------------------------------------------------------------------------- void UpdateCanvases( int verbosityLevel, int max_pulse_order, bool testmode ) { if (verbosityLevel > 3) cout << endl << "...updating canvases" ; for (int pulse_order = 0; pulse_order < max_pulse_order; pulse_order++) { cgpPixelPulses[pulse_order]->Modified(); cgpPixelPulses[pulse_order]->Update(); cgpDistributions[pulse_order]->Modified(); cgpDistributions[pulse_order]->Update(); if ( testmode ) { cgpTestHistos->Modified(); cgpTestHistos->Update(); } } } // end of UpdateCanvases //---------------------------------------------------------------------------- void DeletePixelCanvases( int maxPulseOrder, int verbosityLevel ) { if (verbosityLevel > 2) cout << endl << "...delete pixel Canvases" ; for (int pulse_order = 0; pulse_order < maxPulseOrder; pulse_order++ ) { delete cgpPixelPulses[pulse_order]; delete cgpDistributions[pulse_order]; } delete[] cgpPixelPulses; delete[] cgpDistributions; } //---------------------------------------------------------------------------- //---------------------------------------------------------------------------- bool ProduceDebugHistos( vector* pZXings ) { // 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 < gRegionOfInterest; 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] ); } bool has_negative_slope = false; cFiltered->cd(1); gPad->SetGrid(); debugHistos[Ameas_].Draw(); TF1 *OneEdge; TLine *OneLine; TBox *OneBox; vector MyEdges; vector MyBoxes; vector MyLines; int left = debugHistos[Ameas_].GetXaxis()->GetFirst(); int right = debugHistos[Ameas_].GetXaxis()->GetLast(); vector::iterator reg; for (reg = pZXings->begin() ; reg < pZXings->end() ; reg++) { if (reg->slopeOfRisingEdge < 0) { has_negative_slope = true; } } if (has_negative_slope) { for (unsigned int i=0; isize(); i++){ Region* reg = &pZXings->at(i); cout << "Slope is: " << reg->slopeOfRisingEdge << endl; // xBegin = (reg->halfRisingEdgePos - reg->distanceEdgeToMax); // yBegin = reg->interceptRisingEdge + reg->slopeOfRisingEdge * (xBegin) OneEdge = new TF1( "OneEdge", "[0]*x+[1]", left, right ); OneEdge->SetParameters( reg->slopeOfRisingEdge, reg->interceptRisingEdge ); OneEdge->SetLineColor(kRed); OneEdge->SetLineWidth(1); MyEdges.push_back(OneEdge); OneEdge->Draw("SAME"); // delete OneEdge; OneBox = new TBox( pZXings->at(i).maxPos -10 , pZXings->at(i).maxVal -0.5, pZXings->at(i).maxPos +10 , pZXings->at(i).maxVal +0.5); OneBox->SetLineColor(kBlue); OneBox->SetLineWidth(1); OneBox->SetFillStyle(0); OneBox->SetFillColor(kRed); MyBoxes.push_back(OneBox); OneBox->Draw(); OneLine = new TLine( pZXings->at(i).halfRisingEdgePos, Ameas[pZXings->at(i).halfRisingEdgePos], pZXings->at(i).halfRisingEdgePos, pZXings->at(i).halfRisingEdgePos * pZXings->at(i).slopeOfRisingEdge + pZXings->at(i).interceptRisingEdge ); OneLine->SetLineColor(kBlue); OneLine->SetLineWidth(3); MyLines.push_back(OneBox); OneLine->Draw(); } cFiltered->cd(2); gPad->SetGrid(); debugHistos[Vslide_].Draw(); debugHistos[Vcorr_].Draw("SAME"); for (unsigned int i=0; iDraw("SAME"); } for (unsigned int i=0; iDraw(); } for (unsigned int i=0; iDraw(); } cFiltered->cd(3); gPad->SetGrid(); debugHistos[Vcorr_].SetLineColor(kBlue); debugHistos[Vcorr_].Draw(); debugHistos[Ameas_].Draw("SAME"); debugHistos[Vslide_].SetLineColor(kRed); debugHistos[Vslide_].Draw("SAME"); debugHistos[Vcfd_].SetLineColor(kGreen); debugHistos[Vcfd_].Draw("SAME"); debugHistos[Vcfd2_].SetLineColor(kViolet); debugHistos[Vcfd2_].Draw("SAME"); for (unsigned int i=0; iDraw("SAME"); } for (unsigned int i=0; iDraw(); } for (unsigned int i=0; iDraw(); } // cFiltered->cd(3); // gPad->SetGrid(); // debugHistos[Vcfd2_].Draw(); // TLine *zeroline = new TLine(0, 0, 1024, 0); // zeroline->SetLineColor(kBlue); // zeroline->Draw(); cFiltered->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") { return true; } MyLines.clear(); MyBoxes.clear(); MyEdges.clear(); // for (unsigned int i=0; isize(); i++) // { // delete MyEdges.at(i) // delete MyBoxes.at(i); // } // delete OneBox; // delete OneEdge; // delete OneLine; // MyEdges.erase(); // MyBoxes.erase(); } //TODO!!!!!!!!! // do some Garbage collection ... // all the Objects on the heap should be deleted here. return false; }// end of if(spikeDebug) TString SetHostsPaths(TString filename, bool Input) { ghostname = gSystem->HostName(); TString temp_filename; if ( ghostname.Contains("isdc") ) //IF ONE IS WORKING AT ISDC { gDataDirectory = "/fact/"; gHomeDirectory = "/home_nfs/isdc/jbbuss/"; } if ( ghostname.Contains("hpc") ) //IF ONE IS WORKING AT PHIDO { gDataDirectory = "/fhgfs/groups/app/fact-construction/"; gHomeDirectory = "/home/jbuss/"; } if (Input) { temp_filename = gDataDirectory; temp_filename += filename; filename = temp_filename; } else if (!Input) { temp_filename = gHomeDirectory; temp_filename += filename; filename = temp_filename; } return filename; } int main() { FPulseOverlay(); return 0; }