void calc_phi_t_comp(int dataset,double Eg){
  gStyle->SetOptStat(0);
  double eff[50],deff[50];
  double y[50],dy[50];  
  double eff2[50],deff2[50];
  double y2[50],dy2[50];
  double t[50],dt[50];
  double Nps[42];

  int bin=int((Eg-2.8)/0.2);
    //  int bin=int((Eg-2.8)/0.2);
  TFile *f=NULL;
  switch(dataset){
  case 0:
    f=TFile::Open("/work/halld/home/staylor/recon/xsecs/flux_51384_51457.root");
    cout << "Flux file:  flux_51384_51457.root" << endl;
    break;
  case 1:
    f=TFile::Open("/work/halld/home/staylor/recon/xsecs/flux_40856_41105.root");
    cout << "Flux file:  flux_40856_41105.root" << endl;
   break;
  case 2:
    f=TFile::Open("/work/halld/home/staylor/recon/xsecs/flux_41106_41257.root");
    cout << "Flux file:  flux_41106_41257.root" << endl;
    break;
  default:
    cout << "No data set specified" << endl;
    break;
  }
  if (f!=NULL){  
    f->cd();
    TH1D *tagged_flux=(TH1D*)f->FindObjectAny("tagged_flux");
    double Nps_sum=tagged_flux->GetBinContent(bin+1);
    f->Close();
    
    ifstream in3("phi_t_yields.dat");
    for (int i=0;i<bin;i++){
      string skipline;
      getline(in3,skipline);
      //   cout << skipline << endl;
      
    }
    double dummy;
    in3 >> dummy;
    for (int i=0;i<50;i++){
      in3 >> y[i];
      in3 >> dy[i];
    }
    in3.close();
    ifstream in4("phi_kk_t_acc.dat"); 
    for (int i=0;i<bin;i++){
      string skipline;
      getline(in4,skipline);
      //   cout << skipline << endl;
      
    }
    in4 >> dummy;
    for (int i=0;i<50;i++){
      in4 >> eff[i];
      //eff[i]=0.2;
      in4 >> deff[i];
    }
    in4.close();
    ifstream in1("phi_t_yields2.dat");
    for (int i=0;i<bin;i++){
      string skipline;
      getline(in1,skipline);
      //   cout << skipline << endl;
      
    }
    in1 >> dummy;
    for (int i=0;i<50;i++){
      in1 >> y2[i];
      in1 >> dy2[i];
    }
    in1.close();
    ifstream in2("phi_kk_t_acc2.dat"); 
    for (int i=0;i<bin;i++){
      string skipline;
      getline(in2,skipline);
      //   cout << skipline << endl;
      
    }
    in2 >> dummy;
    for (int i=0;i<50;i++){
      in2 >> eff2[i];
      //	eff2[i]=0.5;
      in2 >> deff2[i];
    }
    in2.close();
    
    unsigned int ind=0;
    double m_phi=1.01946;
    double m_phi_sq=m_phi*m_phi;
    double m_p=0.93827;
    double m_p_sq=m_p*m_p;
    ifstream inclas("/work/halld/home/staylor/recon/xsecs/2018L/clasdb_E63M1.txt");
    double Eclas_old=0.,Eclas=0.;
    vector<double>Eclas_vec;
    vector<vector<double> >xclas_vec;
    vector<vector<double> >dxclas_vec;
    vector<vector<double> >tclas_vec;
    double w_old=0.,w=0.,costheta,xclas,dxclas;
    vector<double>tempt;
    vector<double>tempx;
    vector<double>tempdx;
    for (int i=0;i<2040;i++){
      inclas >> w;
      inclas >> costheta;
      inclas >> xclas;
      inclas >> dxclas;
      
      if (fabs(w_old-w)>0.001){
	xclas_vec.push_back(tempx);
	dxclas_vec.push_back(tempdx);
	tclas_vec.push_back(tempt);
	Eclas_vec.push_back(Eclas_old);
	
	tempx.clear();
	tempdx.clear();
	tempt.clear();
      }
      else{
	double s=w*w;
	double E_phi=(s-m_p_sq+m_phi_sq)/(2*w);
	double p_phi_sq=E_phi*E_phi-m_phi_sq;
	double p_phi=sqrt(p_phi_sq);
	double E_phi_sq=E_phi*E_phi;
	double p_beam=(s-m_p_sq)/2./w;
	Eclas=(s-m_p_sq)/(2.*m_p);      
	double myt=pow(E_phi-p_beam,2)-pow(p_phi-p_beam,2)
	  -2.*p_phi*p_beam*(1.-costheta);
	
	double scale=0.5/(p_beam*p_phi);
	tempt.push_back(-myt);
	tempx.push_back(scale*xclas);
	tempdx.push_back(scale*dxclas);
      }
      
      Eclas_old=Eclas;
      w_old=w;

    }
    
    double diff_min=1e6;
    for (unsigned int i=0;i<Eclas_vec.size();i++){
      if (fabs(Eg-Eclas_vec[i])<diff_min){
	diff_min=fabs(Eg-Eclas_vec[i]);
	ind=i;
      }
    }
    
    double y_acc[50],dy_acc[50]; //with kaon id
    double y2_acc[50],dy2_acc[50];
    for (int i=0;i<50;i++){
      //cout << eff[i] << endl;
      if (eff[i]>0.){
	double scale=1e6/1.293/Nps_sum/0.492/0.04;
	y_acc[i]=scale*y[i]/eff[i];
	// cout << y[i] << endl;
	dy_acc[i]=y_acc[i]*sqrt(deff[i]*deff[i]/(eff[i]*eff[i])+dy[i]*dy[i]/(y[i]*y[i]));  
	y2_acc[i]=scale*y2[i]/eff2[i];
	// cout << y[i] << endl;
	dy2_acc[i]=y2_acc[i]*sqrt(deff[i]*deff[i]/(eff[i]*eff[i])+dy2[i]*dy2[i]/(y2[i]*y2[i]));
      }
      t[i]=0.02+0.04*i;
      dt[i]=0.02;
    }
    
    TCanvas *c2=new TCanvas("c2","c2");
    c2->Draw();
    TH2F *h=new TH2F("h","Differential cross section for #phi",10,0,1.5,10,0,1.5);
    h->SetXTitle("-t [GeV^{2}]");
    h->SetYTitle("d#sigma/dt [#mub/GeV^{2}]");
    h->Draw();

    int n=2;
    TGraphErrors *g1=new TGraphErrors(35-n,&t[n],&y_acc[n],&dt[n],&dy_acc[n]);
    g1->SetLineColor(2);
    g1->SetMarkerColor(2);
    g1->SetMarkerStyle(20);
    //g1->Fit("expo");
    g1->Draw("p"); 
    
    TGraphErrors *g2=new TGraphErrors(35-n,&t[n],&y2_acc[n],&dt[n],&dy2_acc[n]);
    g2->SetLineColor(4);
    g2->SetMarkerColor(4);
    g2->SetMarkerStyle(21);
    TF1 *f1=new TF1("f1","[0]*exp([1]*x)",0.1,1.);
    f1->SetParameters(1,-4);
    //  g2->Fit(f1,"r");
    g2->Draw("p");
    
    
    double ballamx[5]={1.5,0.68,0.23,0.15,0.04};
    double ballamdx[5]={0.3,0.15,0.11,0.10,0.03};
    double ballamt[5]={0.11,0.3,0.5,0.7,1.2};
    double ballamdt[5]={0.09,0.1,0.1,0.1,0.4};
    
    double ballamx_lo[5]={1.22,0.44,0.17,0.24,0.03};
    double ballamdx_lo[5]={0.22,0.12,0.10,0.10,0.02};
    double ballamt_lo[5]={0.1225,0.3,0.5,0.7,1.2};
    double ballamdt_lo[5]={0.0775,0.1,0.1,0.1,0.4};
    
    TGraphErrors *g3=NULL;
    if (Eg>6){
      g3=new TGraphErrors(5,ballamt,ballamx,ballamdt,ballamdx);
    }
    else{
      g3=new TGraphErrors(5,ballamt_lo,ballamx_lo,ballamdt_lo,ballamdx_lo);
    }
    g3->SetLineColor(8);
    g3->SetMarkerColor(8);
    g3->SetMarkerStyle(22);
    g3->Draw("p");
    
    
    double tclas_arr[24],xclas_arr[24],dxclas_arr[24];
    int num_clas=tclas_vec[ind].size();
    for (int i=0;i<num_clas;i++){
      tclas_arr[i]=tclas_vec[ind][i];
      xclas_arr[i]=xclas_vec[ind][i];
      dxclas_arr[i]=dxclas_vec[ind][i];
    }
    
    TGraphErrors *g4=new TGraphErrors(num_clas,tclas_arr,xclas_arr,0,dxclas_arr);
    g4->SetLineColor(6);
    g4->SetMarkerColor(6);
    g4->SetMarkerStyle(25);
    g4->Draw("p");
    
    
    TLegend *legend =new TLegend(0.4,0.5,0.8,0.7);
    char mytitle[80];
    sprintf(mytitle,"Method 1, E#gamma=%3.1f GeV ",Eg);
    legend->AddEntry(g1,mytitle,"p");
    sprintf(mytitle,"Method 2, E#gamma=%3.1f GeV ",Eg);
    legend->AddEntry(g2,mytitle,"p");
    if (Eg>6){
      legend->AddEntry(g3,"Ballam 1973 E#gamma= 9.3 GeV","p");
    }
    else{
      legend->AddEntry(g3,"Ballam 1973 E#gamma= 4.7 GeV","p");
    }
    sprintf(mytitle,"CLAS 2013 E#gamma=%4.2f GeV", Eclas_vec[ind]);
    legend->AddEntry(g4,mytitle,"p");
    legend->Draw();
    
    sprintf(mytitle,"phi_t_xsecs_%3.1f.png",Eg);
    c2->SaveAs(mytitle);
  }
}