#include <iostream>
void plot_sigma(){
  ifstream file("eta_xsec.dat");

  gStyle->SetOptStat(0);

  double t[10000];
  double costheta_cm[10000];
  double theta_lab[10000];
  double dsigma_dt[10000];
  double dsigma_dOmega[10000];
  double Eg=0.;
  double mass=0.;

  file >> Eg;
  file >> mass;
  for (unsigned int i=0;i<10000;i++){
   file >> t[i];
   file >> costheta_cm[i];
   file >> theta_lab[i];
   file >> dsigma_dt[i];
   file >> dsigma_dOmega[i];
  }
  file.close();

  TCanvas *c1 = new TCanvas("c1","",600,500);

  char mytitle[80];
  sprintf(mytitle,"#eta(%f) cross section for E_{#gamma}=%f GeV",mass,Eg);
  TH2F *hr13 = new TH2F("hr13",mytitle,
                        10000,0,2.0,1000,0,0.5);
  hr13->SetXTitle("-t [GeV^{2}]");
  hr13->GetYaxis()->SetTitleOffset(1.25);
  hr13->SetYTitle("d#sigma/dt [#mubarn/GeV^{2}]");
  hr13->Draw();
  
  double sum=0.;
  double myt,t_old=t[0];
  for (unsigned int i=0;i<10000;i++){
    myt=t[i];
    //    if (fabs(myt)>0.09 && fabs(myt)<0.7)
      {
      sum+=dsigma_dt[i]*(myt-t_old);
    }
    t_old=myt;
  }
  cout << "xsec " << sum << " micro-barns" <<endl;

  TGraph *g0=new TGraph(10000,t,dsigma_dt);
  g0->SetMarkerStyle(1);
  g0->Draw("PC");
  
  if (Eg>5.9&&Eg<6.1){
    double t6gev[12]={0.012,0.03,0.09,0.18,0.31,0.39,0.48,0.58,0.68,0.79,0.91,1.18};
    double sig6gev[12]={0.191,0.158,0.232,0.194,0.133,0.142,0.106,0.091,0.055,0.048,0.037,0.024};
    double err6gev[12]={0.037,0.039,0.042,0.034,0.028,0.028,0.026,0.023,0.023,0.019,0.018,0.015};
    
    TGraphErrors *g2=new TGraphErrors(12,t6gev,sig6gev,0,err6gev);
    g2->SetMarkerStyle(20);
    g2->Draw("p");
  }
  else if (Eg>3.9&&Eg<4.1){ 
    double t4gev[12]={0.017,0.028,0.063,0.120,0.21,0.31,0.44,0.59,0.76,0.94,1.12,1.37};
    double sig4gev[12]={0.295,0.270,0.406,0.408,0.432,0.330,0.234,0.153,0.128,0.062,0.041,0.044};
    double err4gev[12]={0.05,0.049,0.062,0.058,0.065,0.050,0.042,0.034,0.030,0.022,0.024,0.022};

 

    TGraphErrors *g2=new TGraphErrors(12,t4gev,sig4gev,0,err4gev);
    g2->SetMarkerStyle(20);
    g2->Draw("p");
  }
  else if (Eg>7.9 && Eg<8.1){
    double t8gev[5]={0.3,0.4,0.5,0.65,0.8};
    double sig8gev[5]={0.0692,0.0612,0.0608,0.0372,0.0258};
    double err8gev[5]={0.0084,0.0076,0.0074,0.0046,0.0035};
    
    
    TGraphErrors *g2=new TGraphErrors(5,t8gev,sig8gev,0,err8gev);
    g2->SetMarkerStyle(20);
    g2->Draw("p");
  }

  TCanvas *c2 = new TCanvas("c2","",600,500);

  TH2F *hr14 = new TH2F("hr14",mytitle,
                        100,0,10.0,100,0,0.3);
  hr14->SetXTitle("Lab #theta [degrees]");
  hr14->GetYaxis()->SetTitleOffset(1.25);
  hr14->SetYTitle("d#sigma/d#Omega [#mubarn/steradian]");
  hr14->Draw();
  
  for (unsigned int i=0;i<10000;i++){ theta_lab[i]*=57.3;}

  TGraph *g3=new TGraph(10000,theta_lab,dsigma_dOmega);
  g3->SetMarkerStyle(1);
  g3->Draw("PC");

  double theta_max=57.3*atan2(118.,1200.);
  TLine *line=new TLine(theta_max,0,theta_max,0.3);
  line->SetLineColor(2);
  line->SetLineWidth(2);
  line->SetLineStyle(2);
  line->Draw();
}