| 1 | void gammarate()
|
|---|
| 2 | {
|
|---|
| 3 |
|
|---|
| 4 | gROOT->Reset();
|
|---|
| 5 | c1 = new TCanvas("c1","Gamma rate",700,300,500,500);
|
|---|
| 6 | c1->SetGridx();
|
|---|
| 7 | c1->SetGridy();
|
|---|
| 8 | c1->SetLogy(1);
|
|---|
| 9 | c1->SetLogx(1);
|
|---|
| 10 |
|
|---|
| 11 | TFile f("area.root");
|
|---|
| 12 | MHMcCollectionArea* mhmc = (MHMcCollectionArea*)f.Get("MHMcCollectionArea");
|
|---|
| 13 | TH1D* h1 = (TH1D*)mhmc->GetHist();
|
|---|
| 14 |
|
|---|
| 15 |
|
|---|
| 16 | h1->GetXaxis();
|
|---|
| 17 | Float_t xmin = h1->GetXaxis()->GetXmin();
|
|---|
| 18 | Float_t xmax = h1->GetXaxis()->GetXmax();
|
|---|
| 19 | Int_t nbin = h1->GetNbinsX();
|
|---|
| 20 |
|
|---|
| 21 | // gamma spectrum:
|
|---|
| 22 | fun1 = new TF1("fun1","10000*2.706e-11*pow(x,-2.47-0.11*log(x))",xmin/1000.,xmax/1000.);
|
|---|
| 23 |
|
|---|
| 24 | TH1F *h3f = (TH1F*)h1->Clone("h3f");
|
|---|
| 25 | h3f->SetTitle("Integrale per bin");
|
|---|
| 26 | h3f->Reset();
|
|---|
| 27 | h3f->SetMinimum(1.e-13);
|
|---|
| 28 | h3f->SetMaximum(1.e-7);
|
|---|
| 29 |
|
|---|
| 30 | TH1F *h4f = (TH1F*)h1->Clone("h4f");
|
|---|
| 31 |
|
|---|
| 32 | h4f->SetTitle("Rate");
|
|---|
| 33 | h4f->Reset();
|
|---|
| 34 | h4f->SetMinimum(1.e-4);
|
|---|
| 35 | h4f->SetMaximum(1.e-1);
|
|---|
| 36 |
|
|---|
| 37 |
|
|---|
| 38 | Float_t Error = 0.;
|
|---|
| 39 |
|
|---|
| 40 | cout << "Bin\tArea\tCenter\tFlux\t\tMolt\tInt\tRate"<<endl;
|
|---|
| 41 | for (Int_t i = 1; i <= nbin; i++)
|
|---|
| 42 | {
|
|---|
| 43 | //area->GetXaxis()->FindBin(i)
|
|---|
| 44 | Float_t valArea = h1->GetBinContent(i);
|
|---|
| 45 | Float_t bincenter = h1->GetBinCenter(i);
|
|---|
| 46 | Float_t valFlusso = fun1->Eval(bincenter/1000., 0., 0.);
|
|---|
| 47 | Float_t Integral = fun1->Integral(h1->GetBinLowEdge(i)/1000.,h1->GetBinLowEdge(i+1)/1000.);
|
|---|
| 48 | h3f->SetBinContent(i,Integral);
|
|---|
| 49 | Float_t Rate = Integral * valArea;
|
|---|
| 50 | Float_t ErArea = h1->GetBinError(i);
|
|---|
| 51 | Float_t ErRate = ErArea * Integral;
|
|---|
| 52 | Error += ErRate*ErRate;
|
|---|
| 53 | h4f->SetBinContent(i,Rate);
|
|---|
| 54 | h4f->SetBinError(i,ErRate);
|
|---|
| 55 | cout << i << "\t" << bincenter << "\t" << Rate << "\t"<< endl;
|
|---|
| 56 | }
|
|---|
| 57 | //h1->DrawCopy();
|
|---|
| 58 | //h3f->DrawCopy();
|
|---|
| 59 | h4f->SetStats(kFALSE);
|
|---|
| 60 | h4f->SetTitle("Gamma rate vs energy");
|
|---|
| 61 | h4f->GetYaxis()->SetTitle("Rate");
|
|---|
| 62 | h4f->DrawCopy("pe");
|
|---|
| 63 |
|
|---|
| 64 | cout << "Rate (Hz):" << h4f->Integral(0,100)<< " +- " <<pow(Error, .5)<<endl;
|
|---|
| 65 |
|
|---|
| 66 | h4f->Delete();
|
|---|
| 67 | h3f->Delete();
|
|---|
| 68 | h1->Delete();
|
|---|
| 69 | }
|
|---|
