/*
 * pxiData.h
 *
 *  Created on: Jan 14, 2013
 *      Author: yqiang
 */

#ifndef PXIDATA_H_
#define PXIDATA_H_

#include <TFile.h>
#include <TTree.h>
#include <TBranch.h>
#include <TH1.h>
#include <TString.h>
#include <TTimeStamp.h>
#include <TROOT.h>
#include <TCollection.h>
#include <TChain.h>
#include <stdio.h>

#define RECORDSIZE 10000
//#define RECORDSIZE 5000

const Int_t arraysize = RECORDSIZE;

class pxiData {

private:
	TChain *chain;
	// Declaration of leaf types
	Long64_t record_tsec;
	Long64_t record_tnsec;
	Float_t record_data[RECORDSIZE];
	// List of branches
	TBranch *b_record;

	virtual Int_t GetEntry(Int_t entry);

public:
	pxiData(TCollection *list, TString treename);
	virtual ~pxiData();
	TH1F *GetHistogram(TTimeStamp time_offset, Int_t length, Double_t freq);
	Double_t GetFreq(Int_t length);
	Int_t GetLength(Double_t &freq);
	TTimeStamp starttime;
	TTimeStamp endtime;
	Int_t entries;
};

pxiData::pxiData(TCollection *list, TString treename) {
	record_tsec = 0;
	record_tnsec = 0;
	// create a chain with treename
	chain = new TChain(treename);
	chain->AddFileInfoList(list, 100);
	if (!chain) {
		printf("Tree %s not found!\n", treename.Data());
		entries = 0;
		starttime = 0;
		endtime = 0;
		return;
	}
	chain->SetMakeClass(1);
	chain->SetBranchAddress("record", &record_tsec, &b_record);
	entries = chain->GetEntries();
	chain->GetEntry(0);
	starttime = TTimeStamp((Time_t) record_tsec, (Int_t) record_tnsec);
	chain->GetEntry(entries - 1);
	endtime = TTimeStamp((Time_t) record_tsec, (Int_t) record_tnsec);
	return;
}

pxiData::~pxiData() {
	if (!chain)
		return;
	delete chain;
	return;
}

Int_t pxiData::GetEntry(Int_t entry) {
	if (!chain)
		return 0;
	return chain->GetEntry(entry);
}

Double_t pxiData::GetFreq(Int_t length = 0) {
	Double_t freq = 0;
	// calculate report frequency, sampling frequency is rounded to 100 Hz
	Double_t t_diff = Double_t(endtime.GetSec() - starttime.GetSec())
			+ Double_t(endtime.GetNanoSec() - starttime.GetNanoSec()) * 1.e-9;
	if (length == 0) {
		length = entries;
		std::cout << "Period of time (s): " << t_diff << std::endl;
	}
	if (t_diff > 0)
		// round the report frequency to 1 Hz
		freq = floor(Double_t(length - 1) / t_diff + 0.5);
	return freq;
}

Int_t pxiData::GetLength(Double_t &freq) {
	Int_t length = 0;
	freq = 0;
	if (entries > 0) {
		Double_t t_diff;
		freq = this->GetFreq();
		length++;
		TTimeStamp itime = starttime;
		TTimeStamp ftime;
		for (int i = 1; i < entries; i++) {
			chain->GetEntry(i);
			ftime = TTimeStamp((Time_t) record_tsec, (Int_t) record_tnsec);
			t_diff = Double_t(ftime.GetSec() - itime.GetSec())
					+ Double_t(ftime.GetNanoSec() - itime.GetNanoSec()) * 1.e-9;
			length += Int_t(floor(t_diff * freq + 0.5));
			itime = ftime;
		}
		freq = this->GetFreq(length);
	}
	std::cout << "Expected data length: " << length << std::endl;
	std::cout << "Received data length: " << entries << std::endl;
	std::cout << "Report frequency (Hz): " << freq << std::endl;
	std::cout << "DAQ frequency (Hz): " << freq * arraysize << std::endl;
	return length;
}

TH1F *pxiData::GetHistogram(TTimeStamp time_offset, Int_t length,
		Double_t freq = 0.) {
	if (!chain)
		return 0;
	TString hname = chain->GetName();
	hname.Replace(0, hname.Last(':') + 1, 0, 0);
	TString htitle = TString("History plot of ") + hname
			+ TString(";Time (s);Value (V)");
	hname.Prepend("hist_");
	TH1F* hh = (TH1F*) gROOT->FindObject(hname);
	if (hh)
		delete hh; // deleted existing object

// get frequency if not passed
	if (freq == 0)
		freq = this->GetFreq();

	hh = new TH1F(hname, htitle, length * arraysize, 0, length / freq);

// fill histogram
	bool fillzero = true;
	TTimeStamp itime = time_offset;
	TTimeStamp ftime;
	Int_t index = 0;
	Int_t gap_count = 0;
	Int_t n_gap = 1;
	Int_t n_miss = 0;

	for (Int_t i = 0; i < length; i++) {
		if (gap_count < 1) {
			if (index >= entries)
				fillzero = true;
			else {
				chain->GetEntry(index);
				ftime = TTimeStamp((Time_t) record_tsec, (Int_t) record_tnsec);
				n_gap = Int_t(
						floor(
								(Double_t(ftime.GetSec() - itime.GetSec())
										+ Double_t(
												ftime.GetNanoSec()
														- itime.GetNanoSec())
												* 1.e-9) * freq + 0.5));
				// fill data
				if ((i == 0 && n_gap == 0) || (n_gap == 1)) {
					// if the first data is delayed by 1, fill zeros first
					if (i == 0 && n_gap == 1) {
						fillzero = true;
						index--;
					} else {
						fillzero = false;
						itime = ftime;
					}
				}
				// missing data start to count gap
				else if (n_gap > 1) {
					fillzero = true;
					std::cout << "Data missing from " << i << " for "
							<< n_gap - 1 << " readings" << std::endl;
					gap_count = 1;
					n_miss += n_gap - 1;
				}
				// data overlapping
				else {
					fillzero = true;
					i--;
					if (i >= 0)
						std::cout << "Data overlapping at " << i << std::endl;
				}
				index++;
			}
		} else {
			gap_count++;
			if (gap_count < n_gap)
				fillzero = true;
			else {
				fillzero = false;
				itime = ftime;
				gap_count = 0;
			}
//			}
		}

		// fill data
		if (fillzero)
			for (Int_t j = 0; j < arraysize; j++)
				hh->SetBinContent(i * arraysize + j + 1, 0);
		else
			for (Int_t j = 0; j < arraysize; j++)
				hh->SetBinContent(i * arraysize + j + 1, record_data[j]);
	}
	std::cout << hname << " was created";
	if (n_miss > 0)
		std::cout << ", missing " << n_miss << " readings";
	std::cout << std::endl;

	return hh;
}

#endif /* PXIDATA_H_ */