#include "MonteCarlo.h"

MonteCarlo::MonteCarlo()
{
    InitVariables();

    return;
}

MonteCarlo::MonteCarlo(
        TString filename
        )
{
    InitVariables();

    mFileName   = filename;

    OpenRootFile();

    LoadEventTree(  "Events");
    LoadHeaderTree( "RunHeaders");

    ReadRunHeader();

    return;
}

MonteCarlo::MonteCarlo(
        TString filename,
        TString evtsTreeName
        )
{
    InitVariables();

    mFileName   = filename;

    OpenRootFile();

    LoadEventTree(  evtsTreeName);
    LoadHeaderTree( "RunHeaders");

    ReadRunHeader();

    return;
}

MonteCarlo::MonteCarlo(
        TString filename,
        TString evtsTreeName,
        TString headerTreeName
        )
{
    InitVariables();

    mFileName   = filename;

    OpenRootFile();

    LoadEventTree(  evtsTreeName);
    LoadHeaderTree( headerTreeName);

    ReadRunHeader();

    return;
}

MonteCarlo::~MonteCarlo(void)
{
    delete[] mpPixel;

    CloseRootFile();

    delete mpRootFile;

    return;
}

void
MonteCarlo::InitVariables()
{
    mFileName       = "";
    mSeparator      = " ";

    // Änderung !!!
    mpRootFile      = NULL;
    // Änderung zuende

    mpEventTree     = NULL;
    mpHeaderTree    = NULL;
    mpPixel         = NULL;

    // Änderung !!!
    mpIntendedPulsePos      = NULL;
    mpMcRunHeader           = NULL;
    mpGeomCam               = NULL;
    mpRawRunHeader          = NULL;
    mpCorsikaRunHeader      = NULL;

    mpElectronicNoise       = NULL;
    mpRawEventData          = NULL;
    mpIncidentAngle         = NULL;
    mpMcEventMetaData       = NULL;
    mpRawEventHeader        = NULL;
    mpCorsikaEvtHeader      = NULL;
    // Änderung zuende

    mpSamples       = NULL;
    mEventNumber    = 0;
    mNumberOfEvents = 2;

    mVerbosityLvl   = 0;

    return;
}

void
MonteCarlo::SetVerbosityLevel(int verbLvl)
{
    mVerbosityLvl = verbLvl;

    return;
}

int
MonteCarlo::GetVerbosityLevel()
{
    return mVerbosityLvl;
}

void
MonteCarlo::OpenRootFile()
{
    if (mVerbosityLvl > 0) cout << "...opening root file: " << mFileName << endl;

    mpRootFile = new TFile(mFileName, "READ");

    if (!mpRootFile->IsOpen())
    {
        cout << "ERROR - Could not find file " << mFileName << endl;
        return ;
    }

    return;
}

void
MonteCarlo::CloseRootFile()
{
    if (mVerbosityLvl > 0) cout << "...closing root file: " << mFileName << endl;
    mpRootFile->Close("R");

    return;
}

void
MonteCarlo::OpenCsvFile(TString fileName)
{
    mCsvFileName    = fileName;

    if (mVerbosityLvl > 0) cout << "...opening csv file: " << mCsvFileName << endl;

    mCsv.open( mCsvFileName );

    return;
}

void
MonteCarlo::CloseCsvFile()
{
    if (mVerbosityLvl > 0) cout << "...closing csv file: " << mCsvFileName << endl;

    mCsv.close();

    return;
}

void
MonteCarlo::LoadHeaderTree(TString treeName)
{
    if (mVerbosityLvl > 0) cout << "...loading run header tree" << endl;

    mpHeaderTree       = (TTree*)mpRootFile->Get(treeName);

    if (mpHeaderTree->IsZombie())
    {
        cout << "...could not load tree" << endl;
        return;
    }



    //Set Adresses to Branches in RunHeader-Tree
    mpHeaderTree->SetBranchAddress("MGeomCam.",         &mpGeomCam);
    mpHeaderTree->SetBranchAddress("IntendedPulsePos.", &mpIntendedPulsePos);
    mpHeaderTree->SetBranchAddress("MMcRunHeader.",     &mpMcRunHeader);
    mpHeaderTree->SetBranchAddress("ElectronicNoise.",  &mpElectronicNoise);
    mpHeaderTree->SetBranchAddress("MRawRunHeader.",    &mpRawRunHeader);
    mpHeaderTree->SetBranchAddress("MCorsikaRunHeader.",&mpCorsikaRunHeader);
    return;
}

void
MonteCarlo::ReadRunHeader()
{
    if (mVerbosityLvl > 0)
        cout << "...reading run header " << mpHeaderTree << endl;

    //Read Values from RunHeader-Tree

    mpHeaderTree->GetEntry();    //This causes problems


    //Get values from "MGeomCam" Branch
        //Getter functions from "MGeomCamFACT.h"
    mCamDist                = mpGeomCam->GetCameraDist();
    mNumberOfPixels         = mpGeomCam->GetNumPixels();
    mNumberOfSectors        = mpGeomCam->GetNumSectors();
    mNumberOfAreas          = mpGeomCam->GetNumAreas();

    //Get values from "IntendedPulsePos" Branch
    mIntendedPulsePos       = mpIntendedPulsePos->GetVal();

    //Get values from "MMcRunHeader" Branch from event Tree
    mNumSimulatedShowers    = mpMcRunHeader->GetNumSimulatedShowers();

    //Get number of Events from event Tree
    mNumberOfEntries        = mpEventTree->GetEntries();
    if (mVerbosityLvl > 0)
        cout << "Event Tree has " << mNumberOfEntries << " entries" << endl;

    //Get values from "MRawRunHeader" Branch
    mNumberOfEvents         = mpRawRunHeader->GetNumEvents(); // couses problems
    mNumberOfEventsRead     = mpRawRunHeader->GetNumEventsRead();
    mSamplingFrequency      = mpRawRunHeader->GetFreqSampling();
    mSourceName             = mpRawRunHeader->GetSourceName();
    mFileNumber             = mpRawRunHeader->GetFileNumber();
    mNumberOfSamples        = mpRawRunHeader->GetNumSamplesHiGain();
    mNumberOfBytes          = mpRawRunHeader->GetNumBytesPerSample();
    mRunNumber              = mpRawRunHeader->GetRunNumber();
    mRunType                = mpRawRunHeader->GetRunType();


    //Get values from "MCorsikaRunHeader" Branch
    mSlopeSpectrum          = mpCorsikaRunHeader->GetSlopeSpectrum();
    mEnergyMin              = mpCorsikaRunHeader->GetEnergyMin();
    mEnergyMax              = mpCorsikaRunHeader->GetEnergyMax();
    mZdMin                  = mpCorsikaRunHeader->GetZdMin();
    mZdMax                  = mpCorsikaRunHeader->GetZdMax();
    mAzMin                  = mpCorsikaRunHeader->GetAzMin();
    mAzMax                  = mpCorsikaRunHeader->GetAzMax();

    return;
}

void
MonteCarlo::LoadEventTree(TString treeName)
{
    if (mVerbosityLvl > 0) cout << "...loading event tree" << endl;

    mpEventTree       = (TTree*)mpRootFile->Get(treeName);

    if (mpEventTree->IsZombie())
    {
        cout << "...could not load tree" << endl;
        return;
    }

    //Set Adresses to Branches in Events-Tree
    mpEventTree->SetBranchAddress("MRawEvtData.",       &mpRawEventData);
    mpEventTree->SetBranchAddress("IncidentAngle.",     &mpIncidentAngle);
    mpEventTree->SetBranchAddress("MMcEvt.",            &mpMcEventMetaData);
    mpEventTree->SetBranchAddress("MRawEvtHeader.",     &mpRawEventHeader);
    mpEventTree->SetBranchAddress("MCorsikaEvtHeader.", &mpCorsikaEvtHeader);

    return;
}

void
MonteCarlo::ReadEventMetaData()
{
    if (mVerbosityLvl > 1) cout << "...reading event header" << endl;

    //Get values from "MGeomCamFACT" Branch
    mIncidentAngle          = mpIncidentAngle->GetVal();

    //Get values from "MMcEvt" Branch
        //The following Getter Functions can be found in MMcEvt.h
    mCorsikaEventNumber     = mpMcEventMetaData->GetEvtNumber();
    mPhotElFromShower       = mpMcEventMetaData->GetPhotElfromShower();
    mPhotElinCamera         = mpMcEventMetaData->GetPhotElinCamera();
        //The following Getter Functions can be found in MMcEvtBasic.h
    mPartId                 = mpMcEventMetaData->GetPartId();
    mPartName               = mpMcEventMetaData->GetParticleName(mPartId);
    mPartSymbol             = mpMcEventMetaData->GetParticleSymbol(mPartId);
    mEnergy                 = mpMcEventMetaData->GetEnergy();
    mImpact                 = mpMcEventMetaData->GetImpact();
    mTelescopePhi           = mpMcEventMetaData->GetTelescopePhi();
    mTelescopeTheta         = mpMcEventMetaData->GetTelescopeTheta();
    mPhi                    = mpMcEventMetaData->GetParticlePhi();
    mTheta                  = mpMcEventMetaData->GetParticleTheta();

    //Get values from "MRawEventHeader" Branch
    mEventNumber            = mpRawEventHeader->GetDAQEvtNumber();
    mNumTriggerLvl1         = mpRawEventHeader->GetNumTrigLvl1();
    mNumTriggerLvl2         = mpRawEventHeader->GetNumTrigLvl2();

    //Get values from "MCorsikaEvtHeader" Branch
    mFirstInteractionHeight = mpCorsikaEvtHeader->GetFirstInteractionHeight();
    mEvtReuse               = mpCorsikaEvtHeader->GetNumReuse();
    mMomentumX              = mpCorsikaEvtHeader->GetMomentum().X();
    mMomentumY              = mpCorsikaEvtHeader->GetMomentum().Y();
    mMomentumZ              = mpCorsikaEvtHeader->GetMomentum().Z();
    mZd                     = mpCorsikaEvtHeader->GetZd();
    mAz                     = mpCorsikaEvtHeader->GetAz();
    mX                      = mpCorsikaEvtHeader->GetX();
    mY                      = mpCorsikaEvtHeader->GetY();

//    mWeightedNumPhotons     = mpCorsikaEvtHeader->Get;
    // no getter function, no idea how to extract information

    return;
}

void
MonteCarlo::ReadEventRawData()
{
    if (mVerbosityLvl > 1) cout << "...reading event raw data" << endl;

    // delete Pixel Array before you set refferences for a new one
    if (mpPixel != NULL)
    {
        delete[] mpPixel;
    }
    mpPixel = new pixel_t[mNumberOfPixels];

    mpRawEventData->InitRead(mpRawRunHeader);

    int pix_first_sample;
//    int pix_last_sample;

    unsigned short* all_raw_data    = NULL;

    all_raw_data            = (unsigned short*) mpRawEventData->GetSamples();

    if (mVerbosityLvl > 1)
        cout << "...pixel progress: ";

    for ( int i = 0; i < mNumberOfPixels; i++ )
    {
        if (mVerbosityLvl > 1){
            if ( !(i%(mNumberOfPixels/10) ))
                cout << i/(mNumberOfPixels*1.0)*100 << "%.." ;
            if ( i == mNumberOfPixels-1)
                cout << "100% " ;
        }


        mpPixel[i].pedestal   = mpElectronicNoise[0][i].GetPedestal();
        mpPixel[i].SoftId     = mpRawEventData->GetPixelIds()[i];

        pix_first_sample        = i*mNumberOfSamples;

        // point beginning of pixel's rawdata array to address
        // of pixel's first sample's adress

        mpPixel[i].rawData      = &(all_raw_data[pix_first_sample]);

//        for (int j = pix_first_sample; j < pix_last_sample; j++)
//        {
//            cout << mpPixel[i].rawData[j] << " ";
//        }
//        cout << endl;

//        // this could be improved by
//         DONE BUT NOT TESTED
//        for (int j = i*mSampleSize; j < (i+1)*mSampleSize; j++)
//        {
//        }
    }

    if (mVerbosityLvl > 1)
        cout << endl;

    return;
}

void
MonteCarlo::ReadEvent(int Event)
{
    if (mVerbosityLvl > 0) cout << endl
                                << "====================" << endl
                                << "...reading Event: " << Event << endl
                                << "====================" << endl;


    //load certain event from tree
    mpEventTree->GetEntry(Event);

    ReadEventMetaData();
    ReadEventRawData();

    return;
}


void
MonteCarlo::WriteMc2Csv(TString filename)
{
    cout << "...writing mc to csv: " << filename << endl;
    cout << "...processing " << mNumberOfEntries << "Events" << endl;

    mCsvFileName = filename;
    OpenCsvFile(mCsvFileName);

    WriteFileHeader2Csv();
    WriteRunHeaderNames2Csv();
    WriteRunHeader2Csv();
    WriteEventHeaderNames2Csv();
    WriteEventDataNames2Csv();

//    for (int evt = 0; evt < mNumberOfEvents; evt++)
    for (int evt = 0; evt < mNumberOfEntries; evt++)
    {
        ReadEvent(evt);
        WriteEventHeader2Csv();
        WriteEventData2Csv();
    }

    cout << endl << "...conversion done " << endl;

    CloseCsvFile();
    return;
}


void
MonteCarlo::WritePixelData2Csv(int pixelID)
{
    if (mVerbosityLvl > 3) cout << "...writing pixel data to csv" << endl;
    mCsv << mEventNumber << mSeparator;
    mCsv << mpPixel[pixelID].SoftId << mSeparator;
    mCsv << mpPixel[pixelID].pedestal << mSeparator;

    for (int i = 0; i < mNumberOfSamples; i++)
    {
         mCsv << mpPixel[pixelID].rawData[i] << mSeparator;
    }
    mCsv << endl;

    return;
}

void
MonteCarlo::WriteEventData2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing event data to csv" << endl;

    for (int i = 0; i < mNumberOfPixels; i++)
    {
        WritePixelData2Csv(i);
    }

    return;
}

void
MonteCarlo::WriteEventDataNames2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing event categories to csv" << endl;

    mCsv << "### [RawData]" << endl;
    mCsv << "# mEventNumber" << mSeparator;
    mCsv << "pixelID" << mSeparator;
     mCsv << "pixelPedestal" << mSeparator;
    for (int i = 0; i < mNumberOfSamples; i++)
    {
        mCsv << "Raw_" << i << mSeparator;
    }
    mCsv << endl;

    return;
}

void
MonteCarlo::WriteEventHeaderNames2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing event header names to csv" << endl;

    mCsv << "### [EventHeader]" << endl;

    mCsv << "# mEventNumber" << mSeparator;
    mCsv << "mNumberOfBytes" << mSeparator;
    mCsv << "mIncidentAngle" << mSeparator;
    mCsv << "mPartId" << mSeparator;
    mCsv << "mEnergy" << mSeparator;
    mCsv << "mImpact" << mSeparator;
    mCsv << "mTelescopePhi" << mSeparator;
    mCsv << "mTelescopeTheta" << mSeparator;
    mCsv << "mPhi" << mSeparator;
    mCsv << "mTheta" << mSeparator;
    mCsv << "mCorsikaEventNumber" << mSeparator;
    mCsv << "mPhotElFromShower" << mSeparator;
    mCsv << "mEvtReuse" << mSeparator;
    mCsv << "mNumTriggerLvl1" << mSeparator;
    mCsv << "mFirstInteractionHeight" << mSeparator;
    mCsv << "mMomentumX" << mSeparator;
    mCsv << "mMomentumY" << mSeparator;
    mCsv << "mMomentumZ" << mSeparator;
    mCsv << "mZd" << mSeparator;
    mCsv << "mAz" << mSeparator;
    mCsv << "mX" << mSeparator;
    mCsv << "mY" << mSeparator;
//    mCsv << "mWeightedNumPhotons" ;
    mCsv << endl;

    return;
}


void
MonteCarlo::WriteEventHeader2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing event header to csv" << endl;

    mCsv << mEventNumber << mSeparator;
    mCsv << mNumberOfBytes << mSeparator;
    mCsv << mIncidentAngle << mSeparator;
    mCsv << mPartId << mSeparator;
    mCsv << mEnergy << mSeparator;
    mCsv << mImpact << mSeparator;
    mCsv << mTelescopePhi << mSeparator;
    mCsv << mTelescopeTheta << mSeparator;
    mCsv << mPhi << mSeparator;
    mCsv << mTheta << mSeparator;
    mCsv << mCorsikaEventNumber << mSeparator;
    mCsv << mPhotElFromShower << mSeparator;
    mCsv << mEvtReuse << mSeparator;
    mCsv << mNumTriggerLvl1 << mSeparator;
    mCsv << mFirstInteractionHeight << mSeparator;
    mCsv << mMomentumX << mSeparator;
    mCsv << mMomentumY << mSeparator;
    mCsv << mMomentumZ << mSeparator;
    mCsv << mZd << mSeparator;
    mCsv << mAz << mSeparator;
    mCsv << mX << mSeparator;
    mCsv << mY << mSeparator;
//    mCsv << mWeightedNumPhotons ;
    mCsv << endl;

    return;
}

void
MonteCarlo::WriteRunHeaderNames2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing run header names to csv" << endl;

    mCsv << "### [RunHeader]" << endl;

    mCsv << "# mNumberOfEntries" << mSeparator;
    mCsv << "mIntendedPulsePos" << mSeparator;
//    Csv << "mPedestalOffset" << mSeparator;
    mCsv << "mNumSimulatedShowers" << mSeparator;
    mCsv << "mNumberOfPixels" << mSeparator;
    mCsv << "mNumberOfSamples" << mSeparator;
//    mCsv << "mSampleSize" << mSeparator;
    mCsv << "mCamDist" << mSeparator;
    mCsv << "mSourceName" << mSeparator;
    mCsv << "mSlopeSpectrum" << mSeparator;
    mCsv << "mEnergyMin" << mSeparator;
    mCsv << "mEnergyMax" << mSeparator;
    mCsv << "mZdMin" << mSeparator;
    mCsv << "mZdMax" << mSeparator;
    mCsv << "mAzMin" << mSeparator;
    mCsv << "mAzMax" << mSeparator;
    mCsv << "mFileNumber" << mSeparator;
    mCsv << "mRunnumber" << mSeparator;
    mCsv << "mRunType" << endl;

    return;
}


void
MonteCarlo::WriteRunHeader2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing run header to csv" << endl;

    mCsv << mNumberOfEntries << mSeparator;
    mCsv << mIntendedPulsePos << mSeparator;
//    mCsv << mPedestalOffset << mSeparator;
    mCsv << mNumSimulatedShowers << mSeparator;
    mCsv << mNumberOfPixels << mSeparator;
    mCsv << mNumberOfSamples << mSeparator;
//    mCsv << mSampleSize << mSeparator;
    mCsv << mCamDist << mSeparator;
    mCsv << mSourceName << mSeparator;
    mCsv << mSlopeSpectrum << mSeparator;
    mCsv << mEnergyMin << mSeparator;
    mCsv << mEnergyMax << mSeparator;
    mCsv << mZdMin << mSeparator;
    mCsv << mZdMax << mSeparator;
    mCsv << mAzMin << mSeparator;
    mCsv << mAzMax << mSeparator;
    mCsv << mFileNumber << mSeparator;
    mCsv << mRunNumber << mSeparator;
    mCsv << mRunType << endl;

    return;
}

void
MonteCarlo::WriteFileHeader2Csv()
{
    if (mVerbosityLvl > 0) cout << "...writing file header to csv" << endl;

    mCsv << "### ==========================================================="
            << endl;
    mCsv << "### =             FACT Monte Carlo                            ="
            << endl;
    mCsv << "### ==========================================================="
            << endl;
    mCsv << "### = FileInfo:                                               "
            << endl;
    mCsv << "### =                                                         "
            << endl;
    mCsv << "### = Source Name:       " << mSourceName << endl;
    mCsv << "### = Number of Entries: " << mNumberOfEntries << endl;
    mCsv << "### = Run Number:        " << mRunNumber << endl;
    mCsv << "### = Run Type  :        " << mRunType << endl;
    mCsv << "### = File Number:       " << mFileNumber << endl;
    mCsv << "### ==========================================================="
         << endl ;
    mCsv << "###"
         << endl ;

    return;
}



//--------------------------------------------------------------------------------
//
//
//


//  FADC samples (hi gain) of all pixels
//  This is an array of Byte_t variables. The value of a FADC sample has a
//  size of n=fNumBytesPerSample bytes. Therefore, a FADC sample value will
//  occupy n consecutive elements of this array (little endian ordering, i.e,
//  less significant bits (and bytes) go first.
//  If m = GetNumHiGainSamples(), the n bytes corresponding to the value of the
//  i-th FADC sample of the j-th pixel are stored in the n consecutive
//  positions of this array, starting from fHiGainFadcSamples[j*n*m+i*n]


/*
int extractMC()
{

    TTree *tree = (TTree*)f.Get("Events");

//    MArrayB *array = NULL;

    const float dconv = 2000/4096;

    unsigned short *short_array;
    TH1F histo("histo", "histo", 150, 0, 150);

    TCanvas *canvas = new TCanvas("canvas", "canvas", 0, 0, 400, 400 );
    canvas->cd();
    tree->SetBranchAddress("MRawEvtData.", &branch1);

    tree->GetEntry(22);

    cout << "Branch1 " << branch1 << endl;

    array = branch1->GetSamples();


    for(int i = 0; i < 150; i++){
        cout << "i_" << i << "\t 1. " << short_array[i] ;
        short_array[i] = le16toh(short_array[i]);

//        uint16_t htobe16(uint16_t host_16bits);
//        uint16_t htole16(uint16_t host_16bits);
//        uint16_t be16toh(uint16_t big_endian_16bits);


        cout << "\t 2. " << short_array[i] ;
        cout << "\t arr. " << (int)array[2*i] << " " << (int)array[2*i+1] << endl;

        histo.SetBinContent(i, short_array[i]-1000);
    }
    histo.DrawCopy();
    canvas->Modified();
    canvas->Update();

//    int m = branch1->GetNumHiGainSamples();
//    int n = branch1->GetNumBytesPerSample();

    cout << "GetNumHiGainSamples " << branch1->GetNumHiGainSamples()
         << " GetNumBytesPerSample " << branch1->GetNumBytesPerSample() << endl;

//    tree->GetEntry(2);

//    array = branch1->GetSamples();

//    cout << "Branch2 " << branch1 << endl;

////    for(int i = 1439; i < 1441; i++){
////        cout << "Pixel: " << i << endl;

//        for(int j = 431950; j < 432020; j++){
//            cout << "Sl: " << j << "\t" ;
//            cout << (unsigned short) array[j] << "\t";
//            if ( j%5 == 0 ) cout <<endl;

//        }
//        cout << endl;
////    }

////    tree->SetBranchAddress("MRawEvtData.fHiGainFadcSamples.fArray", array);

}

*/