void test_macro()//arguements:th1, th2, y axis max 
{
  gROOT->Reset();
  gStyle->SetNdivisions(510,"X");
  gStyle->SetNdivisions(505,"Y");
  gStyle->SetStatBorderSize(1);
  gStyle->SetStatW(0.2);
  gStyle->SetStatH(0.1);
  //gStyle->SetFrameBorderMode(1);
  //gStyle-> SetDrawBorder(0);
  //gStyle->SetOptTitle(0);
  gStyle->SetCanvasBorderMode(0);
  gStyle->SetPadBorderMode(0);
  gStyle->SetCanvasColor(0);
  gStyle->SetFrameBorderMode(0);
  gStyle->SetPadColor(0);
  gStyle->SetStatColor(0);
  gStyle->SetLabelSize(.06,"xy");
  gStyle->SetLabelOffset(.011,"xy");
  gStyle->SetTitleXOffset(1.1);
  gStyle->SetTitleXSize(0.06);
  gStyle->SetTitleYOffset(1.12);
  gStyle->SetTitleYSize(0.06);
  gStyle->SetTitleW(0.5);
  gStyle->SetTitleH(0.08);
  gStyle->SetTitleBorderSize(1);
  gStyle->SetTitleFillColor(0);
  gStyle->SetPadLeftMargin(0.18);
  gStyle->SetPadRightMargin(0.15);
  gStyle->SetPadBottomMargin(.18);
  gStyle->SetTextSize(.03);
  gStyle->SetStatFontSize(.03);
  gStyle->SetOptStat(11);
  gStyle->SetOptTitle(0);
  gStyle->SetPalette(1);
  //___________________________________________________________________  
  Int_t TH1 = 16, TH2 =16;		  
  TGraph *tg1, *tg2, *tg3, *tg4;
  int nn1 = 0, nn2 = 0;
  float kB1[5500], kB2[5500], kB3[5500], kB4[5500], kB5[5500];

  char line[256], cc1[1], kHash1[1], buff1[256], ch1;
  char line[256], cc2[1], kHash2[1], buff2[256], ch2;
  float var1a, var2a, var3a, var4a, var5a;
  float var1b, var2b, var3b, var4b, var5b;
  float sigma1[200], ss1[200], ed1[200], ep1[200], en1[200];
  float sigma2[200], ss2[200], ed2[200], ep2[200], en2[200];
  float ee1, th1a, th2a, phi12a, dummy1a;
  float ee2, th1b, th2b, phi12b, dummy1b;
  int points1, dummy2a, dummy3a;
  int points2, dummy2b, dummy3b;
  //###########
  //read a text file
  //fill vectors
  //plot TGraphError
  //save it as png file
  char file1[500], file2[500];
  char text1[500];
  float kS1, kCounts1, kNorm;
  float kS[200], kS_err[200], kCounts[200];
  float kSQRT_count1;
  float kSQRT_count[200];
  FILE *pfile, *pfile1, *pfile2; 
  TGraphErrors *tg;
  int k1[49], k2[49];
  float k3[49], k4[49];
  int kk1, kk2;
  float kk3, kk4;
  int m1[49], m2[49];
  float m3[49], m4[49];
  int mm1, mm2;
  float mm3, mm4;
  float deltaS = 4.0;
  float tot_effi_d[7], tot_effi_p[7];

  //normalising factor:
  float n_el[7] = {1465730, 793332, 502017, 371621, 279511, 259465, 239448};
  float dOhm_p[7] = {0.0605, 0.0678, 0.0750, 0.0822, 0.0892, 0.0961, 0.1030};
  float sigma_el[7] = {7.5, 3.55, 2.0, 1.38, 1.02, 0.78, 0.62};
  float effi_dE[7]={0.9980, 0.9981, 0.9982, 0.9982, 0.9984, 0.9984, 0.9981};
  float effi_MWPC_D[7]={0.8698, 0.8759, 0.8739, 0.8712, 0.8768, 0.8719, 0.8761};//for deuterons
  float effi_MWPC_P[7]={0.8496, 0.8533, 0.8484, 0.8484, 0.8499, 0.8642, 0.8539};//for protons
  float effi_cross_over_el[7]={0.9929, 0.9938, 0.9898, 0.9915, 0.9937, 0.9936, 0.9902};//from elastic d simulation
  float effi_conf[7][7];
  float effi_cross_over_br[7][7];//for breakup configurations
  int TriggReScale = 32; // to scale back the single trigg events of el d
  int micro2mili = 1000;//convert the elastic cro-sec from micro to mili barn (because theory is in mili)
  float nonrel2rel = 3.0;//approximate conversion from non-reltivistic theory to relativistic (in MeV)
  //##############################
  // open config & cross-over effi 
  //##############################
/*  int n1 = 0;
  pfile1 = fopen("conf.effi","r");
  while(true)
    {
      fscanf(pfile1," %d %d %f %f", &kk1, &kk2, &kk3, &kk4);
      if (feof(pfile1)) break;
      n1++; 
      effi_conf[0.5*kk1-8][0.5*kk2-8] = kk3;
      //  cout<<n1<<"--- "<<kk1<<" "<<kk2<<" "<<kk3<<" "<<kk4<<"\n";
    }
      cout<<"HELLO = "<<effi_conf[0.5*TH1-8][0.5*TH2-8]<<"\n";
  //##############################
  int n2 = 0;
  pfile2 = fopen("cros_over.effi","r");
  while(true)
    {
      fscanf(pfile2," %d %d %f %f", &mm1, &mm2, &mm3, &mm4);
      if (feof(pfile2)) break;
      n2++; 
      effi_cross_over_br[0.5*mm1-8][0.5*mm2-8] = mm3;
    }
  cout<<"HELLO2 = "<<effi_cross_over_br[0.5*TH1-8][0.5*TH2-8]<<"\n";

  //##############################
  tot_effi_d[0.5*TH1-8] = effi_dE[0.5*TH1-8]*effi_MWPC_D[0.5*TH1-8];
  tot_effi_p[0.5*TH2-8] = effi_dE[0.5*TH2-8]*effi_MWPC_P[0.5*TH2-8];
*/
  //   kNorm = (n_el[0.5*TH1-8]*dOhm_p[0.5*TH2-8])/(0.25*sigma_el[0.5*TH1-8]);//old
//  kNorm = (TriggReScale*n_el[0.5*TH1-8]*dOhm_p[0.5*TH2-8]*deltaS*tot_effi_d[0.5*TH1-8]*tot_effi_p[0.5*TH2-8]*effi_cross_over_br[0.5*TH1-8][0.5*TH2-8]*effi_conf[0.5*TH1-8][0.5*TH2-8])/(micro2mili*sigma_el[0.5*TH1-8]);//new
   kNorm = 1.0;
  //cout<<kNorm<<"\n";
  //  sprintf(file1,"./test/S_thd%d_thp%d_fi180.dat", TH1, TH2);// NO BACKGROUND SUB.	
  sprintf(file1,"S_bg_sub_thd%d_thp%d_fi180.dat", TH1, TH2);//WITH BACKGROUND SUB.	
  int n = 0;
  pfile = fopen(file1,"r");
  while(true)
    {
      fscanf(pfile," %f %f %f", &kS1, &kCounts1, &kSQRT_count1);
      if (feof(pfile)) break;
      kS[n] = kS1;
      kS_err[n] = 0;
      kCounts[n] = kCounts1/kNorm;
      kSQRT_count[n] = kSQRT_count1/kNorm;
      n++; 
    }
  tg = new TGraphErrors(n, kS, kCounts, kS_err, kSQRT_count);
  tg->SetMarkerColor(kOrange);
  tg->SetMarkerStyle(21);
  tg->SetMarkerSize(1.0);
  tg->GetXaxis()->SetTitle("S [MeV]");
  tg->GetYaxis()->SetTitle("#frac{d^{5}#sigma}{dSd#Omega_{1}d#Omega_{2}} [#mub#upointsr^{-2}#upointMeV^{-1}]");
  tg->GetXaxis()->CenterTitle();
  tg->GetYaxis()->CenterTitle();
  //  tg->SetTitle("mean value");
  float gunn = 20.0;

  //#############################
  // ### OPEN theory ###
  //#############################
/*
  FILE *infile1 = fopen("./d160_180.dat", "r"); // test lines
  while (  fgets(buff1, 255, infile1) != NULL) 
    {
      ch1 = buff1[0];
      // if (ch1 = ' ')
      // 	continue;
      if ( ch1 == '#') // some comment characters 
	{
	  //#  points  E  th1  0.0  th2  phi12  1   3
	  // sscanf(buff, "%c %f %d %f %f %f %f %d %d", &control, &ener, &trush, &th1, &phi1, &th2, &phi2, &trush, &trush);
	  sscanf(buff1, "%c %f %d %f %f %f %f %d %d", &kHash1, &ee1, &points1, &th1a, &dummy1a, &th2a, &phi12a, &dummy2a, &dummy3a);
	  if (fabs(TH1-th1a)<0.1 && fabs(TH2-th2a)<0.01)
	    fprintf(stderr, "%c %f %d %f %f %f %f %d %d\n",  kHash1[0], ee1, points1, th1a, dummy1a, th2a, phi12a, dummy2a, dummy3a);
	}
      else
	if (fabs(TH1-th1a)<0.1 && fabs(TH2-th2a)<0.01)
	  {
	    //fprintf(stderr, "%f %f %f %f\n", TH1, th1, TH2, th2);
	    sscanf(buff1, "%f %f %f %f %f", &var1a, &var2a, &var3a, &var4a, &var5a);
	    
	    // fprintf(stderr, "%f %f %f %f %f\n", var1, var2, var3, var4, var5);
	    ss1[nn1] = var1a-nonrel2rel;// S value
	    sigma1[nn1] = var2a;// d5sigma value
	    ed1[nn1] = var3a;// Ed
	    ep1[nn1] = var4a;// Ep
	    // en1[nn1] = 0.1*var5a;// En
	    en1[nn1] = var5a/gunn;// En
	    nn1++; 
	  }
    }
      
  tg1 = new TGraph(nn1, ss1, sigma1);
  tg3 = new TGraph(nn1, ss1, en1);
  tg1->SetLineColor(4);
  tg3->SetLineColor(2);
  tg3->SetLineStyle(2);
*/
    //#############################
/*
  FILE *infile2 = fopen("./d160_180_4.dat", "r"); // SECOND version theory of Deltuva
  while (  fgets(buff2, 255, infile2) != NULL) 
    {
      ch2 = buff2[0];
      // if (ch2 = ' ')
      // 	continue;
      if ( ch2 == '#') // some comment characters 
	{
	  //#  points  E  th1  0.0  th2  phi12  1   3
	  // sscanf(buff, "%c %f %d %f %f %f %f %d %d", &control, &ener, &trush, &th1, &phi1, &th2, &phi2, &trush, &trush);
	  sscanf(buff2, "%c %f %d %f %f %f %f %d %d", &kHash2, &ee2, &points2, &th1b, &dummy1b, &th2b, &phi12b, &dummy2b, &dummy3b);
	  if (fabs(TH1-th1b)<0.1 && fabs(TH2-th2b)<0.01)
	    fprintf(stderr, "%c %f %d %f %f %f %f %d %d\n",  kHash2[0], ee2, points2, th1b, dummy1b, th2b, phi12b, dummy2b, dummy3b);
	}
      else
	if (fabs(TH1-th1b)<0.1 && fabs(TH2-th2b)<0.01)
	  {
	    //fprintf(stderr, "%f %f %f %f\n", TH1, th1, TH2, th2);
	    sscanf(buff2, "%f %f %f %f %f", &var1b, &var2b, &var3b, &var4b, &var5b);
	    
	    // fprintf(stderr, "%f %f %f %f %f\n", var1, var2, var3, var4, var5);
	    ss2[nn2] = var1b-nonrel2rel;// S value
	    sigma2[nn2] = var2b;// d5sigma value
	    ed2[nn2] = var3b;// Ed
	    ep2[nn2] = var4b;// Ep
	    en2[nn2] = var5b;// En
	    nn2++; 
	  }
    }
      
  tg2 = new TGraph(nn2, ss2, sigma2);
  tg2->SetLineColor(3);
*/
   
  //#############################
  //##### DRAW ####
  //#############################
  TCanvas *can1 = new TCanvas ("can1","can1",700,450);
    Float_t ymax_my = 1.0;

  //  can1->cd();
  tg->Draw("AP");
  TAxis *axis = tg->GetXaxis();
  axis->SetLimits(20.,190.);                 // along X
  tg->GetHistogram()->SetMaximum(ymax_my);   // along          
  // tg->GetHistogram()->SetMaximum(2.);   // along          
  tg->GetHistogram()->SetMinimum(0);  //   Y

  sprintf(text1,"#it{#theta_{p}}=%d#circ", TH2);	
  
  pt1=new TPaveText(0.2,0.8,0.4,0.85,"NDC");
  pt1->SetFillColor(0);
  pt1->AddText(text1);
  pt1->SetBorderSize(0);
  pt1->SetTextFont(82);
  pt1->SetTextSize(0.05);
  pt1->SetTextColor(1);
  pt1->Draw(); 
 

 // tg1->Draw("same");
 // tg2->Draw("same");
 // tg3->Draw("same");


  //  sprintf(file2,"./test/S_thd%d_thp%d_fi180.png", TH1, TH2);	
    can1->SaveAs(Form("S_thd%d_thp%d_fi180.png", TH1, TH2));  

  Float_t rightmax = ymax_my*gunn;
  // Float_t rightmax = 2.0*gunn;
  //  Float_t rightmax = 20.0;
  TGaxis *right_axis = new TGaxis(gPad->GetUxmax(),gPad->GetUymin(),
			    gPad->GetUxmax(), gPad->GetUymax(),0,rightmax,505,"+L");
  right_axis->SetLineColor(kRed);
  right_axis->SetTextColor(kRed);
  right_axis->SetTitle("neutron energy [MeV]");
  right_axis->SetTitleSize(0.06);
  right_axis->CenterTitle();
  right_axis->SetLabelOffset(0.008);
  right_axis->SetTextFont(82);
  right_axis->SetLabelFont(82);
  right_axis->SetLabelSize(0.06);
  right_axis->SetLabelColor(kRed);
  right_axis->SetTitleOffset(0.9);	
  // right_axis->SetName("right_axis");
  right_axis->Draw("same");
  can1->Modified(); 
  can1->Update();
    //can1->SaveAs(Form("S_thd%d_thp%d_fi180.png", TH1, TH2));  



}