#include #include #include "TFile.h" #include "TH1D.h" #include "TProfile.h" #include "TCanvas.h" #include "TString.h" #include "TStyle.h" #include "TChain.h" #include "TH2.h" #include "TH1.h" #include "TF1.h" #include "TTree.h" #include "TKey.h" #include "Riostream.h" #include "TMath.h" #include "TRandom.h" #include "TRandom3.h" #include "TVirtualFitter.h" #include #include // Simple routine needed to calculate the inverse error function // ------------------------------------------------------------- Double_t ErfInverse (Double_t x); // Macro that computes p-value and Z-value of N observed vs B predicted // Poisson counts // -------------------------------------------------------------------- void Poisson_prob_fix (double B, double N) { int maxN = N*3/2; if (N<20) maxN=2*N; TH1D * Pois = new TH1D ("Pois", "", maxN, -0.5, maxN-1); TH1D * PoisGt = new TH1D ("PoisGt", "", maxN, -0.5, maxN-1); double sum=0.; double fact=1.; for (int i=0; i1) fact*=i; double poisson = exp(-B)*pow(B,i)/fact; if (iSetBinContent(i+1,poisson); if (i>=N) PoisGt->SetBinContent(i+1,poisson); } double P=1-sum; double Z = sqrt(2) * ErfInverse(1-2*P); cout << "P of observing N=" << N << " or more events if B=" << B << " : P= " << 1-sum << endl; cout << "This corresponds to " << Z << " sigma for a Gaussian one-tailed test." << endl; Pois->SetLineWidth(3); PoisGt->SetFillColor(kBlue); TCanvas* T = new TCanvas ("T","Poisson distribution", 500, 500); T->Divide(1,2); T->cd(1); Pois->Draw(); PoisGt->Draw("SAME"); T->cd(2); T->GetPad(2)->SetLogy(); Pois->Draw(); PoisGt->Draw("SAME"); } Double_t ErfInverse(Double_t y) { // returns the inverse error function obtained // y must be <-1