Script for producing abnormal ratio plots for given phi bin.

// Abnormal Ratio Program

// This is a C++ program to be used with ROOT that subracts pedestals (run1060) from two quartz fiber
// runs, one with one fiber illuminated and the others with all 8 fibers illuminated and finds their 
// ratio. The average ratio for a particular tower is calculated and the ratio of one fiber to all fibers
// illuminated is divided by  the average ratio. The graphs show ratio/(averaged ratio) on the vertical
// axis and phi bin on the horizontal axis.

#include 

//====================================
//====================================
//====================================
float meancalc(int cent,float yval[],int offset){
  float meanval=0;
  float N=0;
  int low=50;
  
  if (cent2)
	    {
	      rr[j]=0;
	     printf("ratior > 2: j= %d\n",j);
	   }

	  cnt++;
	}

      if (cnt==13)
	{
	  cnt=1;
	  i++;
	}
    }
}
//====================================
//====================================
//====================================
abnormPhi() {

  gStyle->SetOptStat(0);
  TCanvas *c1=new TCanvas("layer","Laser Data: Single quartz-fiber tests",400,400);

  TH2F *hd[13];
 
  hd[1] = new TH2F("","",1,0,13,1,0,2);
   
  c1->cd(1);
  hd[1]->Draw();

  // Parameters to change: run numbers 
  int runa=1071;
  int runb=1074;
  int peds=1060;

  int absector=0;
  int absubsector=0;

  printf("Enter the sector (1-12): \n");
  cin >> absector;
  cout << "sector is " << absector << endl;
  cout << "Enter the subsector (1-5, where 'A'=1,'B'=2,...): " << endl;
  cin >> absubsector;
  cout << "subsector is " << absubsector << endl;

  
  int begin = ((absector-1)*5+(absubsector-1))*12+1;
  printf("begin =%d\n",begin);
  int end = ((absector-1)*5+(absubsector-1))*12+12;  
  printf("end =%d\n",end);
  char sub = 'Z';

  if (absubsector==1)
    sub ='A';
  else if (absubsector==2)
    sub ='B';
  else if (absubsector==3)
    sub ='C';
  else if (absubsector==4)
    sub ='D';
  else if (absubsector==5)
    sub ='E';

  int fibnum=0;
  if (runb==1072)
    fibnum=1;
  else if (runb==1073)
    fibnum=2;  
  else if (runb==1074)
    fibnum=3;
  else if (runb==1075)
    fibnum=4;
  else if (runb==1076)
    fibnum=5;
  else if (runb==1078)
    fibnum=6;
  else if (runb==1079)
    fibnum=7;
  else if (runb==1074)
    fibnum=8;

  char titleHist[100];
  sprintf(titleHist," ratio/avg ratio %02dT%c-Fiber %d",absector,sub,fibnum);
  printf("title %s\n",titleHist);

  hd[1]->SetTitle(titleHist);
  hd[1]->SetXTitle("eta");
  hd[1]->SetYTitle("ratio/avg ratio");

  c1->SetBorderMode(0);

  TString hFilea="/N/u1/bullinger/outLaser/pi";
  hFilea+=runa; hFilea+=".hist.root";
 
  TString hFileb="/N/u1/bullinger/outLaser/pi";
  hFileb+=runb; hFileb+=".hist.root";
  
  TString hFilec="/N/u1/bullinger/outLaser/pi";
  hFilec+=peds; hFilec+=".hist.root";

  fda=new TFile(hFilea);
  assert(fda->IsOpen());

  fdb=new TFile(hFileb);
  assert(fdb->IsOpen());
  
  fdc=new TFile(hFilec);
  assert(fdc->IsOpen());
  
  int j=1;
  int i=1;
  int k,section;
  int index=1;
  char ch;
  int count=1;
  float centa=0;
  float centb=0;
  float centc=0;
  float meana[800];
  float meanb[800];
  float meanc[800];
  float neta[65];
  float ratio[800];
  float ratior[800];
  float average[65];
  float rb[65];          // array to be used as y-axis in graph
  int bb=0;

  // Create an array filled with eta values-to be used in the graph on the x-axis
    for (int d=0; d<12;d++){
     neta[d]=d+1;
   }
 
  //Loop through all sec values 1 to 12
  for (section=1;section<=12;section++)
    // Loop through all sub values 'A' to 'E'
    for (ch='A';ch<='E';ch++)
	// Loop through all eta values 1 to 12
	for(index=1;index<=12;index++)
	  { 
 
	    int sec=section,eta=index;
	    char sub=ch;
	    char namea[100];
            char nameb[100];
	    
	    sprintf(namea,"%02dT%c%02d",sec,sub,index);
	    sprintf(nameb,"%02dT%c%02d",sec,sub,index);
            
	    TH1F *ha=fda->Get(namea);
	    assert(ha);
  
            TH1F *hb=fdb->Get(nameb);
            assert(hb);
	    
	    TH1F *hc=fdc->Get(namea);
	    assert(hc);

	    centa=ha->GetMaximumBin();
            centb=hb->GetMaximumBin();
	    centc=hc->GetMaximumBin();

	    float *ya=ha->GetArray();
            float *yb=hb->GetArray();
	    float *yc=hc->GetArray();
	  

	    // Calculating means for run with all fibers illuminated and mean for one fiber illuminated	
	    meana[j]=meancalc(centa,ya,50);
            meanb[j]=meancalc(centb,yb,50);
	    
	    // Calculating mean for pedestal run
	    meanc[j]=meancalc(centc,yc,5);

	    // Subtracting pedestals
	    meana[j]=meana[j]-meanc[j]; 
	    meanb[j]=meanb[j]-meanc[j];

	    // Flagging bad towers (with just pedestals). If mean a or b is around +/-10, set ratio to 0
	    if (((meana[j]>-10)&&(meana[j]<10))||((meanb[j]>-10)&&(meanb[j]<10)))
	      {
	        ratio[j]=0;
	      }
	    else
	      ratio[j]=meanb[j]/meana[j];
	    
	    // If ratio is negative or greater than 2 or equal to 1, set ratio to zero
	    if (ratio[j]<0 || ratio[j]>2 || ratio[j]==1)
	      ratio[j]=0;

	    j++;
	  }
 
  // Calculate the average ratio for one phi bin (the average ratio for all 12 eta)
  avgcalc(average,ratio);


  // Calculate the "ratio of ratios": ratior=ratio/[avg ratio * average]
  ratiorcalc(average,ratio,ratior);
  
  // Make array to fill graph for each eta
  for (k=begin;k<=end;k++)
    {
      rb[bb]=ratior[k];
      printf("ratio %f\n",rb[bb]);
      bb++;
    }

 
  // Make graph
  int mc=2;  // marker color
  TGraph *g1 = new TGraph(12,neta,rb); 
  g1->SetMarkerStyle(7);
  g1->SetMarkerColor(mc);
  c1->cd(1);
  g1->Draw("P");

  c1->cd();

   return 0;
}