// $Id$
//
//    File: JEventProcessor_ps_rate.cc
// Created: Wed Mar  8 21:07:40 EST 2017
// Creator: davidl (on Linux gluon112.jlab.org 2.6.32-642.3.1.el6.x86_64 x86_64)
//

#include <iostream>
using namespace std;

#include "JEventProcessor_ps_rate.h"
using namespace jana;


#include <DAQ/DBeamCurrent.h>
#include <PAIR_SPECTROMETER/DPSHit.h>
#include <PAIR_SPECTROMETER/DPSCHit.h>




// Routine used to create our JEventProcessor
#include <JANA/JApplication.h>
#include <JANA/JFactory.h>
extern "C"{
void InitPlugin(JApplication *app){
	InitJANAPlugin(app);
	app->AddProcessor(new JEventProcessor_ps_rate());
}
} // "C"


//------------------
// JEventProcessor_ps_rate (Constructor)
//------------------
JEventProcessor_ps_rate::JEventProcessor_ps_rate()
{
	// This tells how many event times to store before
	// we start taking differences. This accomodates
	// the out-of-order nature of the multi-threaded
	// process.
	max_event_times = 100;

}

//------------------
// ~JEventProcessor_ps_rate (Destructor)
//------------------
JEventProcessor_ps_rate::~JEventProcessor_ps_rate()
{

}

//------------------
// init
//------------------
jerror_t JEventProcessor_ps_rate::init(void)
{
	pshit_tdiff   = new TH1D("pshit_tdiff"  ,  "Time diff. between DPSHit of different arms"    , 400, -50.0, 50.0);
	pschit_tdiff  = new TH1D("pschit_tdiff",   "Time diff. between DPSCHit of different arms"   , 400, -50.0, 50.0);
	psc_ps_tdiff  = new TH1D("psc_ps_tdiff" ,  "Time diff. between DPSHit pair and DPSCHit pair", 400, -50.0, 50.0);
	ps_teventdiff = new TH2D("ps_teventdiff" , "Time diff. between events with PS coinc;time from start of run (sec);tdiff (msec).", 10000, 0.0, 10000.0, 500, 0.0, 5.0);
	Nps_coinc     = new TH1D("Nps_coinc" ,     "Number of PS coincidences per second", 10000, 0.0, 10000.0);
	Nps_coinc_events = new TH1D("Nps_coinc_events" ,     "Number of events with one or more PS coincidences per second", 10000, 0.0, 10000.0);
	tfiducial     = new TH1D("tfiducial" ,     "Map of fiducial time regions", 10000, 0.0, 10000.0);

	return NOERROR;
}

//------------------
// brun
//------------------
jerror_t JEventProcessor_ps_rate::brun(JEventLoop *eventLoop, int32_t runnumber)
{
	// This is called whenever the run number changes
	return NOERROR;
}

//------------------
// evnt
//------------------
jerror_t JEventProcessor_ps_rate::evnt(JEventLoop *loop, uint64_t eventnumber)
{
	vector<const DPSHit*> pshits;
	vector<const DPSCHit*> pschits;
	vector<const DBeamCurrent*> beamcurrents;
	loop->Get(pshits);
	loop->Get(pschits);
	loop->Get(beamcurrents);
	
	// PS coincidences
	vector<double> psdiffs;
	vector<double> psmeans;
	for(uint32_t i=0; i<pshits.size(); i++){
		auto h1 = pshits[i];
		for(uint32_t j=i+1; j<pshits.size(); j++){
			auto h2 = pshits[j];
			if(h1->arm == h2->arm) continue; // choose opposite arms
			
			psdiffs.push_back( h1->t - h2->t);
			psmeans.push_back((h1->t + h2->t)/2.0);
		}
	}

	// PSC coincidences
	vector<double> pscdiffs;
	vector<double> pscmeans;
	for(uint32_t i=0; i<pschits.size(); i++){
		auto h1 = pschits[i];
		for(uint32_t j=i+1; j<pschits.size(); j++){
			auto h2 = pschits[j];
			if(h1->arm == h2->arm) continue; // choose opposite arms
			
			pscdiffs.push_back( h1->t - h2->t);
			pscmeans.push_back((h1->t + h2->t)/2.0);
		}
	}
	
	// Event time
	double tevent = -1.0E6;
	bool is_fiducial = false;
	if(!beamcurrents.empty()){
		tevent = beamcurrents[0]->t;
		is_fiducial = beamcurrents[0]->is_fiducial;
	}

	// Lock ROOT mutex and fill histograms
	japp->RootFillLock(this);

	if(is_fiducial){
		uint32_t ibin = (uint32_t)(floor(tevent)+1.0);
		tfiducial->SetBinContent(ibin, 1.0);
	}

	for(auto tdiff : psdiffs ) pshit_tdiff->Fill(tdiff);
	for(auto tdiff : pscdiffs) pschit_tdiff->Fill(tdiff);
	
	bool pscoinc = false;
	for(auto tps : psmeans ){
		for(auto tpsc : pscmeans){
			psc_ps_tdiff->Fill(tps - tpsc);
			if( fabs(tps-tpsc) < 5.0 ){
				pscoinc = true;
				Nps_coinc->Fill(tevent);
			}
		}
	}
	
	if(pscoinc) Nps_coinc_events->Fill(tevent);
	
	if(pscoinc && tevent > -1.0E3){
		event_times.insert(tevent);
		
		while(event_times.size() > max_event_times){
			auto it = event_times.begin();;
			double t1 = *it++;
			double t2 = *it;
			
			if( floor(t1) == floor(t2) ){
				ps_teventdiff->Fill(floor(t1), (t2-t1)*1.0E3);
			}
			
			event_times.erase(event_times.begin());
		}
	}

	japp->RootFillUnLock(this);


	return NOERROR;
}

//------------------
// erun
//------------------
jerror_t JEventProcessor_ps_rate::erun(void)
{
	return NOERROR;
}

//------------------
// fini
//------------------
jerror_t JEventProcessor_ps_rate::fini(void)
{
	if(event_times.size() > 1 ){
		auto it = event_times.begin();
		while( true ){
			double t1 = *it++;
			if( it == event_times.end() ) break;
			double t2 = *it;
			
			if( floor(t1) == floor(t2) ){
				ps_teventdiff->Fill(floor(t1), (t2-t1)*1.0E3);
			}
		}
	}

	return NOERROR;
}