// $Id$
//
//    File: JEventProcessor_bcal_single_photon.cc
// Created: Fri Jul 13 15:47:34 EDT 2012
// Creator: davidl (on Linux ifarm1102 2.6.18-274.3.1.el5 x86_64)
//

#include <iostream>
#include <set>
using namespace std;

#include <TRACKING/DMCThrown.h>
#include <BCAL/DBCALShower.h>
#include <BCAL/DBCALCluster.h>
#include <BCAL/DBCALTDCHit.h>
#include <BCAL/DBCALPoint.h>

#include "JEventProcessor_bcal_single_photon.h"
using namespace jana;



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


//------------------
// JEventProcessor_bcal_single_photon (Constructor)
//------------------
JEventProcessor_bcal_single_photon::JEventProcessor_bcal_single_photon()
{

}

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

}

//------------------
// init
//------------------
jerror_t JEventProcessor_bcal_single_photon::init(void)
{
	ratio_Erecon_Ethrown = new TH1D("ratio_Ethrown_Erecon", "#frac{E_{thrown}}{E_{recon}}", 10000, 0.0, 2.0E-3);

	Ethrown_vs_Erecon = new TH2D("Ethrown_vs_Erecon", "", 550, 0.0, 5.5, 225, 0.0, 4.5);

	tdc_match = new TTree("tdc_match", "TDC hits match to ADC points");
	tdc_match->Branch("M", &match, "module/I:layer:sector:end:E/F:t_adc:t_tdc:t_mean_adc");

	return NOERROR;
}

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

//------------------
// evnt
//------------------
jerror_t JEventProcessor_bcal_single_photon::evnt(JEventLoop *loop, int eventnumber)
{
	vector<const DMCThrown*> mcthrowns;
	vector<const DBCALShower*> bcalshowers;
	vector<const DBCALTDCHit*> bcaltdchits;
	
	loop->Get(mcthrowns);
	loop->Get(bcalshowers);
	loop->Get(bcaltdchits);

	// Filter out single block showers
	vector<const DBCALShower*> multiblock_bcalshowers;
	for(unsigned int i=0; i<bcalshowers.size(); i++){
		if(bcalshowers[i]->N_cell > 1) multiblock_bcalshowers.push_back(bcalshowers[i]);
	}
	
	if(mcthrowns.size()!=1 || multiblock_bcalshowers.size()!=1)return NOERROR;
	
	double Ethrown = mcthrowns[0]->energy();
	double Erecon = multiblock_bcalshowers[0]->E;
	double ratio = Ethrown/Erecon;
	ratio_Erecon_Ethrown->Fill(ratio);
	
	Ethrown_vs_Erecon->Fill(Erecon, Ethrown);
	
	// Get DBCALPoint objects associated with this reconstructed shower
	// and match them to DBCALTDCHit objects. The DBCALShower object can
	// (in principle) be made of multiple DBCALCluster objects, each with
	// its own list of DBCALPoint objects. (At this point in time, the
	// IU algorithm has a 1 to 1 correspondence between clusters and showers
	// but we build this to handle the more complex case). We want a single
	// list of all DBCALPoint objects from all clusters. I suppose that
	// it's possible, that the clustering algorithm could use the same
	// DBCALPoint in multiple DBCALClusters so we want to make sure it
	// does not show up twice. Hence, the use of an STL set.
	vector<const DBCALCluster*> bcalclusters;
	multiblock_bcalshowers[0]->Get(bcalclusters);
	
	set<const DBCALPoint*> bcalpoints;
	for(unsigned int i=0; i<bcalclusters.size(); i++){
		vector<const DBCALPoint*> tmp;
		bcalclusters[i]->Get(tmp);
		for(unsigned int j=0; j<tmp.size(); j++) bcalpoints.insert(tmp[j]);
	}

	// Loop over DBCALPoint objects and match with DBCALTDCHit objects.
	//tdc_match_t &match = *matchptr;
	set<const DBCALPoint*>::iterator iter = bcalpoints.begin();
	for(; iter!=bcalpoints.end(); iter++){
		const DBCALPoint *adchit = *iter;
		
		// Get DBCALHit for upstream and downstream ends.
		// Note that it's possible there is only one end hit!
		const DBCALHit *uphit=NULL;
		const DBCALHit *dnhit=NULL;
		vector<const DBCALHit*> bcalhits;
		adchit->Get(bcalhits);
		for(unsigned int i=0; i<bcalhits.size(); i++){
			if(bcalhits[i]->end == DBCALGeometry::kUpstream)   uphit = bcalhits[i];
			if(bcalhits[i]->end == DBCALGeometry::kDownstream) dnhit = bcalhits[i];
		}

		for(unsigned int i=0; i<bcaltdchits.size(); i++){
			const DBCALTDCHit *tdchit = bcaltdchits[i];
			
			if(adchit->module() != tdchit->module)continue;
			if(adchit->layer()  != tdchit->layer )continue;
			if(adchit->sector() != tdchit->sector)continue;

			match.module = tdchit->module;
			match.layer = tdchit->layer;
			match.sector = tdchit->sector;
			match.end = tdchit->end;
			match.E = adchit->E();
			match.t_adc = -1000.0; // flag that no hit was found (overwritten below)
			match.t_tdc = tdchit->t;
			match.t_mean_adc = adchit->t();
			
			if(uphit!=NULL && tdchit->end==DBCALGeometry::kUpstream) match.t_adc = uphit->t;
			if(dnhit!=NULL && tdchit->end==DBCALGeometry::kDownstream) match.t_adc = dnhit->t;

			tdc_match->Fill();
		}
	}
	

	return NOERROR;
}

//------------------
// erun
//------------------
jerror_t JEventProcessor_bcal_single_photon::erun(void)
{
	// This is called whenever the run number changes, before it is
	// changed to give you a chance to clean up before processing
	// events from the next run number.
	return NOERROR;
}

//------------------
// fini
//------------------
jerror_t JEventProcessor_bcal_single_photon::fini(void)
{
	// Called before program exit after event processing is finished.
	return NOERROR;
}