source: trunk/MagicSoft/Mars/mtemp/mpisa/macros/AlphaPlot.C@ 4671

void AlphaPlot(TString onsource="hillas_Mrk421.root", TString offsource="hillas_OffMrk421.root")
{
  const Int_t kMinDeg = 25;

  Int_t i;
  Float_t EvtOn = 0;
  Float_t EvtOff = 0;

  TH1F *hon = new TH1F("hon","On Source Histogram",19,0,95);
  TH1F *hoff = new TH1F("hoff","Off Source Histogram",19,0,95);
  TH1F *htemp;

  fon = new TFile(onsource,"READ");

  htemp = (TH1F *)fon->FindObjectAny("hAlpha");
  
  for(i=0;i<19;i++)
    hon->SetBinContent(i,htemp->GetBinContent(i));

  for(i=0;i<hon->GetNbinsX();i++){
    if(hon->GetBinCenter(i) > kMinDeg)
      EvtOn += hon->GetBinContent(i);
  }

  fon->Close();

  foff = new TFile(offsource,"READ");

  htemp = (TH1F *)foff->FindObjectAny("hAlpha");
  for(i=0;i<19;i++)
    hoff->SetBinContent(i,htemp->GetBinContent(i));

  for(i=0;i<hoff->GetNbinsX();i++){
    if(hoff->GetBinCenter(i) > kMinDeg)
      EvtOff += hoff->GetBinContent(i);
  }

  foff->Close();

  Float_t norm = EvtOn/EvtOff;
  cout << "EvtOn/EvtOff = " << norm << endl;

  Float_t eccessi = 0;
  Float_t Noff = 0;

  TH1F *hAlphaPlot = new TH1F("hAlphaPlot","Alpha Plot", 19,0,95);
  for(i=0; i<hon->GetNbinsX(); i++){
    hAlphaPlot->SetBinContent(i,hon->GetBinContent(i)-hoff->GetBinContent(i)*norm);
    if(hAlphaPlot->GetBinCenter(i) < 15){
      eccessi += hAlphaPlot->GetBinContent(i);
      Noff += norm*hoff->GetBinContent(i);
    }
  }

  cout << "Eccessi = " << eccessi << "   Nsigma = " << eccessi/(TMath::Sqrt(Noff)) <<endl;

  new TCanvas("Alpha On","Alpha plot ON source");
  hon->DrawCopy();
  new TCanvas("Alpha Off","Alpha plot OFF source");
  hoff->DrawCopy();
  new TCanvas("Alpha","Resulting alpha plot");
  hAlphaPlot->DrawCopy();
}