// The plot of rho vs time from Greiner and Muller is suspicious.
// Look at raw numbers from other sources.
// 
// Andrei Gaponenko, 2004, 2005

enum{MAXSET=30};

rho_measurements() {
  double labelSize = 0.06;
  double axisTitleSize = 0.08;
  
  gStyle  ->SetLabelFont(132);
  gStyle  ->SetLabelFont(132,"X");
  gStyle  ->SetLabelFont(132,"Y");
  gStyle  ->SetLabelFont(132,"Z");
  gStyle  ->SetLabelFont(132,"T");
  gStyle  ->SetTitleFont(132);
  gStyle  ->SetTitleFont(132,"X");
  gStyle  ->SetTitleFont(132,"Y");
  gStyle  ->SetTitleFont(132,"Z");
  gStyle  ->SetTitleFont(132,"T");
  gStyle  ->SetTextFont(132);
  gStyle  ->SetStatFont(132);
  gStyle  ->SetPadColor(10);
  gStyle  ->SetCanvasColor(10);


  double rhov[MAXSET];
  double rhoe[MAXSET];
  double year[MAXSET];
  int i=-1;
  double eps = 0.18;

  // Table II from 
  // Lawrence Rosenson, "Momentum Spectrum of Positrons 
  // from the Decay of mu+ Mesons",
  // Phys Rev 109 (1958) 958.
  i++; rhov[i] = 0.075; rhoe[i] = 0.20;  year[i] = 1949;
  i++; rhov[i] = 0.19 ; rhoe[i] = 0.12 ; year[i] = 1951;
  i++; rhov[i] = 0.41 ; rhoe[i] = 0.13 ; year[i] = 1952;
  i++; rhov[i] = 0.26 ; rhoe[i] = 0.26 ; year[i] = 1952+eps;
  i++; rhov[i] = 0.50 ; rhoe[i] = 0.13 ; year[i] = 1954;
  i++; rhov[i] = 0.68 ; rhoe[i] = 0.11 ; year[i] = 1955;
  i++; rhov[i] = 0.57 ; rhoe[i] = 0.14 ; year[i] = 1956-eps;
  i++; rhov[i] = 0.62 ; rhoe[i] = 0.031; year[i] = 1956;
  i++; rhov[i] = 0.72 ; rhoe[i] = 0.058; year[i] = 1956+eps;
  // The last entry from this table is in the next source
  // (more consistent year?).
  //i++; rhov[i] =      ; rhoe[i] =      ; year[i] =     ;


// Table I from 
// Bruce Arne Sherwood, "Measurement of the Muon-Decay Spectrum 
// with a Wire Spark-Chamber Spectrometer",
// Phys Rev 156 (1967) 1475
  i++; rhov[i] = 0.67 ; rhoe[i] = 0.05 ; year[i] = 1958;
  i++; rhov[i] = 0.741; rhoe[i] = 0.027; year[i] = 1959;
  i++; rhov[i] = 0.780; rhoe[i] = 0.025; year[i] = 1960;
  //  i++; rhov[i] = 0.774; rhoe[i] = 0.042; year[i] = 1961;
  i++; rhov[i] = 0.751; rhoe[i] = 0.034; year[i] = 1962; // mu-
  i++; rhov[i] = 0.661; rhoe[i] = 0.016; year[i] = 1964;
  i++; rhov[i] = 0.867; rhoe[i] = 0.035; year[i] = 1964+eps; // mu-
  //i++; rhov[i] = 0.750; rhoe[i] = 0.003; year[i] = 1965;
  //i++; rhov[i] = 0.760; rhoe[i] = 0.009; year[i] = 1965+eps;



// PDG, including unused data
  i++; rhov[i] = 0.7503; rhoe[i] = 0.0026; year[i] = 1966;
  i++; rhov[i] = 0.760; rhoe[i] = 0.009; year[i] = 1967;
  i++; rhov[i] = 0.762; rhoe[i] = 0.008; year[i] = 1968;
  i++; rhov[i] = 0.7518; rhoe[i] = 0.0026; year[i] = 1969;
  i++; rhov[i] = 0.7508; rhoe[i] = 0.0011; year[i] = 2005;
  const int n = i+1;

  double yearerr[MAXSET];
  for(i=0; i<MAXSET; i++) { yearerr[i] = 0; }

  TGraphErrors *g = new TGraphErrors(n, year, rhov, yearerr, rhoe);
  g->SetTitle("rho vs year");
  TCanvas *c = new TCanvas("rho_vs_year");
  c->SetLeftMargin(0.18);
  c->SetBottomMargin(0.15);
  c->SetTopMargin(0.05);
  c->SetRightMargin(0.05);
  // c->Divide(1,2)
  // TPad* p1 = c->cd(1);

  g->SetMarkerColor(4);
  g->SetMarkerStyle(21);
  g->SetMarkerSize(1.4);
  g->SetLineWidth(2);
  g->SetMinimum(0.);
  
  TAxis* yaxis = g->GetYaxis();
  yaxis->SetRangeUser(0.,1.);
  yaxis->SetNdivisions(-504);
  yaxis->SetLabelSize(labelSize);
  yaxis->SetLabelOffset(0.008);
  yaxis->SetTitle("Value of #rho");
  yaxis->SetTitleOffset(1.05);
  yaxis->SetTitleSize(axisTitleSize);
  yaxis->CenterTitle();

  TAxis* xaxis = g->GetXaxis();
  xaxis->SetLabelSize(labelSize);
  xaxis->SetTitle("Year");
  xaxis->SetTitleSize(axisTitleSize);
  xaxis->SetTitleOffset(0.8);
  xaxis->CenterTitle();
  
  g->SetMaximum(1.);


  gStyle->SetOptTitle(0);
  
  g->Draw("AP");

  //----------------------------------------------------------------
  double vadate = 1958.;

  TF1 *va = new TF1("va", "[0]", vadate, 2000.);
  va->SetParameter(0, 0.75);
  va->SetLineColor(2);
  va->Draw("same");

//  TLine *vayear = new TLine(vadate, 0., vadate, 1.);
//  vayear->SetLineColor(6);
//  vayear->Draw("same");

  c->Print("rho_vs_year.eps");

  g->GetHistogram()->GetXaxis()->SetRangeUser(1960, 1972);
  g->GetHistogram()->GetYaxis()->SetRangeUser(0.6,0.95);
  
  c->Print("rho_vs_year_1960.eps");

}