/********************** FPulseOverlay ***********************
* read FACT raw data
* remove spikes
* calculate baseline
* find peaks (CFD and normal discriminator)
+ seagRCh 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 <TROOT.h>
#include <TCanvas.h>
#include <TProfile.h>
#include <TTimer.h>
#include <TH1F.h>
#include <TH2F.h>
#include <Getline.h>
#include <TLine.h>
#include <TBox.h>
#include <TMath.h>
#include <TFile.h>
#include <TStyle.h>
#include <TString.h>
#include <TObjArray.h>
#include <TF1.h>
#include<TSystem.h>

#include <stdio.h>
#include <stdint.h>
#include <cstdio>
#include <iostream>
#include <fstream>
using namespace std;

#define NPIX  1440
#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.h"
# 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"

#include "configfile.h"


//----------------------------------------------------------------------------
// Decleration of global variables
//----------------------------------------------------------------------------

bool            gBreakout=false;

vector<int16_t> gAllPixelDataVector;
vector<int16_t> gStartCellVector;
unsigned int    gCurrentEventID;
size_t          gPXLxROI;
UInt_t          gRegionOfInterest;
UInt_t          gNumberOfPixels;
TString         ghistotitle;

size_t          gTriggerOffsetROI;
size_t          gRC;
size_t          gDrs_n;
vector<float>   gDrsBaseMean;
vector<float>   gDrsGainMean;
vector<float>   gDrsTriggeroffsetMean;

vector<float>   gAmeas(FAD_MAX_SAMPLES);   // copy of the data (measured amplitude
vector<float>   gVcorr(FAD_MAX_SAMPLES);   // corrected Values
vector<float>   gVslide(FAD_MAX_SAMPLES);  // sliding average result
vector<float>   gVcfd(FAD_MAX_SAMPLES);    // CDF result
vector<float>   gVcfd2(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
        gAmeas_ 	= 0,
        N1mean_ = 1,
        gVcorr_  = 2,
        Vtest_  = 3,
        gVslide_ = 4,
        gVcfd_   = 5,
        gVcfd2_  = 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};
    TList*  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<Region>*, int, int, int , int );
void DrawTestHistograms( int);
bool ProduceDebugHistos( vector<Region> *pZXings);
bool  UseThisPulse( int, int, int, float, int, int);
void UpdateCanvases( int, int, bool);
void DeletePixelCanvases( int, int );

int HandleFitsFile(
        fits *      datafile,
        TString     datafilename,
        TString     drsfilename,
        int         nevents,
        int         npixel,
        int         verbosityLevel
        );

void ReadSequenzFile(
        TString fileName,
        vector<TString> *fitsFiles,
        vector<TString> *drsFiles
        );

TString gDataDirectory;
TString gHomeDirectory;
TString ghostname;

//----------------------------------------------------------------------------
//----------------------------------------------------------------------------
// MAIN - Funtion
//----------------------------------------------------------------------------
//----------------------------------------------------------------------------
int FPulseOverlay_filelist(
        TString     inputPath           = "raw/2012/03/09/",
        TString     sequenzFileName = "/home/jbuss/macros/isdc_rootmacros/PulseTemplates/seq",
//        TString     dataFileSuffix      = ".fits.gz",
//        int         firstRunId          = 18,
//        int         nRuns               = 50,
//        TString     drsfilename         = "20120309_012.drs.fits.gz",
        TString     OutputPath          = "analysis/FPulseTemplate/20120309_018/Overlay/",
        TString     OutRootFileName     = "20120309_018.root",
        int         firstevent          = 1,        // Hast to be between 1 and infinity NOT 0!!!!
        int         nevents             = 100,
        int         firstpixel          = 0,
        int         npixel              = 2,
        int         pixelSetSize       = 150,
        int         maxPulseOrder       = 1,
        int         verbosityLevel      = 6,        // different verbosity levels can be implemented here
        TString     histoOptions        = "SRM",
        int         AmplWindowWidth     = 14,       //Width of Window for selection of pulses to histograms
        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,
        bool        ProduceGraphic      = true,
        bool        spikeDebug          = false,
        bool        debugPixel          = false,
        bool        stats               = false,
        bool        saveResults         = false,
        bool        testmode            = false,
        int         refresh_rate        = 500      //refresh rate for canvases
        )
{
//-----------------------------------------------------------------------------
// prepare datafiles
//-----------------------------------------------------------------------------

    inputPath       = SetHostsPaths(false, inputPath );

    // declare 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.

    //vector of fits files to read
    vector<TString> fitsFileNames;
    vector<TString> drsFileNames;

    ReadSequenzFile(
            sequenzFileName,
            &fitsFileNames,
            &drsFileNames
            );

    cout << "filelist\n" ;

    //fill vector of fits files to read
    cout << "generating file list:" << endl;
    for (unsigned int fileNr = 0; fileNr < fitsFileNames.size(); fileNr++)
    {
        if ( fileNr >= drsFileNames.size() )
        {
            cout << "too less drs-files in vector, skipping rest" << endl;
            break;

            drsFileNames.resize(fileNr);
            fitsFileNames.resize(fileNr);
        }

        cout << "...checking raw files";
//        cout << "\t fits:\t" << fitsFileNames.at(fileNr);
//        cout << "\n\t\t drs:\t" << drsFileNames.at(fileNr) << endl;

        //filename/path construction
        TString fits_filename = inputPath;
        fits_filename.Append(fitsFileNames.at(fileNr));
        TString drs_filename = inputPath;
        drs_filename.Append(drsFileNames.at(fileNr));
        cout << "\n\t fits:\t" << fits_filename;
        cout << "\n\t drs:\t" << drs_filename << endl;

        //check if files are readable
        ifstream fitsFile(fits_filename);
        ifstream drsFile(drs_filename);
        if (fitsFile.good() && drsFile.good())
        {
            //fill vector
            cout << "...added\n"  << endl;
            fitsFileNames.at(fileNr) = fits_filename;
            drsFileNames.at(fileNr) = drs_filename;
        }
        else {
            cout << "...FAILED:";
            cout << "\t deleting line"
//                 << fitsFileNames.at(fileNr)
//                 << " and "
//                 << drsFileNames.at(fileNr)
                 << " from list\n" << endl;

            fitsFileNames.erase(fitsFileNames.begin()+fileNr);
            drsFileNames.erase(drsFileNames.begin()+fileNr);
            fileNr--;
//            continue;
        }
    }
    cout << endl << "size of file vectors: fits: "
         << fitsFileNames.size()
         << " drs: "
         << drsFileNames.size() << endl;

//----------------------------------------------------------------------------
//	Save-Root-File Settings
//----------------------------------------------------------------------------

    OutputPath      = SetHostsPaths(true, OutputPath );
    OutRootFileName.Prepend(OutputPath);
//    drsfilename         = SetHostsPaths(false,   drsfilename );
//    OutRootFileName     = SetHostsPaths(true,  OutRootFileName );
    cout <<  endl;
//    cout << "drs filename:\t" << drsfilename << endl;
//    cout << "data filename:\t" << datafilename << endl;
    cout << "out filename:\t" << OutRootFileName << endl;


    if (saveResults)
    {
        CreateRootFile( OutRootFileName, true, 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<double> a_cfd(k_cfd, 0);
    double b_cfd = 1.;
    a_cfd[0]=-0.75;
    a_cfd[k_cfd-1]=1.;

//-----------------------------------------------------------------------------
// Book the histograms
//-----------------------------------------------------------------------------



    if (testmode)
    {
        BookTestHistos( verbosityLevel );
    }
//----------------------------------------------------------------------------
// Initialize Pixel
//----------------------------------------------------------------------------
    Pixel** pixel = new Pixel*[NPIX];

    for (int i = 0 ; i < NPIX; i++)
    {
        pixel[i] = NULL;
    }


//-----------------------------------------------------------------------------
// Main Cycle
//-----------------------------------------------------------------------------

    if ( pixelSetSize == -1 )
    {
        pixelSetSize = firstpixel +npixel;
    }
    int lastPixelOfSet = 0;
//-------------------------------------
// Loop over Pixel Sets
//-------------------------------------
    for ( int firstPixelOfSet = firstpixel;
          firstPixelOfSet < firstpixel + npixel;
          firstPixelOfSet += pixelSetSize )
    {
        lastPixelOfSet = firstPixelOfSet + pixelSetSize-1;
        if (lastPixelOfSet > firstpixel + npixel)
        {
            lastPixelOfSet = firstpixel + npixel;
        }
        if (lastPixelOfSet > 1439)
        {
            lastPixelOfSet = 1439;
        }
        if (verbosityLevel == 0)
        {
            cout << "------------------------------------------------" << endl
                 << "...processing Set from Pixel "
                 << firstPixelOfSet
                 << " to Pixel "
                 << lastPixelOfSet << endl;
        }

        //-------------------------------------
        // Loop over Fits Files
        //-------------------------------------

        for( unsigned int fileId = 0; fileId < fitsFileNames.size(); fileId++ )
        {
            datafile = NULL;
//            int return_val = 0;
//            return_val = HandleFitsFile(
//                        datafile,
//                        datafiles.at(fileId),
//                        drsfilename,
//                        nevents,
//                        npixel,
//                        verbosityLevel
//                        );
//            if (return_val)
//            {
//                cerr << "ERROR" << endl;
//                return return_val;
//            }
//==================================================

                // 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.
                cout << "\nopening " << fitsFileNames.at(fileId) << endl;

                gNEvents = openDataFits(
                                fitsFileNames.at(fileId),
                                &datafile,
                                gAllPixelDataVector,
                                gStartCellVector,
                                gCurrentEventID,
                                gRegionOfInterest,
                                gNumberOfPixels,
                                gPXLxROI,
                                verbosityLevel
                            );

                if (gNEvents == 0)
                {
                    cout << "return code of OpenDataFile:" << fitsFileNames.at(fileId)<< 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: "<< gNumberOfPixels << "\t";
                }

                if ( npixel == -1 || npixel > (int)gNumberOfPixels  ) npixel = (int)gNumberOfPixels; // -1 means all!

                if (verbosityLevel > 0)
                {
                    cout <<"of, which "<< npixel << " will be processed"<< endl;
                }

                //Get the DRS calibration
                gRC = openCalibFits(
                            drsFileNames.at(fileId),
                            gDrsBaseMean,
                            gDrsGainMean,
                            gDrsTriggeroffsetMean,
                            gTriggerOffsetROI
                            );

                if (gRC == 0)
                {
                 cout << "return code of openCalibFits:" << drsFileNames.at(fileId) << endl;
                 cout << "is zero -> aborting." << endl;
                 return 1;
                }

            //==================================================================

//--------------------------------------------------------------------
// 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: " << gCurrentEventID
                      << "/" << nevents << endl;
                }

                cout << "\t\t success!" << 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: " << gCurrentEventID
                          << "/" << nevents << endl
                          << " Pixel: " << pixelID
                          << "/" << firstpixel + npixel -1 << endl;
                    }

    //-------------------------------------
    // Create individual Pixel
    //-------------------------------------
                    if (ev == firstevent)
                    {
                        pixel[pixelID] = new Pixel(
                                    pixelID,
                                    maxPulseOrder,
                                    verbosityLevel,
                                    stats,
                                    histoOptions,
                                    gPixelOverlayXaxisLeft,
                                    gPixelOverlayXaxisRight,
                                    gBSLMean,
                                    gGainMean
                                    );
                    }

    //-------------------------------------
    // Apply Calibration
    //-------------------------------------
                    if (verbosityLevel > 2) cout << "...applying DrsCalibration";
                    applyDrsCalibration(
                                gAmeas,
                                pixelID,
                                12,
                                12,
                                gDrsBaseMean,
                                gDrsGainMean,
                                gDrsTriggeroffsetMean,
                                gRegionOfInterest,
                                gAllPixelDataVector,
                                gStartCellVector
                                );
                    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 (gAmeas, gVcorr);
                    if (verbosityLevel > 2) cout << "...done " << endl;

                    // filter gVcorr with sliding average using FIR filter function
                    if (verbosityLevel > 2) cout << "...applying sliding average filter";
                    sliding_avg(gVcorr, gVslide, avg1);
                    if (verbosityLevel > 2) cout << "...done " << endl;

                    // filter gVslide with CFD using FIR filter function
                    if (verbosityLevel > 2) cout << "...apllying factfir filter";
                    factfir(b_cfd , a_cfd, k_cfd, gVslide, gVcfd);
                    if (verbosityLevel > 2) cout << "...done " << endl;

                    // filter gVcfd with sliding average using FIR filter function
                    if (verbosityLevel > 2) cout << "...applying 2nd sliding average filter";
                    sliding_avg(gVcfd, gVcfd2, avg2);
                    if (verbosityLevel > 2) cout << "...done " << endl;

    //-------------------------------------
    // SeagRCh vor Zero crossings
    //-------------------------------------
                    if (verbosityLevel > 2) cout << endl << "...seagRChing zero crossings" ;
                    // peaks in gAmeas[] are found by seagRChing for zero crossings
                    // in gVcfd2
                    // first Argument 1 means ... *rising* edge
                    // second Argument 1 means ... seagRCh with stepsize 1 ... 10 is okay as well
                    vector<Region>* pZXings = zerosearch( gVcfd2 , 1 , 8);
                    // pZXings means "zero cross ings"
                    EnlargeRegion(*pZXings, 10, 10);
                    findAbsMaxInRegions(*pZXings, gVslide);
                    removeMaximaBelow( *pZXings, 3.0);
                    removeRegionWithMaxOnEdge( *pZXings, 2);
                    removeRegionOnFallingEdge( *pZXings, 100);
                    findTimeOfHalfMaxLeft(*pZXings, gVslide, 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 );
                    gBreakout = ProduceDebugHistos( pZXings );
                    }// end of if(spikeDebug)

                    delete pZXings;
                    if (gBreakout)   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 (gBreakout) 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
    //-------------------------------------

        }//End of Loop over Files
    //-------------------------------------
    // end of Loops over Files
    //-------------------------------------


//-------------------------------------
// 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, <return> 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
//-------------------------------------
    TString test = "root";
    SaveHistograms(     //save histograms of generell results into output root file
                OutRootFileName,
                test,
                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 <return> 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 TList();
    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<Region>* 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<Region>::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)gAmeas.size() )         Right       =   gAmeas.size();
        if (EdgeRight > (int)gAmeas.size() )     EdgeRight   =   gAmeas.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), gAmeas[pos]) ;
                CurrentPixel->hMaxProfile[order_of_pulse]->Fill( pos - (reg->maxPos), gAmeas[pos]) ;
            }
            //Histograms for Edge Overlay
            for ( int pos = EdgeLeft; pos < EdgeRight; pos++)
            {
                CurrentPixel->hEdgeOverlay[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), gAmeas[pos]) ;
                CurrentPixel->hEdgeProfile[order_of_pulse]->Fill( pos - (reg->halfRisingEdgePos), gAmeas[pos]) ;
            }

            if (verbosityLevel > 2) cout << "...done" << endl;
        }

        //gBreakout option
        if (gBreakout)
        {
            break;
        }
   }
    if (verbosityLevel > 2) cout << "...done" << endl;
}
// end of FillHistograms
//----------------------------------------------------------------------------

bool UseThisPulse(
        int             maxPos,
        int             order,
        int             AmplWindowWidth,
        float           gainMean,
        int             maxPulseOrder,
        int             verbosityLevel
        )
{
    //determine if pulse is of given order
bool use_this_peak = false;
        if (gAmeas[maxPos] >= ((gainMean*(order+1)) - (AmplWindowWidth/2))
            && gAmeas[maxPos] < ((gainMean*(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<Region>* pZXings
        )
{
    // TODO do this correct. The vectors should be the rigt ones... this is just luck
    debugHistos[gAmeas_].GetXaxis()->Set(gAmeas.size() , -0.5 , gAmeas.size()-0.5);
    debugHistos[gVcorr_].GetXaxis()->Set(gAmeas.size() , -0.5 , gAmeas.size()-0.5);
    debugHistos[gVslide_].GetXaxis()->Set(gAmeas.size() , -0.5 , gAmeas.size()-0.5);
    debugHistos[gVcfd_].GetXaxis()->Set(gAmeas.size() , -0.5 , gAmeas.size()-0.5);
    debugHistos[gVcfd2_].GetXaxis()->Set(gAmeas.size() , -0.5 , gAmeas.size()-0.5);

    for ( unsigned int sl = 0; sl < gRegionOfInterest; sl++)
    {
       debugHistos[gAmeas_].SetBinContent(sl, gAmeas[sl]);
       debugHistos[gVcorr_].SetBinContent(sl, gVcorr[sl]);
       debugHistos[gVslide_].SetBinContent( sl, gVslide[sl] );
       debugHistos[gVcfd_].SetBinContent( sl, gVcfd[sl] );
       debugHistos[gVcfd2_].SetBinContent( sl, gVcfd2[sl] );
    }

    bool has_negative_slope = false;
    cFiltered->cd(1);
    gPad->SetGrid();
    debugHistos[gAmeas_].Draw();

    TF1             *OneEdge;
    TLine           *OneLine;
    TBox            *OneBox;
    vector<TF1*>    MyEdges;
    vector<TBox*>   MyBoxes;
    vector<TBox*>   MyLines;

    int left    = debugHistos[gAmeas_].GetXaxis()->GetFirst();
    int right   = debugHistos[gAmeas_].GetXaxis()->GetLast();
    vector<Region>::iterator region;
    for (region = pZXings->begin() ; region < pZXings->end() ; region++)
    {
        if (region->slopeOfRisingEdge < 0)
        {
            has_negative_slope = true;
        }
    }

    if (has_negative_slope)
    {
    for (unsigned int i=0; i<pZXings->size(); i++){
        Region* reg = &pZXings->at(i);
        cout << "Slope is: " << reg->slopeOfRisingEdge << endl;

//            xBegin = (reg->halfRisingEdgePos - reg->distanceEdgeToMax);
//            yBegin = reg->integRCeptRisingEdge + 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,
                    gAmeas[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[gVslide_].Draw();
    debugHistos[gVcorr_].Draw("SAME");

    for (unsigned int i=0; i<MyEdges.size(); i++){
        MyEdges.at(i)->Draw("SAME");
    }
    for (unsigned int i=0; i<MyBoxes.size(); i++){
        MyBoxes.at(i)->Draw();
    }
    for (unsigned int i=0; i<MyLines.size(); i++){
        MyLines.at(i)->Draw();
    }

    cFiltered->cd(3);
    gPad->SetGrid();
    debugHistos[gVcorr_].SetLineColor(kBlue);
    debugHistos[gVcorr_].Draw();
    debugHistos[gAmeas_].Draw("SAME");
    debugHistos[gVslide_].SetLineColor(kRed);
    debugHistos[gVslide_].Draw("SAME");
    debugHistos[gVcfd_].SetLineColor(kGreen);
    debugHistos[gVcfd_].Draw("SAME");
    debugHistos[gVcfd2_].SetLineColor(kViolet);
    debugHistos[gVcfd2_].Draw("SAME");
    for (unsigned int i=0; i<MyEdges.size(); i++){
        MyEdges.at(i)->Draw("SAME");
    }
    for (unsigned int i=0; i<MyBoxes.size(); i++){
        MyBoxes.at(i)->Draw();
    }
    for (unsigned int i=0; i<MyLines.size(); i++){
        MyLines.at(i)->Draw();
    }


//            cFiltered->cd(3);
//            gPad->SetGrid();
//            debugHistos[gVcfd2_].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, <return> to go on: ");
        timer.TurnOff();
        if (input=="q\n") {
                return true;
        }

        MyLines.clear();
        MyBoxes.clear();
        MyEdges.clear();

//        for (unsigned int i=0; i<pZXings->size(); 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)

int HandleFitsFile(
        fits *      datafile,
        TString     datafilename,
        TString     drsfilename,
        int         nevents,
        int         npixel,
        int         verbosityLevel
        )
{
    // 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.
    cout << "\nopening " << datafilename << endl;

    gNEvents = openDataFits(
                    datafilename,
                    &datafile,
                    gAllPixelDataVector,
                    gStartCellVector,
                    gCurrentEventID,
                    gRegionOfInterest,
                    gNumberOfPixels,
                    gPXLxROI,
                    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: "<< gNumberOfPixels << "\t";
    }

    if ( npixel == -1 || npixel > (int)gNumberOfPixels  ) npixel = (int)gNumberOfPixels; // -1 means all!

    if (verbosityLevel > 0)
    {
        cout <<"of, which "<< npixel << " will be processed"<< endl;
    }

    //Get the DRS calibration
    gRC = openCalibFits(
                drsfilename,
                gDrsBaseMean,
                gDrsGainMean,
                gDrsTriggeroffsetMean,
                gTriggerOffsetROI
                );

    if (gRC == 0)
    {
     cout << "return code of openCalibFits:" << drsfilename << endl;
     cout << "is zero -> aborting." << endl;
     return 1;
    }
    return 0;
}

void ReadSequenzFile(TString fileName, vector<TString> *fitsFiles, vector<TString> *drsFiles)
{
    ifstream seqFile;
    seqFile.open(fileName,ios::in);

    TString fitsName, drsName;

    if(!seqFile)
    {
        cerr << "***** no config-file in the user defined directory.\n"
             << "Check the path .\n" << endl;
        return;
    }

    if(!seqFile.good()) {
        cerr << "***** not able to read config-file." << endl;
        return;
    }

    while(!seqFile.eof())
    {
        seqFile >> fitsName >> drsName;

        if (!fitsName.IsWhitespace() && !fitsName.IsNull())
        {
            fitsFiles->push_back(fitsName);

            if (!drsName.IsWhitespace() && !drsName.IsNull())
            {
                drsFiles->push_back(drsName);
            }
            else drsFiles->push_back("0");
        }
    }

    seqFile.close();

    return;
}


//int main(int argc,char *argv[])
//{

//    TString     test;
//    TString     gRCFileName;
//    TString     processType             = "overlay";
//    bool        gRCFileNameCmdSet        = false;
//    int         verbLevel               = 0; // different verbosity levels can be implemented here
//    bool        verbLevelCmdSet         = false;
//    bool        save                    = false;
//    bool        produceGraphic          = false;

//    // decode arguments
//    if(argc < 2)
//    {
//        printf("no arguements given, using standard arguments\n");
//    }

//    // set conditions for functions arguments
//    for (int i=1;i<argc;i++)
//    {
//        test = argv[i];

//        if (test.Contains("--config_file") || test.Contains("-c"))
//        {
//            cout << "gRC-File: \"" << argv[i + 1] << "\"" << endl;
//            gRCFileName = argv[i + 1];
//            gRCFileNameCmdSet = true;
//            continue;
//        }

//        if (test.Contains("--verbosity") || test.Contains("-v"))
//        {
//            cout << "Verbosity Level: \"" << argv[i + 1] << "\"" << endl;
//            verbLevel = atoi(argv[i + 1]);
//            continue;
//        }

//        if (test.Contains("--save") || test.Contains("-s"))
//        {
//            cout << "will save results" << endl;
//            save = true;
//            continue;
//        }

//        if (test.Contains("--graphics") || test.Contains("-g"))
//        {
//            cout << "will produce graphics" << endl;
//            produceGraphic = true;
//            continue;
//        }
//    }

//    // reading gRC-File:
//    configfile gRCfile( gRCFileName, processType );

//    if (save)
//    {
//        gRCfile.mSave = save;
//    }

//    if (verbLevelCmdSet)
//    {
//        gRCfile.mVerbLevel = verbLevel;
//    }

//    if (produceGraphic)
//    {
//        gRCfile.mProduceGraphic = produceGraphic;
//    }

//    FPulseOverlay(
//                gRCfile.mInputPath,
//                gRCfile.mDataFileName,
//                gRCfile.mDrsFileName,
//                gRCfile.mOutputPath,
//                gRCfile.mOutputFile,
//                gRCfile.mFirstEvent,
//                gRCfile.mNumEvents,
//                gRCfile.mFirstPixel,
//                gRCfile.mNumPixel,
//                gRCfile.mPixelSetSize,
//                gRCfile.mMaxOrder,
//                gRCfile.mVerbLevel,
//                gRCfile.mHistoOptions,
//                gRCfile.mAmplWindowWidth,
//                gRCfile.mGainMean,
//                gRCfile.mBSLMean,
//                gRCfile.mAvg1,
//                gRCfile.mAvg2,
//                gRCfile.mOverlayWindowLeft,
//                gRCfile.mOverlayWindowRight,
//                gRCfile.mProduceGraphic,
//                gRCfile.mSpikeDebug,
//                gRCfile.mDbgPixel,
//                gRCfile.mPrintStats,
//                gRCfile.mSave,
//                gRCfile.mTestMode,
//                gRCfile.mRefreshRate
//                );
//    return 0;

//}