//Jacobian.C
Double_t SetMtetPhi(Double_t *x, Double_t *par)
{
  Float_t xx =x[0];
  Double_t norm = 3./4. * 1/(acos(-1)/2 +atan(par[0]/par[1])); 
  Double_t integrandS = par[0]*par[1] / (pow((par[4] - par[0]*par[0]),2) + pow(par[0]*par[1],2)) * par[3] / pow(par[4]*(par[4] - par[3]*par[3]),1/2);
  Double_t gamma = pow(par[2]*par[2] + par[0]*par[0],1/2) / par[0];
  Double_t alpha = pow(gamma*gamma - 1,1/2);
  Double_t mu = par[3]/par[0]; 
  Double_t sincarreT =  mu*mu * (mu*mu + alpha*alpha) * pow( mu*mu + mu*mu*alpha*alpha*cos(xx)*cos(xx) + alpha*alpha*sin(xx)*sin(xx),-1);
  Double_t I = (pow(mu,4) + mu*mu*pow(mu*alpha*cos(xx),2) + 2*pow(mu*alpha*sin(xx),2) + mu*mu*pow(mu*sin(xx),2)) * pow( mu*mu + alpha*alpha*sin(xx)*sin(xx),-1/2) * pow( mu*mu + mu*mu*alpha*alpha*cos(xx)*cos(xx) + alpha*alpha*sin(xx)*sin(xx),-3/2);
  Double_t f = 1/(2*acos(-1))* norm * integrandS * (2 - sincarreT)*I;
  return f;
}

Double_t IntegPhi(Double_t *x , Double_t *par)
{
  Double_t s = x[0];
  Double_t mt = par[0];
  TF1 *g1 = new TF1("g1",SetMtetPhi,0,8,5);
  g1->SetParameter(0,80);
  g1->SetParameter(1,2.14);
  g1->SetParameter(2,0);
  g1->SetParameter(3,mt);
  g1->SetParameter(4,s);
  ROOT::Math::GSLIntegrator ig2(1.E-6,1.E-6,1000); 
  ROOT::Math::WrappedTF1 wf2(*g1); 
  ig2.SetFunction(wf2);
  Double_t f = ig2.Integral(0,2*acos(-1));
  return f;
}

Double_t IntegS(Double_t *mt , Double_t *par)
{
  Double_t tampon = mt[0];
  TF1 *f1 = new TF1("f1",IntegPhi,0,10000000,1);
  f1->SetParameter(0,tampon);
  Double_t borne = tampon * tampon;
  ROOT::Math::GSLIntegrator ig(1.E-6,1.E-6,1000); 
  ROOT::Math::WrappedTF1 wf(*f1); 
  ig.SetFunction(wf);
  Double_t f = ig.IntegralUp(borne);
  return f;
}

Double_t fin(Double_t *x, Double_t *par)
{
  return SetMt(x,par) * IntegPhi(x,par);
}

void newIdea()
{ 
 Double_t mt=0;
 TF1 *final = new TF1("final",IntegS,0,100,0);
 TH2F *h1 =new TH2F("h1","h1",150,-0.5,149.5,1000,0,100);
 final->Draw();
}
  //  final2->Draw("same");
  
 //  int i=0;
//   while(i < 150)
//     {
//       //cout<<"entre dans la boucle"<<endl;
//       cout << "valeur de i "<< i<<endl;      
//       mt = i;
//       //cout << "valeur de mt "<< mt << endl;
//       Double_t f = final->Eval(mt);
//       //cout <<"valeur de f pincipal = " << f << endl;      
//       h1->Fill(i,f);
//       //cout << "fin d'un calcul" << endl;
//       i=i+1;
//       //mt = i;
//     }
//   h1->Draw();