#include <Likelihood.h>

//-- TOP scan configuration from "ScanConfig.cc" --//
//int _pmt_ch;
//int _tdc_ch_num;
//int _tdc_full_count;
//int _tdc_overflow_count;
//float _pdf_min;
////std::string _tree_name;
//std::vector<int> _veto_ch;
//std::vector<int> _t0_offset;


Likelihood::Likelihood( void ){ //{{{

   pdf = new TDC( _tdc_ch_num , _tdc_overflow_count );
   pdf->SetPdfMin( (float)_pdf_min );
   pdf->SetTdcOffset( _t0_offset );

   tr = new TTree( "L" , "Likelihood" );
   tr->Branch( "lnL" , &tr_lnL , "lnL/F" );
   tr->Branch( "nhit" , &tr_n , "nhit/I" );

   evt_min = -1;
   evt_max = -1;

   tdc_min = 0;
   tdc_max = _tdc_overflow_count;

   veto_ch_flag.resize( _tdc_ch_num );
   h_l = new TH1F();
   h_nhit = new TH1F();
   h_lch.resize( _tdc_ch_num );
   for( int ch = 0;ch < _tdc_ch_num;ch++){
      norm.push_back( 1. );
      veto_ch_flag[ch] = false;

      h_lch[ch] = new TH1F();
   }
   for( int i = 0;i < _veto_ch.size();i++)veto_ch_flag[ _veto_ch[i] ] = true;

   tree_name = "top";

} //}}}


void Likelihood::SetNormFactor( float n , int ch ){ //{{{
   if( ch < 0){ 
      for( int ch = 0;ch < _tdc_ch_num;ch++) norm[ch] = n;
   }else{
      norm[ch]=n;
   }
} //}}}


float Likelihood::GetSumPdf( int ch , int start_bin , int end_bin ){ //{{{ ch=-1 : all ch

   float sum = 0.;
   if( ch < 0){
      for(int c = 0;c < _tdc_ch_num;c++){
         for(int tdc = start_bin;tdc < end_bin;tdc++)sum += pdf->GetPdf( c , tdc ) * norm[c];
      }
   }else{
      for(int tdc = start_bin;tdc < end_bin;tdc++)sum += pdf->GetPdf( ch , tdc ) * norm[ch];
   }

   return sum;
} //}}}


void Likelihood::Clear( void ){ //{{{
   for( int ch = 0;ch < _tdc_ch_num;ch++)norm[ch] = 0.;
   pdf->Clear();
   rootfile->Clear();

   pdf->SetPdfMin( (float)_pdf_min );
   pdf->SetTdcOffset( _t0_offset );

   tr->Clear();

} //}}}


void Likelihood::Close( void ){ //{{{
   delete rootfile;
   delete tr;
   delete pdf;
   delete h_l;
   delete h_nhit;
}//}}}


void Likelihood::Write( std::string fname ){ //{{{

   TFile *f = new TFile( (TString)fname , "RECREATE" ); 

   tr->Write();

   f->Close();

}//}}}


void Likelihood::Calc( void ){ //{{{

   std::cout <<  "[Likelihood][Calc]Sum of PDF:" << GetSumPdf( -1 , 0 , 4000 ) << std::endl;

   //define the variable of TTree
   const int const_ch_num = _tdc_ch_num;
   //double d_tdc[const_ch_num];
   Double_t d_tdc[36];

   rootfile->cd();
   TTree *tr_root = (TTree*)rootfile->Get( tree_name.c_str() );
   tr_root->SetBranchAddress( "tdc" , &d_tdc );


   //TDC data loop start
   int lp_st;
   if( evt_min < 0 ){
      lp_st = 0;
   }else if( evt_max <= 1 ){
      lp_st = tr_root->GetEntries() * evt_min;
   }else{
      lp_st = evt_min;
   }

   int lp_end;
   if( evt_max < 0 ){
      lp_end = tr_root->GetEntries();
   }else if( evt_max <= 1 ){
      lp_end = tr_root->GetEntries() * evt_max;
   }else{
      lp_end = evt_max;
   }


   //Event loop
   for(Int_t i = lp_st;i < lp_end;i++){
      if(i % 100000 == 0)std::cerr << ".";
      tr_root->GetEntry( i );

      //read out simulation rootfile or not 
      float lnL = 0;
      int i_nhit = 0;

      for( Int_t ch = 0; ch < _tdc_ch_num; ch++ ){
         int tdc_count = (int)d_tdc[ch];

         tdc_count += _t0_offset[ch];

         if( tdc_count > tdc_min && tdc_count < tdc_max && !veto_ch_flag[ch] ){
            float p = pdf->GetPdf( ch , tdc_count );

            if( p > 0){
               lnL += log( p * norm[ch] );
               i_nhit++;

               h_lch[ch]->Fill( log( p * norm[ch] ) );
            }

         }
      }

      tr_lnL = lnL;
      tr_n = i_nhit;

      h_l->Fill( lnL );
      h_nhit->Fill( i_nhit );

      tr->Fill();

   }


} //}}}