#include <stdio.h>
#include <math.h>
#include <iostream>
#include <fstream>
#include <gsl/gsl_integration.h>

using namespace std;

double f(double u, void * params){
double *par =  (double *) params;
double f; 
double A = 1-par[7]*u*u;
double B = par[7]*par[7]*u*u*u*u;
double C = A*A;
double D = pow ((C-B),0.5);//el momentum
double E = pow(u,2*par[3]+1);//el primer termino con A2
double N = -2*log(u);//logaritmo natural 
double G = pow(N,par[4]);//lo elevo a la gamm
double H = pow(u,2*par[3]-1);//el segundo termino con A2
if (u!=0.0){
       f = D*(4*par[1]*0.5*par[7]*E+par[6]*u*G*(2*par[9]-2*par[7])-2*par[1]*H)/
         ((par[8]*par[9])*(par[7]*u*u-1)*(u+pow(par[9],-0.5))*(1-2*par[7]*u*u+par[9]*u*u));
  }
  else f=0.0;
return f;
}

int main  () {
double mp = 0.938272;
double u = 0.0;
double s = 3.0;
double par [10] = {42.5, 25.5, 0.45, 0.565, 2.25, 30.0, 0.10, 2*mp*mp, mp, s};

gsl_integration_workspace * w = gsl_integration_workspace_alloc (1000);

int niter = 1000;
double liminf = 0.0;
double limsup = par[8]/par[7];
double errabs = 1e-11;
double errrel = 1e-11;
double result, error;

gsl_function F;
  F.function = &f;
  F.params = &par;

 for(int i=1; i<=1000; i++){
    s = s + 1;
    par[9] = s;
    gsl_integration_qawc (&F, liminf,limsup, sqrt(1/s), errabs, errrel , niter, w, &result, &error);
    cout<<s<<"   "<<result<<endl;

    

  }
  
  
  return 0;

}