// code for plotting LED intensity

// Marieke Navin - March 2004

{
  gROOT->Reset();                   // reset global variables

  //if (!TClassTable::GetDict("TGraphIntegrator")) {
  //gSystem->Load("tgint/source/libTGrint.so");
  //cout << "loading integrator\n";
  //}
  //else
  //  cout << "integrator already loaded\n";
  
  gStyle->SetCanvasColor(10);
  gStyle->SetTitleColor(10);
  gStyle->SetFrameFillColor(10);
  gStyle->SetOptStat(0);
  
  ifstream data;                        // declare a file to input from
  // data.open("led1new.data");   // open file containing values
  //data.open("led2.data");
  data.open("led3.data");
  //data.open("led4.data");
  const Int_t n = 100;                        // number of points to be plotted
  Int_t in = 0;

  Float_t x[n], y1[n], y2[n], y3[n], y4[n];
  Float_t x1;

  char line[1024];
  while (data >> x1)
    {
      x[in] = x1;
      data >> y1[in];
      data >> y2[in];
      data >> y3[in];
      data >> y4[in];
      in++;
    }


  for (Int_t i=0; i<in; i++)
    {
      cout << x[i] << " " << y1[i] << " " << y2[i] << " " << y3[i] << " " << y4[i]<< endl;
    }
 
  TGraph *inputSpec = new TGraph(in, x, y1);
  inputSpec->SetMinimum(10);   //inputSpec is the graph
  TGraph *inputSpec2 = new TGraph(in, x, y2);
  //inputSpec->SetTitle("Intensity of LED 1");
  // inputSpec->SetTitle("Intensity of LED 2");
  inputSpec->SetTitle("Intensity of LED 3");
  //inputSpec->SetTitle("Intensity of LED 4");

  TGraph *inputSpec3 = new TGraph(in, x, y3);
  TGraph *inputSpec4 = new TGraph(in, x, y4);

  inputSpec->GetHistogram()->SetXTitle("Cos Theta");
  inputSpec->GetHistogram()->SetTitleOffset(0.5);
  inputSpec->GetHistogram()->SetYTitle("pJ/pulse/sr");
  
  //TCanvas *c1 = new TCanvas("c1", "Intensity of LED 1",0,  0, 800, 600);  
  //TCanvas *c1 = new TCanvas("c1", "Intensity of LED 2",0,  0, 800, 600); 
  TCanvas *c1 = new TCanvas("c1", "Intensity of LED 3",0,  0, 800, 600); 
  //TCanvas *c1 = new TCanvas("c1", "Intensity of LED 4",0,  0, 800, 600); 


  c1->SetLogy(); 
  		    
  inputSpec->SetName("LED1positivetheta");
  inputSpec->Draw("AC");
  gPad->Update();

  inputSpec2->SetLineStyle(kDashDotted);
  inputSpec2->Draw("CP");
  inputSpec2->SetName("LED1negativetheta");
  inputSpec3->SetName("LED1cleavedpositivetheta");
  inputSpec4->SetName("LED1cleavednegativetheta");
  //inputSpec3->SetMarkerStyle(21); put box points on
  inputSpec3->SetLineStyle(kDotted);
  inputSpec3->Draw("CP");
  inputSpec4->SetLineStyle(kDashed);
  inputSpec4->Draw("CP");

  //draw an axis on the right hand side
  TGaxis *axis = new TGaxis(gPad->GetUxmax(),
			    gPad->GetUymin(),
			    gPad->GetUxmax(),
			    gPad->GetUymax(),
			    2.5e7,6e10,510,"+LG"); //G for log scale
  axis->SetTitle("Number of photons/pulse/sr");
  axis->Draw();

  //creating a legend
  TLegend *legend = new TLegend(0.4,0.2,0.65,0.44);
  legend->AddEntry(inputSpec,"Uncleaved Positive Theta", "l");
  legend->AddEntry(inputSpec2, "Uncleave Negative Theta", "l");
  legend->AddEntry(inputSpec3, "Cleaved Positive Theta", "l");
  legend->AddEntry(inputSpec4, "Cleaved Negative Theta", "l");
  legend->Draw();

//Create a new TGraphIntegrator and assign the graph to it
///TGraphIntegrator *theIntegrator = new TGraphIntegrator();
//theIntegrator->SetGraph(inputSpec);

//Integrate graph between x=0 and x=1. Note that the last 0 in the call is
//the function identifier (see the source code)
///Double_t integral = theIntegrator->Integrate(0,1,0);

/theIntegrator->SetGraph(inputSpec2);
//Double_t integral2 = theIntegrator->Integrate(0,1,0);
//theIntegrator->SetGraph(inputSpec3);
//Double_t integral3 = theIntegrator->Integrate(0,1,0);
//theIntegrator->SetGraph(inputSpec4);
//Double_t integral4 = theIntegrator->Integrate(0,1,0);

//Output results
//cout<<"Integral of Graph 1   = "<<integral<<endl;
//cout<<"Integral of Graph 2   = "<<integral2<<endl;
//cout<<"Integral of Graph 3   = "<<integral3<<endl;
//cout<<"Integral of Graph 4   = "<<integral4<<endl;

// cout<<"Mean uncleaved = "<<(integral + integral2) / 2.0 << " pJ/pulse/sr"
// << endl;
// cout <<"Mean cleaved = " <<(integral3 + integral4) /2.0 << " pJ/pulse/sr"
//      <<endl;

// cout << "Reduction = " <<((((integral + integral2) / 2.0) - ((integral3 + integral4) /2.0))/ ((integral + //integral2) / 2.0)) * 100 << " %" <<endl;

//Cleanup new'd pointers
//delete theIntegrator;

}