#include <iostream>
#include <stdlib.h>
#include <list>
#include <set>

#include <TROOT.h>
#include <TFile.h>
#include <TTree.h>
#include <TGeoManager.h>
#include <TMath.h>

#include "IOARuntimeParameters.hxx"

#include "IReconFGDHacked.hxx"
#include "fgdUtils.hxx"

/// Initialization.  This method checks the DummyDatabase in
/// order to decide which methods to use for matching and 
/// for FGD isolated reconstruction.
IReconFGDHacked::IReconFGDHacked() : COMET::IAlgorithm("TFGDRecon") 
{

  matcher = new IFgdTpcNoMatcher();          // No matching

  // Choice of which iso-recon algorithm to use is controlled by 
  // ReconFGD parameters file.
  fgdIsoRecon = 0;
  if ((bool)COMET::IOARuntimeParameters::Get().GetParameterI("oaAlignTools.FGDIsoReconP0DCluster"))
    {
      fgdIsoRecon = new IFgdIsoReconP0DCluster();  // P0D Cluster iso recon
      std::cout<<"Using P0D clustering..."<<std::endl;
    }
  else if ((bool)COMET::IOARuntimeParameters::Get().GetParameterI("oaAlignTools.FGDIsoReconRECONSBCAT"))
    {
      fgdIsoRecon = new IFgdIsoReconSbcat();
      std::cout<<"Using RECONSBCAT..."<<std::endl;
    }
  else if ((bool)COMET::IOARuntimeParameters::Get().GetParameterI("oaAlignTools.FGDIsoReconRadon"))
    {
      fgdIsoRecon = new IFgdIsoReconRadon();
      std::cout<<"Using radon transform..."<<std::endl;
    }
  else
    {
      fgdIsoRecon = new IFgdIsoReconP0DCluster();
      std::cout<<"Using P0D clustering by default..."<<std::endl;
    }

  timeBinner = new IFGDTimeBinner();
  if (COMET::IOARuntimeParameters::Get().GetParameterI("ReconFGD.PID.Enable"))
    pid = new IFGDPid();
  else
    pid = 0;
}

IReconFGDHacked::~IReconFGDHacked(){
  delete timeBinner;
  delete fgdIsoRecon;
  delete matcher;
  if (pid)
    delete pid;
};

/// Do some finishing up analysis.
void IReconFGDHacked::FPlot(){    
  fgdIsoRecon->FinishUp();
}

void AddTimeBinBoundary(COMET::IHandle<COMET::IAlgorithmResult> Result, COMET::IHandle<COMET::IHitSelection> hits){

  // Get the interval boundaries
  std::pair<double,double> time_bin_boundary = fgdUtils::GetTimeBinBoundary(*hits);

  // Save boundaries using IRealDatum
  COMET::IRealDatum* binStartTime = new COMET::IRealDatum("binStartTime","binStartTime");
  binStartTime->SetValue(time_bin_boundary.first);
  Result->AddDatum(binStartTime);    

  COMET::IRealDatum* binEndTime = new COMET::IRealDatum("binEndTime","binEndTime");
  binEndTime->SetValue(time_bin_boundary.second);
  Result->AddDatum(binEndTime);    

}
/// This is the main FGD reconstruction program.
/// It basically calls two distinct reconstruction steps.  
///
/// -> Step 1 involves the matching of TPC tracks to FGD hits.  
///    There are several ways this can be done (including 
///    not being done at all).
/// -> Step 2 involves fitting the remaining FGD hits. 
///    Currently this is just implemented with a Hough transform
///    pattern recognition.  Call this "Iso" reconstruction.
///
COMET::IHandle<COMET::IAlgorithmResult> IReconFGDHacked::Process(const COMET::IAlgorithmResult& in){
  COMET::ICOMETEvent& event = GetEvent(); 
  COMET::IHandle<COMET::IHitSelection> fgd_raw = event.GetHitSelection("fgd");
    if (!fgd_raw) 
        return CreateResult(); // Return empty, but valid result.

    // time bin hits (binner does hit filtering too now!)
    TFGDHitBins bins = timeBinner->BinHits(fgd_raw);  

    COMET::IHandle<COMET::IAlgorithmResult> CumulativeResult(CreateResult());

    // Add TPC tracks

    int nbin = 0;
    for (TFGDHitBins::iterator bin = bins.begin(); bin != bins.end();
         bin++) {
        COMET::IAlgorithmResult input;
        input.AddHitSelection(new COMET::IHitSelection(**bin));
 
        // Do the matching reconstruction
        COMET::IHandle<COMET::IAlgorithmResult> MatcherResult = matcher->Process(input);  
        if(!MatcherResult) continue;
        
        // Add the full set of hits to result, so that they can be used 
        // in next step. (do we really need to do this?)
        MatcherResult->AddHitSelection(new COMET::IHitSelection(**bin));    
        
        // Do the fgd-iso reconstruction (this will always return a valid handle)
        COMET::IHandle<COMET::IAlgorithmResult> ReconResult = fgdIsoRecon->Process(*MatcherResult);
        
        // Copy results of TPC/FGD matching and iso-recon into final container.
        COMET::IHandle<COMET::IAlgorithmResult> FinalResult(CreateResult());
        fgdUtils::MergeResults(FinalResult, MatcherResult);
        fgdUtils::MergeResults(FinalResult, ReconResult);

  // Now that we have lists of used hits from both the isorecon
  // and matcher, merge them into a total ReconFGD hit used
  // list.

  COMET::IHandle<COMET::IHitSelection> matchedHits = MatcherResult->GetHitSelection("matchedFgdHits");
  COMET::IHandle<COMET::IHitSelection> isoreconedHits = ReconResult->GetHitSelection("usedFgdIsoreconHits");
  COMET::IHitSelection *usedHits = fgdUtils::combineHitSelections(matchedHits, isoreconedHits);
  usedHits->SetName("used");

  FinalResult->AddHitSelection(usedHits);
  
  // Make sure we have a container of unused hits, even if it's empty
  // This is a bit of a HACK
  COMET::IHandle<COMET::IHitSelection> unusedHits = FinalResult->GetHitSelection("unused");
  if (!unusedHits) 
    FinalResult->AddHitSelection(new COMET::IHitSelection("unused"));

        // Execute FGD iso-recon PID and save results
        if (pid) {
          COMET::IHandle<COMET::IAlgorithmResult> PidResult = pid->Process(*FinalResult);
          if (PidResult) {
            COMET::IHandle<COMET::IReconObjectContainer> fgdPid = 
              PidResult->GetResultsContainer("FGDPids");
            if (fgdPid)
              FinalResult->AddResultsContainer(new COMET::IReconObjectContainer(*fgdPid));
          }
        }

        // Save the information about the boundaries of this time bin
        AddTimeBinBoundary(FinalResult, *bin);
        
        
        // Merge the results from this time bin into global result
        fgdUtils::MergeResults(CumulativeResult, FinalResult);

        // Save the results for this particular time bin as separate algorithm result
        char result_name[100];
        sprintf(result_name,"ReconFGD_timebin%i",nbin);
        event.AddFit(FinalResult, result_name);
        nbin++;
    }
    return CumulativeResult;

}

/// This is the main FGD reconstruction program.
/// It basically calls two distinct reconstruction steps.  
///
/// -> Step 1 involves the matching of TPC tracks to FGD hits.  
///    There are several ways this can be done (including 
///    not being done at all).
/// -> Step 2 involves fitting the remaining FGD hits. 
///    Currently this is just implemented with a Hough transform
///    pattern recognition.  Call this "Iso" reconstruction.
///
COMET::IHandle<COMET::IAlgorithmResult> IReconFGDHacked::HackedProcess(COMET::IHandle<COMET::IHitSelection> fgd_raw){
  COMET::ICOMETEvent& event = GetEvent(); 
    if (!fgd_raw) 
        return CreateResult(); // Return empty, but valid result.

    // time bin hits (binner does hit filtering too now!)
    TFGDHitBins bins = timeBinner->BinHits(fgd_raw);  

    COMET::IHandle<COMET::IAlgorithmResult> CumulativeResult(CreateResult());

    // Add TPC tracks

    int nbin = 0;
    for (TFGDHitBins::iterator bin = bins.begin(); bin != bins.end();
         bin++) {
        COMET::IAlgorithmResult input;
        input.AddHitSelection(new COMET::IHitSelection(**bin));
 
        // Do the matching reconstruction
        COMET::IHandle<COMET::IAlgorithmResult> MatcherResult = matcher->Process(input);  
        if(!MatcherResult) continue;
        
        // Add the full set of hits to result, so that they can be used 
        // in next step. (do we really need to do this?)
        MatcherResult->AddHitSelection(new COMET::IHitSelection(**bin));    
        
        // Do the fgd-iso reconstruction (this will always return a valid handle)

        COMET::IHandle<COMET::IAlgorithmResult> ReconResult = fgdIsoRecon->Process(*MatcherResult);
        
        // Copy results of TPC/FGD matching and iso-recon into final container.
        COMET::IHandle<COMET::IAlgorithmResult> FinalResult(CreateResult());
        fgdUtils::MergeResults(FinalResult, MatcherResult);
        fgdUtils::MergeResults(FinalResult, ReconResult);

  // Now that we have lists of used hits from both the isorecon
  // and matcher, merge them into a total ReconFGD hit used
  // list.

  COMET::IHandle<COMET::IHitSelection> matchedHits = MatcherResult->GetHitSelection("matchedFgdHits");
  COMET::IHandle<COMET::IHitSelection> isoreconedHits = ReconResult->GetHitSelection("usedFgdIsoreconHits");
  COMET::IHitSelection *usedHits = fgdUtils::combineHitSelections(matchedHits, isoreconedHits);
  usedHits->SetName("used");

  FinalResult->AddHitSelection(usedHits);
  
  // Make sure we have a container of unused hits, even if it's empty
  // This is a bit of a HACK
  COMET::IHandle<COMET::IHitSelection> unusedHits = FinalResult->GetHitSelection("unused");
  if (!unusedHits) 
    FinalResult->AddHitSelection(new COMET::IHitSelection("unused"));

        // Execute FGD iso-recon PID and save results
        if (pid) {
          COMET::IHandle<COMET::IAlgorithmResult> PidResult = pid->Process(*FinalResult);
          if (PidResult) {
            COMET::IHandle<COMET::IReconObjectContainer> fgdPid = 
              PidResult->GetResultsContainer("FGDPids");
            if (fgdPid)
              FinalResult->AddResultsContainer(new COMET::IReconObjectContainer(*fgdPid));
          }
        }

        // Save the information about the boundaries of this time bin
        AddTimeBinBoundary(FinalResult, *bin);
        
  // Merge the results from this time bin into global result
        fgdUtils::MergeResults(CumulativeResult, FinalResult);

        // Save the results for this particular time bin as separate algorithm result
        char result_name[100];
        sprintf(result_name,"ReconFGD_timebin%i",nbin);
        event.AddFit(FinalResult, result_name);
        nbin++;
    }
    return CumulativeResult;

}