#include "Trk.hh"
#include "Hit.hh"
#include "Vec.hh"
#include "Evt.hh"
#include "EventFile.hh"
#include "JppShowerPdf.hh"
#include "JppShowerFit.hh"
#include "JppShowerPdf_simple.hh"
#include "JppShowerFit_simple.hh"
#include "rec2/CascadePdf.hh"
#include "Det.hh"
#include "rec2/JppGrad.hh"


int main( int argc, int** argv ){

  EventFile f( "/sps/km3net/users/jseneca/aanet/data/km3_v5/mcv5.1.genhen_nueCC.km3_AAv1.jte.jchain.aashower.94.root" );
  f.set_index( 2147 );
  Det det( "/sps/km3net/users/jseneca/aanet/data/km3_v5/KM3NeT_-00000001_20171212.detx" );
  
  Trk trk = f.evt.mc_trks[0];
  
  vector<Hit> pre_hits = f.evt.hits;
  det.apply( pre_hits );
  vector<Hit> hits;
  
  for( auto h = pre_hits.begin(); h != pre_hits.end(); h++ ){
    Hit hit = *h;
    hit.t = hit.t - ( f.evt.mc_t - f.evt.t.AsDouble()*1e9 );
    double D = ( trk.pos - hit.pos ).len();
    double dt = abs( hit.t - D/v_light );
    if( D < 500 && dt < 100 ){
      hits.push_back( hit );
    }
  };
  
  cout << "len hits: " << hits.size() << endl;
  
  assert( hits.size() );
  
  CascadePdf cascade_pdf( "/pbs/throng/km3net/src/Jpp/master/data/J14p.dat",
			  "/pbs/throng/km3net/src/Jpp/master/data/J13p.dat", 
			  0, 
			  0, //0 : do nothing, 1 = move to shower max, > 1 : sample
			  5000, //buggy behavior when no background
			  -100,
			  900,
			  0.0 );
   
  print( "Made cascadepdf" );
  
  JppShower_simple gradpdf( cascade_pdf, 20 );
  
  print( "Made JppShower_simple" );
  
  gradpdf.init_event( det, hits, trk.pos, 300, true);
  JppShowerFit_simple fitter_gradpdf( "5d_gradpdf",
				      &gradpdf,
				      "Minuit2",
				      "MIGRAD",
				      1 );
  
  print( "Made JppShowerFit_simple" );
  
  JppShower pdf( "/pbs/throng/km3net/src/Jpp/master/data/J14p.dat", //Scattered PDF must be first.
		 "/pbs/throng/km3net/src/Jpp/master/data/J13p.dat",
		 0.0, // blurring of the PDF
		 -100,	
		 900,
		 false,
		 20,
		 7e3,
		 0);

  pdf.init_event( det, hits, trk.pos, 300, true);
  JppShowerFit fitter_pdf( "5d_pdf",
			   &pdf,
			   "Minuit2",
			   "MIGRAD",
			   1 );

 print("JppShowerFit ok");
 
 aa_legacy_timer::Timer t1( " GRAD PDF " );
 Trk fit_gradpdf = fitter_gradpdf.fit( trk, det , hits, true, true );
 aa_legacy_timer::Timer t0( " PDF " );
 Trk fit_pdf = fitter_pdf.fit( trk, det , hits, true, true );

 print( "MC neutrino: ", trk );
 print( "Fit pdf:     ", fit_pdf );
 print( "Fit gradpdf: ", fit_gradpdf );
 
 t1.report();
};