From 6403cc2131e0bd553874b340893359d40e7162ab Mon Sep 17 00:00:00 2001
From: Carlos Yero <cyero002@fiu.edu>
Date: Thu, 29 Jun 2017 10:41:51 -0400
Subject: [PATCH] Add calibration code for the SHMS aerogel detector.
---
CALIBRATION/shms_aero_calib/paero_calib.C | 107 ++++++++++++++++++++++
1 file changed, 107 insertions(+)
create mode 100644 CALIBRATION/shms_aero_calib/paero_calib.C
diff --git a/CALIBRATION/shms_aero_calib/paero_calib.C b/CALIBRATION/shms_aero_calib/paero_calib.C
new file mode 100644
index 00000000..6a7e864e
--- /dev/null
+++ b/CALIBRATION/shms_aero_calib/paero_calib.C
@@ -0,0 +1,107 @@
+#include <TFile.h>
+#include <TTree.h>
+#include <TStyle.h>
+#include <TCanvas.h>
+#include <TH1.h>
+#include <TF1.h>
+#include <iostream>
+#include <fstream>
+#include <iomanip>
+
+//
+// Calibrate SHMS aerogel detector by localizing SPE peaks.
+// The SPE peak positions are found from Gaussian fits to the peaks.
+// The fit ranges can be adjusted by changing fit range limits
+// flo_* and fhi_* in the body of code.
+
+// The code works on root output from hcana, the name of which is the
+// single input parameter of the program.
+
+#define NPMT 7
+#define MaxAdc 2000.
+#define NBin 1000
+
+void paero_calib(string fname) {
+
+ TFile *f = new TFile(fname.c_str());
+ TTree* tree;
+ f->GetObject("T",tree);
+
+ int nentries = tree->GetEntries();
+ cout << "nentries= " << nentries << endl;
+
+ double adc_pos[NPMT];
+ double adc_neg[NPMT];
+
+ TBranch* b_adc_pos;
+ TBranch* b_adc_neg;
+
+ tree->SetBranchAddress("P.aero.goodPosAdcPulseInt",adc_pos,&b_adc_pos);
+ tree->SetBranchAddress("P.aero.goodNegAdcPulseInt",adc_neg,&b_adc_neg);
+
+ TH1D* hpos[NPMT];
+ TH1D* hneg[NPMT];
+
+ for (int i=0; i<NPMT; i++) {
+ hpos[i] = new TH1D(Form("hpos%d",i+1), Form("ADC+ %i",i+1), NBin, 0.001, MaxAdc);
+ hneg[i] = new TH1D(Form("hneg%d",i+1), Form("ADC- %i",i+1), NBin, 0.001, MaxAdc);
+ }
+
+ for (int ientry=0; ientry<nentries; ientry++) {
+
+ tree->GetEntry(ientry);
+
+ for (int i=0; i<NPMT; i++) {
+ hpos[i]->Fill(adc_pos[i], 1.);
+ hneg[i]->Fill(adc_neg[i], 1.);
+ }
+
+ }
+
+ TCanvas* c1 = new TCanvas("adc_spec", "fADC spectra", 600, 800);
+ c1->Divide(2,NPMT);
+
+ double flo_neg[NPMT] {400., 400., 400., 400., 350., 400., 400.};
+ double fhi_neg[NPMT] {700., 650., 800., 700., 750., 650., 900.};
+
+ double flo_pos[NPMT] {250., 300., 400., 550., 250., 300., 250.};
+ double fhi_pos[NPMT] {500., 600., 700., 800., 500., 600., 550.};
+
+ float gain_pos[NPMT] {NPMT*0.};
+ float gain_neg[NPMT] {NPMT*0.};
+
+ int ip=0;
+ for (int i=0; i<NPMT; i++) {
+ c1->cd(++ip);
+ if (hpos[i]->GetSumOfWeights() > 0) {
+ hpos[i]->Fit("gaus","","",flo_pos[i],fhi_pos[i]);
+ hpos[i]->GetFunction("gaus")->SetLineColor(2);
+ hpos[i]->GetFunction("gaus")->SetLineWidth(2);
+ gain_pos[i] = hpos[i]->GetFunction("gaus")->GetParameter(1);
+ }
+ else
+ hpos[i]->Draw();
+ c1->cd(++ip);
+ if (hneg[i]->GetSumOfWeights() > 0) {
+ hneg[i]->Fit("gaus","","",flo_neg[i],fhi_neg[i]);
+ hneg[i]->GetFunction("gaus")->SetLineColor(2);
+ hneg[i]->GetFunction("gaus")->SetLineWidth(2);
+ gain_neg[i] = hneg[i]->GetFunction("gaus")->GetParameter(1);
+ }
+ else
+ hneg[i]->Draw();
+ }
+
+ ofstream of;
+ of.open("gain.r",ios::out);
+ of << "paero_neg_gain = ";
+ for (int i=0; i<NPMT; i++)
+ of << "1./" << gain_neg[i] << ", ";
+ of << endl;
+ of << "paero_pos_gain = ";
+ for (int i=0; i<NPMT; i++)
+ of << "1./" << gain_pos[i] << ", ";
+ of << endl;
+ of.close();
+
+}
--
GitLab