const double W2Min = 4.0;
const double MpSq=0.9383*0.9383;
const int nxB=3, nQ2=4;
const double QSqBinLow[nQ2+1]={1.5,2.5,4.0 ,6.0 ,10.0};
const int nsimu = 10000;
double LogLi[nQ2][nxB];

double xBjMin(double Q2){
    double result = Q2*0.75/12.0;
    return result;
}
double xBjMax(double Q2){
    double result = 1.0;
    result = Q2/(W2Min+Q2-MpSq);
    return result;
}


Double_t Likelihood(Double_t *val,Double_t *par){
    Float_t x = val[0];
    Float_t y = val[1];
    Double_t result=0;
    Double_t func=0;
    if (x>0.0) {
        func= pow(x,par[0]) * pow((1-x),par[1])*pow((1+ (y*y)/ par[3]),-par[2]);
        result+=func;
    }
    return result;
}

Double_t xBQ2binning_tree::Loop(double norm_j_i, int jQ2Bin, int ixBin, TF2 *TFLikeli)
{
    double logli=0.0;
    int nbin= 0;
    static TH2D * h2_Q2_vs_xB = new TH2D("h2_Q2_vs_xB","Q^{2} vs x_{B} ; x_{B}  ; Q^{2} (GeV^{2})  ",50,0.0,1.0, 50,0.0,10.0);
    
    
    
    
//       To read only selected branches, Insert statements like:
// METHOD1:
//    fChain->SetBranchStatus("*",0);  // disable all branches
//    fChain->SetBranchStatus("branchname",1);  // activate branchname
// METHOD2: replace line
//    fChain->GetEntry(jentry);       //read all branches
//by  b_branchname->GetEntry(ientry); //read only this branch
   if (fChain == 0) return 0;
    Long64_t nentries = fChain->GetEntriesFast();
    
    Long64_t nbytes = 0, nb = 0;
    for (Long64_t jentry=0; jentry<nentries;jentry++) {
        Long64_t ientry = LoadTree(jentry);
        if (ientry < 0) break;
        nb = fChain->GetEntry(jentry);   nbytes += nb;
   
      
      // Select events from bin (ixBin, jQ2Bin)
       if ((ixB==ixBin)&&(jQ2==jQ2Bin)) {
           logli += log(TFLikeli->Eval(xBj,QSq))-log(norm_j_i);
           nbin++;
       }
   }
    logli *=1.0/nbin;
    printf("LogLi[jQ2=%02d,ixB=%02d]= %10.3f \n", jQ2Bin, ixBin, logli);
    
    
    return logli;
    
}