//
// This is a c++ script that is meant to be compiled.
// Compilation can happen in two ways.
//
// 1) Run this with: root $AADIR/aa.C residuals.cc+  (the + is important)
//    This will compile residuals.cc on the fly and execute the function
//    'residuals'; i.e. ROOT will call the function that has the same name
//    as the file.
//    A main function is defined, but it is ignored in this case.
//    The script $AADIR/aa.C sets the include path for compilation to incldue
//    the relavant aanet directories.
// 2) This file can also be compiled into a standalone program. The make.py 
//    in this diretory will do that (producing the 'residuals' executable).
//    In this case the main() function is called.
//
// - see residuals.py for same example in python.
// = see residuals.C for a slightly different version that runs in the
//   CINT interpreter.              

#include "Evt.hh"
#include "Det.hh"
#include "EventFile.hh"
#include "TH1D.h"
#include "TCanvas.h"

void fill_residuals( vector<Hit>& hits, 
		     Trk& nu,
		     TH1D& hist, double offset = 0 )
{
  foreach( hit, hits ) 
    {
      double d = (hit.pos - nu.pos).len();
      double r = hit.t - ( nu.t + d / v_light ) - offset;
      hist.Fill(r);
    }
}

TH1D hres1("hres1","hres1;residual(ns)",125,-50,200);
TH1D hres2("hres2","hres2",125,-50,200);

void residuals()
{

  // nb: path relative to the examples directory
  EventFile f( "../data/mc5.1.numuCC.ps10.aa.root" );
  string detfile = "../data/KM3NeT_-00000001_20171212.detx";
  Det det(detfile);

  foreach( evt, f ) // in compiled c++, can use foreach macro (does not work in cint)
    {
      det.apply(evt);
      
      if (f.index > 300) break;
      if (evt.mc_hits.size() < 1000 )  continue;

      double offset = evt.mc_t - evt.t.AsDouble()*1e9;

      // shift tracks to match detector

      foreach( t, evt.mc_trks )
	{
	  t.pos += f.header.coord_origin();
	}

      Trk& nu = evt.mc_trks[0];
      
      fill_residuals( evt.mc_hits, nu, hres1 );
      fill_residuals( evt.hits,    nu, hres2 , offset);  
    }

  hres1.Scale( 1.0 / hres1.Integral() );
  hres2.Scale( 1.0 / hres2.Integral() );
  hres1.SetLineColor(2);
  hres2.SetLineColor(4);
  hres1.Draw();
  hres2.Draw("same"); 
}


// main() is ignored when running the compiled script via root
// but of course, when compiling into an exectuable, we need it.
// Note that no graphics are shown (you need to make a TApplication
// object if you want to show graphics in a compiled executable),
// but we save the plot to a png file.
int main()
{
  cout << "Hi there, this is your main() function talking " << endl;
  residuals();
  gPad->SaveAs("residuals.png");
}