diff --git a/CALIBRATION/hms_hodo_calib/timeWalkCalib.C b/CALIBRATION/hms_hodo_calib/timeWalkCalib.C
new file mode 100644
index 0000000000000000000000000000000000000000..27ccf0dcc8c5175028dd34782b8439a8a65d55be
--- /dev/null
+++ b/CALIBRATION/hms_hodo_calib/timeWalkCalib.C
@@ -0,0 +1,296 @@
+// Macro to perform time-walk fits and extract the calibration parameters
+// Author: Eric Pooser, pooser@jlab.org
+
+// Declare ROOT files
+TFile *histoFile, *outFile;
+
+// Declare constants
+static const UInt_t nPlanes = 4;
+static const UInt_t nSides = 2;
+static const UInt_t nBarsMax = 21;
+static const UInt_t nTwFitPars = 3;
+
+static const Double_t tdcThresh = 120.0; // 30 mV in units of FADC channels
+static const Double_t twFitRangeLow = 10.0;
+static const Double_t twFitRangeHigh = 600.0;
+static const Double_t c0twParInit = -15.0;
+static const Double_t c1twParInit = 5.0;
+static const Double_t c2twParInit = 0.25;
+
+static const Double_t fontSize = 0.05;
+static const Double_t yTitleOffset = 0.75;
+static const Double_t markerSize = 2.0;
+static const Double_t minScale = 0.75;
+static const Double_t maxScale = 0.75;
+
+static const UInt_t lineWidth = 4;
+static const UInt_t lineStyle = 7;
+
+static const UInt_t nbars[nPlanes] = {16, 10, 16, 10};
+static const TString planeNames[nPlanes] = {"1x", "1y", "2x", "2y"};
+static const TString sideNames[nSides] = {"pos", "neg"};
+static const TString twFitParNames[nTwFitPars] = {"c_{0}", "c_{1}", "c_{2}"};
+
+// Declare directories
+TDirectory *dataDir, *planeDir[nPlanes], *sideDir[nPlanes][nSides];
+TDirectory *twDir[nPlanes][nSides];
+// Declare fits
+TF1 *twFit[nPlanes][nSides][nBarsMax], *avgParFit[nPlanes][nSides][nTwFitPars];
+// Declare arrays
+Double_t paddleIndex[nPlanes][nSides][nBarsMax];
+Double_t twFitPar[nPlanes][nSides][nTwFitPars][nBarsMax], twFitParErr[nPlanes][nSides][nTwFitPars][nBarsMax];
+Double_t avgPar[nPlanes][nSides][nTwFitPars];
+Double_t minPar[nPlanes][nSides][nTwFitPars], maxPar[nPlanes][nSides][nTwFitPars];
+// Declare canvases
+TCanvas *twFitCan[nPlanes][nSides], *twFitParCan[nTwFitPars];
+// Declare histos
+// 2D histos
+TH2F *h2_adcTdcTimeDiffWalk[nPlanes][nSides][nBarsMax];
+// Declare graphs
+TGraph *twFitParGraph[nPlanes][nSides][nTwFitPars];
+// Declare multi-graphs
+TMultiGraph *twFitParMultiGraph[nPlanes][nTwFitPars];
+// Declare legends
+TLegend *twFitParLegend[nPlanes][nTwFitPars], *avgParLegend[nPlanes][nTwFitPars];
+// Declare text objects
+TPaveText *twFitParText[nPlanes][nSides][nBarsMax];
+
+// Functions to extract calibration parameters
+//=:=:=:=:=:=:
+//=: Level 1
+//=:=:=:=:=:=:
+
+// Add color to fit lines
+void addColorToFitLine(UInt_t style, UInt_t width, UInt_t color, TF1 *fit1D) {
+ fit1D->SetLineStyle(lineStyle);
+ fit1D->SetLineWidth(lineWidth);
+ fit1D->SetLineColor(color);
+ return;
+} // addColorToFitLine()
+
+// Add color to fit graphs
+void addColorToGraph(UInt_t style, UInt_t size, UInt_t color, TGraph *graph) {
+ graph->SetMarkerStyle(style);
+ graph->SetMarkerSize(size);
+ graph->SetMarkerColor(color);
+ return;
+} // addColorToFitGraph()
+
+// Add centered titles to multigraphs
+void modifyMultiGraph(TMultiGraph *mg, TString xtitle, TString ytitle) {
+ mg->GetYaxis()->SetTitle(ytitle);
+ mg->GetYaxis()->CenterTitle();
+ mg->GetXaxis()->SetTitle(xtitle);
+ mg->GetXaxis()->CenterTitle();
+ return;
+} // modifyMultiGraph()
+
+// Create canvases
+TCanvas *makeCan(UInt_t numColumns, UInt_t numRows, UInt_t winWidth, UInt_t winHeight, TCanvas *can, TString name, TString title) {
+ can = new TCanvas(name, title, winWidth, winHeight);
+ can->Divide(numColumns, numRows);
+ return can;
+} // makeCan()
+
+// Define the time-walk fit function
+Double_t twFitFunc(Double_t *a, Double_t *c) {
+ Double_t twFitVal = c[0] + c[1]/(TMath::Power((a[0]/tdcThresh), c[2]));
+ return twFitVal;
+} // twFitFunc()
+
+// Locate min or max value from input array
+Double_t calcMinOrMax(Double_t *array, UInt_t iplane, TString minOrmax) {
+ auto result = minmax_element(array, array+nbars[iplane]);
+ if (minOrmax == "min") return *result.first;
+ else if (minOrmax == "max") return *result.second;
+ else return 0.0;
+} // calcMinOrMax()
+
+//=:=:=:=:=:=:
+//=: Level 2
+//=:=:=:=:=:=:
+
+// Perform the timw-walk fits
+void doTwFits(UInt_t iplane, UInt_t iside, UInt_t ipaddle) {
+ // Draw fits on canvas
+ twFitCan[iplane][iside]->cd(ipaddle+1);
+ gPad->SetLogz();
+ h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]->Draw("COLZ");
+ gPad->Modified(); gPad->Update();
+ // Declare and initialize the fits
+ twFit[iplane][iside][ipaddle] = new TF1("twFit", twFitFunc, twFitRangeLow, twFitRangeHigh, nTwFitPars);
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++)
+ twFit[iplane][iside][ipaddle]->SetParName(ipar, twFitParNames[ipar]);
+ twFit[iplane][iside][ipaddle]->SetParameters(c0twParInit, c1twParInit, c2twParInit);
+ // Perform the fits and scream if it failed
+ Int_t twFitStatus = h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]->Fit("twFit", "REQ");
+ if (twFitStatus != 0)
+ cout << "ERROR: Time Walk Fit Failed!!! " << "Status = " << twFitStatus << " For Plane: " << planeNames[iplane] << " Side: " << sideNames[iside] << " Paddle: " << ipaddle+1 << endl;
+ // Create text box to display fir parameters
+ twFitParText[iplane][iside][ipaddle] = new TPaveText(0.4, 0.6, 0.895, 0.895, "NBNDC");
+ twFitParText[iplane][iside][ipaddle]->AddText(Form("Entries = %.0f", h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]->GetEntries()));
+ // Obtain the fit parameters and associated errors
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ twFitPar[iplane][iside][ipar][ipaddle] = twFit[iplane][iside][ipaddle]->GetParameter(ipar);
+ twFitParErr[iplane][iside][ipar][ipaddle] = twFit[iplane][iside][ipaddle]->GetParError(ipar);
+ twFitParText[iplane][iside][ipaddle]->AddText(Form(twFitParNames[ipar]+" = %.2f #pm %.2f", twFitPar[iplane][iside][ipar][ipaddle], twFitParErr[iplane][iside][ipar][ipaddle]));
+ } // Parameter loop
+ // Draw the fit parameter text
+ twFitParText[iplane][iside][ipaddle]->SetFillColor(kWhite);
+ twFitParText[iplane][iside][ipaddle]->SetTextAlign(12);
+ twFitParText[iplane][iside][ipaddle]->SetFillStyle(3003);
+ twFitParText[iplane][iside][ipaddle]->Draw();
+ gPad->Modified(); gPad->Update();
+ return;
+} // doTwFits()
+
+// Calculate the averege of the time-walk fit parameters
+void calcParAvg(UInt_t iplane, UInt_t iside) {
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ // Calculate the weighted average while ignoring fit errors provided by Minuit
+ avgParFit[iplane][iside][ipar] = new TF1("avgParFit", "pol0", 1, nbars[iplane]);
+ avgParFit[iplane][iside][ipar]->SetParName(0, "#bar{"+twFitParNames[ipar]+"}");
+ // Add color to fit lines
+ if (iside == 0) addColorToFitLine(lineStyle, lineWidth, kRed, avgParFit[iplane][iside][ipar]);
+ if (iside == 1) addColorToFitLine(lineStyle, lineWidth, kBlue, avgParFit[iplane][iside][ipar]);
+ // Initialize the parameters
+ if (ipar == 0) avgParFit[iplane][iside][ipar]->SetParameter(0, c0twParInit);
+ if (ipar == 1) avgParFit[iplane][iside][ipar]->SetParameter(0, c1twParInit);
+ if (ipar == 2) avgParFit[iplane][iside][ipar]->SetParameter(0, c2twParInit);
+ // Perform the least squares fit
+ Int_t avgParFitStatus = twFitParGraph[iplane][iside][ipar]->Fit("avgParFit", "REQ");
+ if (avgParFitStatus != 0)
+ cout << "ERROR: Parameter Fit Failed!!! " << "Status = " << avgParFitStatus << " For Plane: " << planeNames[iplane] << " Side: " << sideNames[iside] << endl;
+ // Store the fit parameters
+ avgPar[iplane][iside][ipar] = avgParFit[iplane][iside][ipar]->GetParameter(0);
+ // Add graphs to multi graph
+ twFitParMultiGraph[iplane][ipar]->Add(twFitParGraph[iplane][iside][ipar]);
+ // Calculate the min and max value of each parameter
+ minPar[iplane][iside][ipar] = calcMinOrMax(twFitPar[iplane][iside][ipar], iplane, "min");
+ maxPar[iplane][iside][ipar] = calcMinOrMax(twFitPar[iplane][iside][ipar], iplane, "max");
+ } // Parameter loop
+ return;
+} // calcParAvg()
+
+// Draw the time-walk fit parameter multi-graph
+void drawParams(UInt_t iplane) {
+ // Populate the parameter canvases
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ twFitParCan[ipar]->cd(iplane+1);
+ twFitParMultiGraph[iplane][ipar]->Draw("AP");
+ gPad->Modified(); gPad->Update();
+ // Add aesthetics
+ modifyMultiGraph(twFitParMultiGraph[iplane][ipar], "Paddle Number", "Fit Parameter "+twFitParNames[ipar]);
+ twFitParMultiGraph[iplane][ipar]->GetXaxis()->SetNdivisions(21);
+ twFitParMultiGraph[iplane][ipar]->Draw("AP");
+ // Draw legends
+ twFitParLegend[iplane][ipar] = new TLegend(0.1, 0.7, 0.3, 0.9);
+ avgParLegend[iplane][ipar] = new TLegend(0.6, 0.7, 0.9, 0.9);
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ if (iside == 0) twFitParLegend[iplane][ipar]->AddEntry(twFitParGraph[iplane][iside][ipar], "Positive Side", "p");
+ if (iside == 1) twFitParLegend[iplane][ipar]->AddEntry(twFitParGraph[iplane][iside][ipar], "Negative Side", "p");
+ avgParLegend[iplane][ipar]->AddEntry(avgParFit[iplane][iside][ipar], "#bar{"+twFitParNames[ipar]+"}"+Form(" = %.2f", avgPar[iplane][iside][ipar]),"l");
+ // Define the time-walk fit parameter multigraph y-axis range
+ if (iside == 0) {
+ // Set y-axis range min
+ if (minPar[iplane][iside][ipar] < 0.0 || minPar[iplane][iside+1][ipar] < 0.0) {
+ if (minPar[iplane][iside][ipar] < minPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside][ipar] + minScale*minPar[iplane][iside][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside+1][ipar] + minScale*minPar[iplane][iside+1][ipar]);
+ }
+ if (minPar[iplane][iside][ipar] > 0.0 || minPar[iplane][iside+1][ipar] > 0.0) {
+ if (minPar[iplane][iside][ipar] < minPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside][ipar] - minScale*minPar[iplane][iside][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside+1][ipar] - minScale*minPar[iplane][iside+1][ipar]);
+ }
+ // Set y-axis range max
+ if (maxPar[iplane][iside][ipar] < 0.0 || maxPar[iplane][iside+1][ipar] < 0.0) {
+ if (maxPar[iplane][iside][ipar] < maxPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside+1][ipar] - maxScale*maxPar[iplane][iside+1][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside][ipar] - maxScale*maxPar[iplane][iside][ipar]);
+ }
+ if (maxPar[iplane][iside][ipar] > 0.0 || maxPar[iplane][iside+1][ipar] > 0.0) {
+ if (maxPar[iplane][iside][ipar] < maxPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside+1][ipar] + maxScale*maxPar[iplane][iside+1][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside][ipar] + maxScale*maxPar[iplane][iside][ipar]);
+ }
+ } // Side index = 0 condition
+ } // Side loop
+ twFitParLegend[iplane][ipar]->Draw();
+ avgParLegend[iplane][ipar]->Draw();
+ } // Parameter loop
+ return;
+} // drawParams
+
+//=:=:=:=:=
+//=: Main
+//=:=:=:=:=
+
+void timeWalkCalib() {
+
+ // ROOT settings
+ gStyle->SetTitleFontSize(fontSize);
+ gStyle->SetLabelSize(fontSize, "XY");
+ gStyle->SetTitleSize(fontSize, "XY");
+ gStyle->SetTitleOffset(yTitleOffset, "Y");
+ gStyle->SetStatFormat(".2f");
+ gStyle->SetOptFit(0);
+ gStyle->SetOptStat(0);
+
+ // Read the ROOT file containing the time-walk histos
+ histoFile = new TFile("timeWalkHistos.root", "READ");
+ // Obtain the top level directory
+ dataDir = dynamic_cast <TDirectory*> (histoFile->FindObjectAny("hodoUncalib"));
+ // Create the parameter canvases
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++)
+ twFitParCan[ipar] = makeCan(2, 2, 1600, 800, twFitParCan[ipar], twFitParNames[ipar]+"FitParCan", "Parameter "+twFitParNames[ipar]+" Canvas");
+ // Loop over the planes
+ for(UInt_t iplane = 0; iplane < nPlanes; iplane++) {
+ //for(UInt_t iplane = 0; iplane < 1; iplane++) {
+ // Obtain the plane directory
+ planeDir[iplane] = dynamic_cast <TDirectory*> (dataDir->FindObjectAny(planeNames[iplane]));
+ // Create multigraphs
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++)
+ twFitParMultiGraph[iplane][ipar] = new TMultiGraph(planeNames[iplane]+"_"+twFitParNames[ipar]+"_Multigraph", "Plane "+planeNames[iplane]+" Parameter "+twFitParNames[ipar]);
+ // Loop over the sides
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ // Obtain the side and time walk directories
+ sideDir[iplane][iside] = dynamic_cast <TDirectory*> (planeDir[iplane]->FindObjectAny(sideNames[iside]));
+ twDir[iplane][iside] = dynamic_cast <TDirectory*> (sideDir[iplane][iside]->FindObjectAny("adcTdcTimeDiffWalk"));
+ // Create the time-walk histo and fit canvases
+ if (planeNames[iplane] != "1y" || planeNames[iplane] != "2y") twFitCan[iplane][iside] = makeCan(4, 4, 1600, 800, twFitCan[iplane][iside], planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan", planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan");
+ if (planeNames[iplane] == "1y" || planeNames[iplane] == "2y") twFitCan[iplane][iside] = makeCan(4, 3, 1600, 800, twFitCan[iplane][iside], planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan", planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan");
+ // Loop over the paddles
+ for(UInt_t ipaddle = 0; ipaddle < nbars[iplane]; ipaddle++) {
+ // Populate the paddle index arrays
+ paddleIndex[iplane][iside][ipaddle] = Double_t (ipaddle + 1);
+ // Obtain the time-walk histos
+ h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle] = dynamic_cast <TH2F*> (twDir[iplane][iside]->FindObjectAny(Form("h2_adcTdcTimeDiffWalk_paddle_%d", ipaddle+1)));
+ // Perform the time-walk fits while ignoring the S2Y plane (for now)
+ // if (planeNames[iplane] == "2y") continue;
+ doTwFits(iplane, iside, ipaddle);
+ } // Paddle loop
+ // Produce the time-walk fit parameter graphs
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ // Populate graphs and multi-graphs
+ twFitParGraph[iplane][iside][ipar] = new TGraph(nbars[iplane], paddleIndex[iplane][iside], twFitPar[iplane][iside][ipar]);
+ if (iside == 0) addColorToGraph(22, markerSize, kRed, twFitParGraph[iplane][iside][ipar]);
+ if (iside == 1) addColorToGraph(23, markerSize, kBlue, twFitParGraph[iplane][iside][ipar]);
+ // // Add graphs to multi graph
+ // twFitParMultiGraph[iplane][ipar]->Add(twFitParGraph[iplane][iside][ipar]);
+ } // Parameter loop
+ // Calculate the average of the time-walk fit parameters
+ calcParAvg(iplane, iside);
+ } // Side loop
+ // Draw the time-walk parameter graphs
+ drawParams(iplane);
+ } // Plane loop
+ return;
+} // timeWalkFits()
+
+
+
+
+
+
+
diff --git a/CALIBRATION/hms_hodo_calib/timeWalkHistos.C b/CALIBRATION/hms_hodo_calib/timeWalkHistos.C
new file mode 100644
index 0000000000000000000000000000000000000000..5255c4fae7e08e2d822c3339b80f4be01e02c4a4
--- /dev/null
+++ b/CALIBRATION/hms_hodo_calib/timeWalkHistos.C
@@ -0,0 +1,427 @@
+// Macro to generate time-walk histograms for the hodoscopes
+// Author: Eric Pooser, pooser@jlab.org
+
+#include <time.h>
+
+// Declare replay data file and output file
+TFile *replayFile, *outFile;
+
+// Declare data tree
+TTree *rawDataTree;
+
+// Declare constants
+static const UInt_t nPlanes = 4;
+static const UInt_t nSides = 2;
+static const UInt_t nAdcSignals = 3;
+static const UInt_t nTdcSignals = 1;
+static const UInt_t nSignals = nAdcSignals + nTdcSignals;
+static const UInt_t maxTdcHits = 128;
+static const UInt_t maxAdcHits = 4;
+static const UInt_t nMaxBars = 16;
+
+static const TString spec = "H";
+static const TString detec = "hod";
+static const TString trig = "T";
+static const TString hms = "hms";
+
+static const UInt_t nbars[nPlanes] = {16, 10, 16, 10};
+static const TString planeNames[nPlanes] = {"1x", "1y", "2x", "2y"};
+static const TString sideNames[nSides] = {"pos", "neg"};
+static const TString signalNames[nSignals] = {"Adc", "Tdc"};
+static const TString adcData[nAdcSignals] = {"ErrorFlag", "PulseTimeRaw", "PulseAmp"};
+static const TString tdcData[nTdcSignals] = {"TimeRaw"};
+
+static const Double_t tdcChanToTime = 0.100; // Units of ns
+static const Double_t adcChanToTime = 0.0625; // Units of ns
+static const Double_t adcDynamicRange = 1000.0; // Units of mV
+static const Double_t nAdcChan = 4096.0; // Units of ADC channels
+static const Double_t adcChanTomV = adcDynamicRange/nAdcChan; // Units of mV/ADC Chan
+
+static const Double_t hodoPulseAmpCutLow = 10.0; // Units of mV
+static const Double_t hodoPulseAmpCutHigh = 1000.0; // Units of mV
+static const Double_t refAdcPulseAmpCutLow = 50.0; // Units of mV
+static const Double_t refAdcPulseAmpCutHigh = 60.0; // Units of mV
+static const Double_t refAdcPulseTimeCutLow = 210.0; // Units of ns
+static const Double_t refAdcPulseTimeCutHigh = 225.0; // Units of ns
+static const Double_t adcTdcTimeDiffCutLow = -100.0; // Units of ns
+static const Double_t adcTdcTimeDiffCutHigh = 100.0; // Units of ns
+static const Double_t calEtotCutVal = 0.100; // Units of GeV
+static const Double_t cerNpeSumCutVal = 1.5; // Units of NPE
+// static const Double_t adcTdcTimeDiffCutLow = -6000.0; // Units of ns
+// static const Double_t adcTdcTimeDiffCutHigh = 1000.0; // Units of ns
+
+static const Double_t fontSize = 0.04;
+
+// Declare variables to obtain from data tree
+// Number of ADC & TDC hits
+Int_t adcHits[nPlanes][nSides][nSignals], tdcHits[nPlanes][nSides][nSignals];
+// Paddle number which fired
+Double_t adcPaddle[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+Double_t tdcPaddle[nPlanes][nSides][nSignals][nMaxBars*maxTdcHits];
+// FADC data
+Double_t hodoAdcPulseTimeRaw[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+Double_t hodoAdcPulseAmp[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+Double_t hodoAdcErrorFlag[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+// TDC Data
+Double_t hodoTdcTimeRaw[nPlanes][nSides][nSignals][nMaxBars*maxTdcHits];
+// Trigger apparatus data
+Double_t refAdcPulseTimeRaw, refAdcPulseAmp, refAdcMultiplicity;
+Double_t refT1TdcTimeRaw, refT2TdcTimeRaw;
+// Declare variables to steer data to arrays
+TString *adcBaseName, *adcNdataName, *adcPaddleName;
+TString *adcErrorFlagName, *adcPulseTimeRawName, *adcPulseAmpName;
+TString *tdcBaseName, *tdcNdataName, *tdcPaddleName;
+TString *tdcTimeRawName;
+// Declare directories for output ROOT file
+TDirectory *trigRawDir, *hodoRawDir, *hodoUncalibDir;
+TDirectory *planeRawDir[nPlanes], *sideRawDir[nPlanes][nSides], *planeUncalibDir[nPlanes], *sideUncalibDir[nPlanes][nSides];
+TDirectory *adcTimeWalkDir[nPlanes][nSides], *tdcTimeWalkDir[nPlanes][nSides], *adcTdcTimeDiffWalkDir[nPlanes][nSides];
+// Declare histos
+// 1D histos
+TH1F *h1_refAdcPulseTimeRaw, *h1_refAdcPulseAmp, *h1_refAdcMultiplicity;
+TH1F *h1_refT1TdcTimeRaw, *h1_refT2TdcTimeRaw;
+// 2D histos
+TH2F *h2_adcErrorFlag[nPlanes][nSides];
+TH2F *h2_adcPulseTimeRaw[nPlanes][nSides];
+TH2F *h2_adcPulseAmp[nPlanes][nSides];
+TH2F *h2_tdcTimeRaw[nPlanes][nSides];
+
+TH2F *h2_adcPulseTime[nPlanes][nSides];
+TH2F *h2_tdcTime[nPlanes][nSides];
+TH2F *h2_adcTdcTimeDiff[nPlanes][nSides];
+TH2F *h2_adcPulseAmpCuts[nPlanes][nSides];
+TH2F *h2_adcTimeWalk[nPlanes][nSides][nMaxBars];
+TH2F *h2_tdcTimeWalk[nPlanes][nSides][nMaxBars];
+TH2F *h2_adcTdcTimeDiffWalk[nPlanes][nSides][nMaxBars];
+
+// Declare event data objects
+UInt_t nentries, ievent, adcPaddleNum, tdcPaddleNum;
+
+Int_t numAdcHits, numTdcHits;
+
+Double_t adcErrorFlag, adcPulseTimeRaw, adcPulseTime, adcPulseAmp;
+Double_t tdcTimeRaw, tdcTime, adcTdcTimeDiff;
+Double_t calEtot, cerNpeSum;
+
+Bool_t adcRefMultiplicityCut, adcRefPulseAmpCut, adcRefPulseTimeCut;
+Bool_t edtmCut, adcErrorFlagCut, adcAndTdcHitCut;
+Bool_t adcPulseAmpCut, adcTdcTimeDiffCut;
+Bool_t calEtotCut, cerNpeSumCut;
+
+void generatePlots(UInt_t iplane, UInt_t iside, UInt_t ipaddle) {
+
+ // Create trigger aparatus directory
+ trigRawDir = dynamic_cast <TDirectory*> (outFile->Get("trigAppRaw"));
+ if (!trigRawDir) {trigRawDir = outFile->mkdir("trigAppRaw"); trigRawDir->cd();}
+ else outFile->cd("trigAppRaw");
+ // FADC reference
+ if (!h1_refAdcPulseTimeRaw) h1_refAdcPulseTimeRaw = new TH1F("h1_refAdcPulseTimeRaw", "ROC1 Raw FADC Reference Pulse Time; Raw FADC Pulse Time (ns); Number of Entries / 100 ps", 4000, 0, 400);
+ if (!h1_refAdcPulseAmp) h1_refAdcPulseAmp = new TH1F("h1_refAdcPulseAmp", "ROC1 FADC Reference Pulse Amplitude; FADC Pulse Amplitude (mV); Number of Entries / 1 mV", 1000, 0, 1000);
+ if (!h1_refAdcMultiplicity) h1_refAdcMultiplicity = new TH1F("h1_refAdcMultiplicity", "ROC1 FADC Reference Multiplicity; Raw FADC Multiplicity; Number of Entries", 5, -0.5, 4.5);
+ // TDC reference
+ if (!h1_refT1TdcTimeRaw) h1_refT1TdcTimeRaw = new TH1F("h1_refT1TdcTimeRaw", "ROC1 T1 Raw TDC Reference TDC Time (Slot 20, Channel 15); Raw TDC Time (ns); Number of Entries / 100 ps", 6000, 0, 600);
+ if (!h1_refT2TdcTimeRaw) h1_refT2TdcTimeRaw = new TH1F("h1_refT2TdcTimeRaw", "ROC1 T2 Raw TDC Reference TDC Time (Slot 19, Channel 31); Raw TDC Time (ns); Number of Entries / 100 ps", 6000, 0, 600);
+
+ // Create plane directory for raw hodoscope quantities
+ hodoRawDir = dynamic_cast <TDirectory*> (outFile->Get("hodoRaw"));
+ if (!hodoRawDir) {hodoRawDir = outFile->mkdir("hodoRaw"); hodoRawDir->cd();}
+ else outFile->cd("hodoRaw");
+ planeRawDir[iplane] = dynamic_cast <TDirectory*> (hodoRawDir->Get(planeNames[iplane]));
+ if (!planeRawDir[iplane]) {planeRawDir[iplane] = hodoRawDir->mkdir(planeNames[iplane]); planeRawDir[iplane]->cd();}
+ else hodoRawDir->cd(planeNames[iplane]);
+ // Create side directory
+ sideRawDir[iplane][iside] = dynamic_cast <TDirectory*> (planeRawDir[iplane]->Get(sideNames[iside]));
+ if (!sideRawDir[iplane][iside]) {sideRawDir[iplane][iside] = planeRawDir[iplane]->mkdir(sideNames[iside]); sideRawDir[iplane][iside]->cd();}
+ else outFile->cd("hodoRaw/"+planeNames[iplane]+"/"+sideNames[iside]);
+ // Book histos
+ if (!h2_adcErrorFlag[iplane][iside]) h2_adcErrorFlag[iplane][iside] = new TH2F("h2_adcErrorFlag", "FADC Error Flag Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Error Flag", nbars[iplane], 0.5, nbars[iplane]+0.5, 2, -0.5, 1.5);
+ if (!h2_adcPulseTimeRaw[iplane][iside]) h2_adcPulseTimeRaw[iplane][iside] = new TH2F("h2_adcPulseTimeRaw", "Raw FADC Pulse Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; Raw FADC Pulse Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 4000, 0, 400);
+ if (!h2_adcPulseAmp[iplane][iside]) h2_adcPulseAmp[iplane][iside] = new TH2F("h2_adcPulseAmp", "FADC Pulse Amplitude Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Pulse Amplitude (mV) / 1 mV", nbars[iplane], 0.5, nbars[iplane]+0.5, 1000, 0, 1000);
+ if (!h2_tdcTimeRaw[iplane][iside]) h2_tdcTimeRaw[iplane][iside] = new TH2F("h2_tdcTimeRaw", "Raw TDC Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; Raw TDC Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 4000, 0, 400);
+
+ // Create plane directory for uncalibrated hodoscope quantities
+ hodoUncalibDir = dynamic_cast <TDirectory*> (outFile->Get("hodoUncalib"));
+ if (!hodoUncalibDir) {hodoUncalibDir = outFile->mkdir("hodoUncalib"); hodoUncalibDir->cd();}
+ else outFile->cd("hodoUncalib");
+ planeUncalibDir[iplane] = dynamic_cast <TDirectory*> (hodoUncalibDir->Get(planeNames[iplane]));
+ if (!planeUncalibDir[iplane]) {planeUncalibDir[iplane] = hodoUncalibDir->mkdir(planeNames[iplane]); planeUncalibDir[iplane]->cd();}
+ else hodoUncalibDir->cd(planeNames[iplane]);
+ // Create side directory
+ sideUncalibDir[iplane][iside] = dynamic_cast <TDirectory*> (planeUncalibDir[iplane]->Get(sideNames[iside]));
+ if (!sideUncalibDir[iplane][iside]) {sideUncalibDir[iplane][iside] = planeUncalibDir[iplane]->mkdir(sideNames[iside]); sideUncalibDir[iplane][iside]->cd();}
+ else outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]);
+ // Book histos
+ if (!h2_adcPulseTime[iplane][iside]) h2_adcPulseTime[iplane][iside] = new TH2F("h2_adcPulseTime", "FADC Pulse Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Pulse Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 8000, -400, 400);
+ if (!h2_tdcTime[iplane][iside]) h2_tdcTime[iplane][iside] = new TH2F("h2_tdcTime", "TDC Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; TDC Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 8000, -400, 400);
+ if (!h2_adcTdcTimeDiff[iplane][iside]) h2_adcTdcTimeDiff[iplane][iside] = new TH2F("h2_adcTdcTimeDiff", "TDC-ADC Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; TDC - ADC Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 16000, -800, 800);
+ if (!h2_adcPulseAmpCuts[iplane][iside]) h2_adcPulseAmpCuts[iplane][iside] = new TH2F("h2_adcPulseAmpCuts", "FADC Pulse Amplitude Cuts Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Pulse Amplitude (mV) / 1 mV", nbars[iplane], 0.5, nbars[iplane]+0.5, 1000, 0, 1000);
+ // Create ADC time-walk directory
+ adcTimeWalkDir[iplane][iside] = dynamic_cast <TDirectory*> (sideUncalibDir[iplane][iside]->Get("adcTimeWalk"));
+ if (!adcTimeWalkDir[iplane][iside]) {adcTimeWalkDir[iplane][iside] = sideUncalibDir[iplane][iside]->mkdir("adcTimeWalk"); adcTimeWalkDir[iplane][iside]->cd();}
+ else (outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]+"/adcTimeWalk"));
+ // Book histos
+ if (!h2_adcTimeWalk[iplane][iside][ipaddle]) h2_adcTimeWalk[iplane][iside][ipaddle] = new TH2F(Form("h2_adcTimeWalk_paddle_%d", ipaddle+1), "ADC Time vs. Pulse Amp Plane "+planeNames[iplane]+" Side "+sideNames[iside]+Form(" Paddle %d", ipaddle+1)+"; Pulse Amplitude (mV) / 1 mV; ADC Time (ns) / 100 ps", 1000, 0, 1000, 4000, -200, 200);
+ // Create TDC time-walk directory
+ tdcTimeWalkDir[iplane][iside] = dynamic_cast <TDirectory*> (sideUncalibDir[iplane][iside]->Get("tdcTimeWalk"));
+ if (!tdcTimeWalkDir[iplane][iside]) {tdcTimeWalkDir[iplane][iside] = sideUncalibDir[iplane][iside]->mkdir("tdcTimeWalk"); tdcTimeWalkDir[iplane][iside]->cd();}
+ else (outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]+"/tdcTimeWalk"));
+ // Book histos
+ if (!h2_tdcTimeWalk[iplane][iside][ipaddle]) h2_tdcTimeWalk[iplane][iside][ipaddle] = new TH2F(Form("h2_tdcTimeWalk_paddle_%d", ipaddle+1), "TDC Time vs. Pulse Amp Plane "+planeNames[iplane]+" Side "+sideNames[iside]+Form(" Paddle %d", ipaddle+1)+"; Pulse Amplitude (mV) / 1 mV; TDC Time (ns) / 100 ps", 1000, 0, 1000, 4000, -300, 100);
+ // Create TDC-ADC time-walk directory
+ adcTdcTimeDiffWalkDir[iplane][iside] = dynamic_cast <TDirectory*> (sideUncalibDir[iplane][iside]->Get("adcTdcTimeDiffWalk"));
+ if (!adcTdcTimeDiffWalkDir[iplane][iside]) {adcTdcTimeDiffWalkDir[iplane][iside] = sideUncalibDir[iplane][iside]->mkdir("adcTdcTimeDiffWalk"); adcTdcTimeDiffWalkDir[iplane][iside]->cd();}
+ else (outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]+"/adcTdcTimeDiffWalk"));
+ // Book histos
+ if (!h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]) h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle] = new TH2F(Form("h2_adcTdcTimeDiffWalk_paddle_%d", ipaddle+1), "TDC-ADC Time vs. Pulse Amp Plane "+planeNames[iplane]+" Side "+sideNames[iside]+Form(" Paddle %d", ipaddle+1)+"; Pulse Amplitude (mV) / 1 mV; TDC-ADC Time (ns) / 100 ps", 1000, 0, 1000, 150, -20, -5);
+
+} // generatePlots()
+
+void timeWalkHistos(Int_t runNum) {
+
+ // Global ROOT settings
+ gStyle->SetOptFit();
+ gStyle->SetStatFormat(".2f");
+ gStyle->SetTitleFontSize(fontSize);
+ gStyle->SetLabelSize(fontSize, "XY");
+ gStyle->SetTitleSize(fontSize, "XY");
+
+ // Initialize the analysis clock
+ clock_t t;
+ t = clock();
+
+ // Obtain the replay data file and create new output ROOT file
+ // replayFile = new TFile("ROOTfiles/hms_replay_production_all_1267_-1.root", "READ");
+ // replayFile = new TFile("ROOTfiles/hms_replay_production_all_1268_-1.root", "READ");
+ replayFile = new TFile("ROOTfiles/hms_replay_production_all_1267_1268_-1.root", "READ");
+
+ // replayFile = new TFile(Form("ROOTfiles/hms_replay_production_all_%d_-1.root", runNum), "READ");
+ // replayFile = new TFile(Form("ROOTfiles/hms_coin_replay_production_%d_-1.root", runNum), "READ");
+
+ outFile = new TFile("timeWalkHistos_temp.root", "RECREATE");
+ // Obtain the tree
+ rawDataTree = dynamic_cast <TTree*> (replayFile->Get("T"));
+ // Acquire the trigger apparatus data
+ rawDataTree->SetBranchAddress("T.hms.hFADC_TREF_ROC1_adcPulseTimeRaw", &refAdcPulseTimeRaw);
+ rawDataTree->SetBranchAddress("T.hms.hFADC_TREF_ROC1_adcPulseAmp", &refAdcPulseAmp);
+ rawDataTree->SetBranchAddress("T.hms.hFADC_TREF_ROC1_adcMultiplicity", &refAdcMultiplicity);
+ rawDataTree->SetBranchAddress("T.hms.hT1_tdcTimeRaw", &refT1TdcTimeRaw);
+ rawDataTree->SetBranchAddress("T.hms.hT2_tdcTimeRaw", &refT2TdcTimeRaw);
+ rawDataTree->SetBranchAddress("H.cal.etot", &calEtot);
+ rawDataTree->SetBranchAddress("H.cer.npeSum", &cerNpeSum);
+ // Loop over the planes, sides, signals, leafs, and fill data arrays
+ for(UInt_t iplane = 0; iplane < nPlanes; iplane++) {
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ // Generate directory structure and histograms
+ for(UInt_t ipaddle = 0; ipaddle < nbars[iplane]; ipaddle++)
+ generatePlots(iplane, iside, ipaddle);
+ for(UInt_t isignal = 0; isignal < nSignals; isignal++) {
+
+ // Acquire the hodoscope ADC data objects
+ if(signalNames[isignal] == "Adc") {
+ // Define the leaf strings of interest
+ adcBaseName = new TString(spec+"."+detec+"."+planeNames[iplane]+"."+
+ sideNames[iside]+signalNames[isignal]);
+ adcNdataName = new TString("Ndata."+*adcBaseName+"Counter");
+ adcPaddleName = new TString(*adcBaseName+"Counter");
+ // Fill the leafs
+ rawDataTree->SetBranchAddress(*adcNdataName, &adcHits[iplane][iside][isignal]);
+ rawDataTree->SetBranchAddress(*adcPaddleName, &adcPaddle[iplane][iside][isignal]);
+ // Loop over the adc data objects
+ for(UInt_t iadcdata = 0; iadcdata < nAdcSignals; iadcdata++) {
+ if (adcData[iadcdata] == "ErrorFlag") {
+ adcErrorFlagName = new TString(*adcBaseName+"ErrorFlag");
+ rawDataTree->SetBranchAddress(*adcErrorFlagName, &hodoAdcErrorFlag[iplane][iside][isignal]);
+ } // Error flag leaf
+ if (adcData[iadcdata] == "PulseTimeRaw") {
+ adcPulseTimeRawName = new TString(*adcBaseName+"PulseTimeRaw");
+ rawDataTree->SetBranchAddress(*adcPulseTimeRawName, &hodoAdcPulseTimeRaw[iplane][iside][isignal]);
+ } // Raw pulse time leaf
+ if (adcData[iadcdata] == "PulseAmp") {
+ adcPulseAmpName = new TString(*adcBaseName+"PulseAmp");
+ rawDataTree->SetBranchAddress(*adcPulseAmpName, &hodoAdcPulseAmp[iplane][iside][isignal]);
+ } // Pedestal subtracted pulse amplitude leaf
+ } // ADC signal data loop
+ } // ADC signal
+
+ // Acquire the hodoscope TDC data objects
+ if(signalNames[isignal] == "Tdc") {
+ // Define the leaf strings of interest
+ tdcBaseName = new TString(spec+"."+detec+"."+planeNames[iplane]+"."+
+ sideNames[iside]+signalNames[isignal]);
+ tdcNdataName = new TString("Ndata."+*tdcBaseName+"Counter");
+ tdcPaddleName = new TString(*tdcBaseName+"Counter");
+ // Fill the leafs
+ rawDataTree->SetBranchAddress(*tdcNdataName, &tdcHits[iplane][iside][isignal]);
+ rawDataTree->SetBranchAddress(*tdcPaddleName, &tdcPaddle[iplane][iside][isignal]);
+ // Loop over the TDC data objects
+ for(UInt_t itdcdata = 0; itdcdata < nTdcSignals; itdcdata++) {
+ if (tdcData[itdcdata] == "TimeRaw") {
+ tdcTimeRawName = new TString(*tdcBaseName+"TimeRaw");
+ rawDataTree->SetBranchAddress(*tdcTimeRawName, &hodoTdcTimeRaw[iplane][iside][isignal]);
+ } // Raw TDC time leaf
+ } // TDC signal data loop
+ } // TDC signal
+
+ } // Signal loop
+ } // Side loop
+ } // Plane loop
+
+ // Loop over the events and fill histograms
+ nentries = rawDataTree->GetEntries();
+ // nentries = 100000;
+ cout << "\n******************************************" << endl;
+ cout << nentries << " Events Will Be Processed" << endl;
+ cout << "******************************************\n" << endl;
+ for (ievent = 0; ievent < nentries; ievent++) {
+ // Acquire the event from the data tree
+ rawDataTree->GetEntry(ievent);
+ // Fiducial PID cuts
+ calEtotCut = (calEtot < calEtotCutVal);
+ cerNpeSumCut = (cerNpeSum < cerNpeSumCutVal);
+ if (calEtotCut || cerNpeSumCut) continue;
+ // Fill trigger apparatus histos
+ h1_refAdcPulseTimeRaw->Fill(refAdcPulseTimeRaw*adcChanToTime);
+ h1_refAdcPulseAmp->Fill(refAdcPulseAmp);
+ h1_refAdcMultiplicity->Fill(refAdcMultiplicity);
+ h1_refT1TdcTimeRaw->Fill(refT1TdcTimeRaw*tdcChanToTime);
+ h1_refT2TdcTimeRaw->Fill(refT2TdcTimeRaw*tdcChanToTime);
+ // Loop over the planes, sides, signals, data arrays, and fill hodoscope histograms
+ for(UInt_t iplane = 0; iplane < nPlanes; iplane++) {
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ for(UInt_t isignal = 0; isignal < nSignals; isignal++) {
+
+ // Acquire the hodoscope ADC data objects
+ if(signalNames[isignal] == "Adc") {
+ // Loop over the ADC data objects
+ for(UInt_t iadcdata = 0; iadcdata < nAdcSignals; iadcdata++) {
+ numAdcHits = adcHits[iplane][iside][isignal];
+ for (Int_t iadchit = 0; iadchit < numAdcHits; iadchit++) {
+ // Obtain variables
+ adcPaddleNum = adcPaddle[iplane][iside][isignal][iadchit];
+ adcErrorFlag = hodoAdcErrorFlag[iplane][iside][isignal][iadchit];
+ adcPulseTimeRaw = hodoAdcPulseTimeRaw[iplane][iside][isignal][iadchit]*adcChanToTime;
+ adcPulseTime = adcPulseTimeRaw - refAdcPulseTimeRaw*adcChanToTime;
+ adcPulseAmp = hodoAdcPulseAmp[iplane][iside][isignal][iadchit];
+ // Fill histos
+ if (adcData[iadcdata] == "ErrorFlag")
+ h2_adcErrorFlag[iplane][iside]->Fill(adcPaddleNum, adcErrorFlag);
+ if (adcData[iadcdata] == "PulseTimeRaw") {
+ h2_adcPulseTimeRaw[iplane][iside]->Fill(adcPaddleNum, adcPulseTimeRaw);
+ h2_adcPulseTime[iplane][iside]->Fill(adcPaddleNum, adcPulseTime);
+ }
+ if (adcData[iadcdata] == "PulseAmp")
+ h2_adcPulseAmp[iplane][iside]->Fill(adcPaddleNum, adcPulseAmp);
+ } // ADC hit loop
+ } // ADC signal data loop
+ } // ADC signal
+
+ // Acquire the hodoscope TDC data objects
+ if(signalNames[isignal] == "Tdc") {
+ // Loop over the TDC data objects
+ for(UInt_t itdcdata = 0; itdcdata < nTdcSignals; itdcdata++) {
+ numTdcHits = tdcHits[iplane][iside][isignal];
+ // Define cuts
+ // edtmCut = (numTdcHits == nbars[0] || numTdcHits == nbars[2] || numTdcHits == nbars[3]);
+ // Implement cuts
+ // if (edtmCut) continue;
+ for (Int_t itdchit = 0; itdchit < numTdcHits; itdchit++) {
+ // Obtain variables
+ tdcPaddleNum = tdcPaddle[iplane][iside][isignal][itdchit];
+ tdcTimeRaw = hodoTdcTimeRaw[iplane][iside][isignal][itdchit]*tdcChanToTime;
+ tdcTime = tdcTimeRaw - refT2TdcTimeRaw*tdcChanToTime;
+ if (tdcData[itdcdata] == "TimeRaw") {
+ h2_tdcTimeRaw[iplane][iside]->Fill(tdcPaddleNum, tdcTimeRaw);
+ h2_tdcTime[iplane][iside]->Fill(tdcPaddleNum, tdcTime);
+ } // RDC raw time signal
+ } // TDC hit loop
+ } // TDC signal data loop
+ } // TDC signal
+
+ // Define cuts
+ adcRefMultiplicityCut = (refAdcMultiplicity != 1.0);
+ adcRefPulseAmpCut = (refAdcPulseAmp < refAdcPulseAmpCutLow || refAdcPulseAmp > refAdcPulseAmpCutHigh);
+ adcRefPulseTimeCut = (refAdcPulseTimeRaw*adcChanToTime < refAdcPulseTimeCutLow || refAdcPulseTimeRaw*adcChanToTime > refAdcPulseTimeCutHigh);
+ // Implement cuts
+ // if (adcRefMultiplicityCut || adcRefPulseAmpCut || adcRefPulseTimeCut) continue;
+ // Acquire the hodoscope ADC data objects
+ if(signalNames[isignal] == "Adc") {
+ // Loop over the signals again
+ for(UInt_t jsignal = 0; jsignal < nSignals; jsignal++) {
+ // Acquire othe hodoscope TDC data objects
+ if(signalNames[jsignal] == "Tdc") {
+ // Loop over the ADC data objects
+ for(UInt_t iadcdata = 0; iadcdata < nAdcSignals; iadcdata++) {
+ if (adcData[iadcdata] != "PulseTimeRaw") continue;
+ // Acquire the number of ADC hits
+ numAdcHits = adcHits[iplane][iside][isignal];
+ for (Int_t iadchit = 0; iadchit < numAdcHits; iadchit++) {
+ // Obtain variables
+ adcPaddleNum = UInt_t (adcPaddle[iplane][iside][isignal][iadchit]);
+ adcErrorFlag = hodoAdcErrorFlag[iplane][iside][isignal][iadchit];
+ adcPulseTimeRaw = hodoAdcPulseTimeRaw[iplane][iside][isignal][iadchit]*adcChanToTime;
+ adcPulseTime = adcPulseTimeRaw - refAdcPulseTimeRaw*adcChanToTime;
+ adcPulseAmp = hodoAdcPulseAmp[iplane][iside][isignal][iadchit];
+ // Define cuts
+ adcErrorFlagCut = (adcErrorFlag != 0);
+ adcPulseAmpCut = (adcPulseAmp < hodoPulseAmpCutLow || adcPulseAmp > hodoPulseAmpCutHigh);
+ // Implement cuts
+ if (adcErrorFlagCut || adcPulseAmpCut) continue;
+ // Loop over the TDC data objects
+ for(UInt_t itdcdata = 0; itdcdata < nTdcSignals; itdcdata++) {
+ if (tdcData[itdcdata] != "TimeRaw") continue;
+ // Acquire the number of TDC hits
+ numTdcHits = tdcHits[iplane][iside][jsignal];
+ // Define cuts
+ //edtmCut = (numTdcHits == nbars[0] || numTdcHits == nbars[2] || numTdcHits == nbars[3]);
+ // Implement cuts
+ //if (edtmCut) continue;
+ for (Int_t itdchit = 0; itdchit < numTdcHits; itdchit++) {
+ // Obtain variables
+ tdcPaddleNum = UInt_t (tdcPaddle[iplane][iside][jsignal][itdchit]);
+ tdcTimeRaw = hodoTdcTimeRaw[iplane][iside][jsignal][itdchit]*tdcChanToTime;
+ tdcTime = tdcTimeRaw - refT2TdcTimeRaw*tdcChanToTime;
+ adcTdcTimeDiff = tdcTime - adcPulseTime;
+ // Define cuts
+ adcAndTdcHitCut = (adcPaddleNum != tdcPaddleNum);
+ adcTdcTimeDiffCut = (adcTdcTimeDiff < adcTdcTimeDiffCutLow || adcTdcTimeDiff > adcTdcTimeDiffCutHigh);
+ // Implement cuts
+ if (adcAndTdcHitCut || adcTdcTimeDiffCut) continue;
+ h2_adcPulseAmpCuts[iplane][iside]->Fill(tdcPaddleNum, adcPulseAmp);
+ h2_adcTdcTimeDiff[iplane][iside]->Fill(tdcPaddleNum, adcTdcTimeDiff);
+ h2_adcTimeWalk[iplane][iside][tdcPaddleNum-1]->Fill(adcPulseAmp, adcPulseTime);
+ h2_tdcTimeWalk[iplane][iside][tdcPaddleNum-1]->Fill(adcPulseAmp, tdcTime);
+ h2_adcTdcTimeDiffWalk[iplane][iside][tdcPaddleNum-1]->Fill(adcPulseAmp, adcTdcTimeDiff);
+ h2_adcTdcTimeDiffWalk[iplane][iside][tdcPaddleNum-1]->GetXaxis()->CenterTitle();
+ h2_adcTdcTimeDiffWalk[iplane][iside][tdcPaddleNum-1]->GetYaxis()->CenterTitle();
+ } // TDC hit loop
+ } // TDC signal loop
+ } // ADC hit loop
+ } // ADC signal loop
+ } // TDC signal
+ } // Nested signal loop
+ } // ADC signal
+
+ } // Signal loop
+ } // Side loop
+ } // Plane loop
+
+ if (ievent % 100000 == 0 && ievent != 0)
+ cout << ievent << " Events Have Been Processed..." << endl;
+
+ } // rawDataTree event loop
+
+ cout << "\n***************************************" << endl;
+ cout << ievent << " Events Were Processed" << endl;
+ cout << "***************************************\n" << endl;
+
+ // Calculate the analysis rate
+ t = clock() - t;
+ printf ("The Analysis Took %.1f seconds \n", ((float) t) / CLOCKS_PER_SEC);
+ printf ("The Analysis Event Rate Was %.3f kHz \n", (ievent + 1) / (((float) t) / CLOCKS_PER_SEC*1000.));
+
+ outFile->Write();
+ //outFile->Close();
+
+ return 0;
+
+} // time_walk_calib()
diff --git a/CALIBRATION/shms_hodo_calib/timeWalkCalib.C b/CALIBRATION/shms_hodo_calib/timeWalkCalib.C
new file mode 100644
index 0000000000000000000000000000000000000000..7bd0fa234dff6224558113897395a7d96d87aab2
--- /dev/null
+++ b/CALIBRATION/shms_hodo_calib/timeWalkCalib.C
@@ -0,0 +1,296 @@
+// Macro to perform time-walk fits and extract the calibration parameters
+// Author: Eric Pooser, pooser@jlab.org
+
+// Declare ROOT files
+TFile *histoFile, *outFile;
+
+// Declare constants
+static const UInt_t nPlanes = 4;
+static const UInt_t nSides = 2;
+static const UInt_t nBarsMax = 21;
+static const UInt_t nTwFitPars = 3;
+
+static const Double_t tdcThresh = 120.0; // 30 mV in units of FADC channels
+static const Double_t twFitRangeLow = 20.0;
+static const Double_t twFitRangeHigh = 400.0;
+static const Double_t c0twParInit = -35.0;
+static const Double_t c1twParInit = 5.0;
+static const Double_t c2twParInit = 0.25;
+
+static const Double_t fontSize = 0.05;
+static const Double_t yTitleOffset = 0.75;
+static const Double_t markerSize = 2.0;
+static const Double_t minScale = 0.75;
+static const Double_t maxScale = 0.75;
+
+static const UInt_t lineWidth = 4;
+static const UInt_t lineStyle = 7;
+
+static const UInt_t nbars[nPlanes] = {13, 13, 14, 21};
+static const TString planeNames[nPlanes] = {"1x", "1y", "2x", "2y"};
+static const TString sideNames[nSides] = {"pos", "neg"};
+static const TString twFitParNames[nTwFitPars] = {"c_{0}", "c_{1}", "c_{2}"};
+
+// Declare directories
+TDirectory *dataDir, *planeDir[nPlanes], *sideDir[nPlanes][nSides];
+TDirectory *twDir[nPlanes][nSides];
+// Declare fits
+TF1 *twFit[nPlanes][nSides][nBarsMax], *avgParFit[nPlanes][nSides][nTwFitPars];
+// Declare arrays
+Double_t paddleIndex[nPlanes][nSides][nBarsMax];
+Double_t twFitPar[nPlanes][nSides][nTwFitPars][nBarsMax], twFitParErr[nPlanes][nSides][nTwFitPars][nBarsMax];
+Double_t avgPar[nPlanes][nSides][nTwFitPars];
+Double_t minPar[nPlanes][nSides][nTwFitPars], maxPar[nPlanes][nSides][nTwFitPars];
+// Declare canvases
+TCanvas *twFitCan[nPlanes][nSides], *twFitParCan[nTwFitPars];
+// Declare histos
+// 2D histos
+TH2F *h2_adcTdcTimeDiffWalk[nPlanes][nSides][nBarsMax];
+// Declare graphs
+TGraph *twFitParGraph[nPlanes][nSides][nTwFitPars];
+// Declare multi-graphs
+TMultiGraph *twFitParMultiGraph[nPlanes][nTwFitPars];
+// Declare legends
+TLegend *twFitParLegend[nPlanes][nTwFitPars], *avgParLegend[nPlanes][nTwFitPars];
+// Declare text objects
+TPaveText *twFitParText[nPlanes][nSides][nBarsMax];
+
+// Functions to extract calibration parameters
+//=:=:=:=:=:=:
+//=: Level 1
+//=:=:=:=:=:=:
+
+// Add color to fit lines
+void addColorToFitLine(UInt_t style, UInt_t width, UInt_t color, TF1 *fit1D) {
+ fit1D->SetLineStyle(lineStyle);
+ fit1D->SetLineWidth(lineWidth);
+ fit1D->SetLineColor(color);
+ return;
+} // addColorToFitLine()
+
+// Add color to fit graphs
+void addColorToGraph(UInt_t style, UInt_t size, UInt_t color, TGraph *graph) {
+ graph->SetMarkerStyle(style);
+ graph->SetMarkerSize(size);
+ graph->SetMarkerColor(color);
+ return;
+} // addColorToFitGraph()
+
+// Add centered titles to multigraphs
+void modifyMultiGraph(TMultiGraph *mg, TString xtitle, TString ytitle) {
+ mg->GetYaxis()->SetTitle(ytitle);
+ mg->GetYaxis()->CenterTitle();
+ mg->GetXaxis()->SetTitle(xtitle);
+ mg->GetXaxis()->CenterTitle();
+ return;
+} // modifyMultiGraph()
+
+// Create canvases
+TCanvas *makeCan(UInt_t numColumns, UInt_t numRows, UInt_t winWidth, UInt_t winHeight, TCanvas *can, TString name, TString title) {
+ can = new TCanvas(name, title, winWidth, winHeight);
+ can->Divide(numColumns, numRows);
+ return can;
+} // makeCan()
+
+// Define the time-walk fit function
+Double_t twFitFunc(Double_t *a, Double_t *c) {
+ Double_t twFitVal = c[0] + c[1]/(TMath::Power((a[0]/tdcThresh), c[2]));
+ return twFitVal;
+} // twFitFunc()
+
+// Locate min or max value from input array
+Double_t calcMinOrMax(Double_t *array, UInt_t iplane, TString minOrmax) {
+ auto result = minmax_element(array, array+nbars[iplane]);
+ if (minOrmax == "min") return *result.first;
+ else if (minOrmax == "max") return *result.second;
+ else return 0.0;
+} // calcMinOrMax()
+
+//=:=:=:=:=:=:
+//=: Level 2
+//=:=:=:=:=:=:
+
+// Perform the timw-walk fits
+void doTwFits(UInt_t iplane, UInt_t iside, UInt_t ipaddle) {
+ // Draw fits on canvas
+ twFitCan[iplane][iside]->cd(ipaddle+1);
+ gPad->SetLogz();
+ h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]->Draw("COLZ");
+ gPad->Modified(); gPad->Update();
+ // Declare and initialize the fits
+ twFit[iplane][iside][ipaddle] = new TF1("twFit", twFitFunc, twFitRangeLow, twFitRangeHigh, nTwFitPars);
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++)
+ twFit[iplane][iside][ipaddle]->SetParName(ipar, twFitParNames[ipar]);
+ twFit[iplane][iside][ipaddle]->SetParameters(c0twParInit, c1twParInit, c2twParInit);
+ // Perform the fits and scream if it failed
+ Int_t twFitStatus = h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]->Fit("twFit", "REQ");
+ if (twFitStatus != 0)
+ cout << "ERROR: Time Walk Fit Failed!!! " << "Status = " << twFitStatus << " For Plane: " << planeNames[iplane] << " Side: " << sideNames[iside] << " Paddle: " << ipaddle+1 << endl;
+ // Create text box to display fir parameters
+ twFitParText[iplane][iside][ipaddle] = new TPaveText(0.4, 0.6, 0.895, 0.895, "NBNDC");
+ twFitParText[iplane][iside][ipaddle]->AddText(Form("Entries = %.0f", h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]->GetEntries()));
+ // Obtain the fit parameters and associated errors
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ twFitPar[iplane][iside][ipar][ipaddle] = twFit[iplane][iside][ipaddle]->GetParameter(ipar);
+ twFitParErr[iplane][iside][ipar][ipaddle] = twFit[iplane][iside][ipaddle]->GetParError(ipar);
+ twFitParText[iplane][iside][ipaddle]->AddText(Form(twFitParNames[ipar]+" = %.2f #pm %.2f", twFitPar[iplane][iside][ipar][ipaddle], twFitParErr[iplane][iside][ipar][ipaddle]));
+ } // Parameter loop
+ // Draw the fit parameter text
+ twFitParText[iplane][iside][ipaddle]->SetFillColor(kWhite);
+ twFitParText[iplane][iside][ipaddle]->SetTextAlign(12);
+ twFitParText[iplane][iside][ipaddle]->SetFillStyle(3003);
+ twFitParText[iplane][iside][ipaddle]->Draw();
+ gPad->Modified(); gPad->Update();
+ return;
+} // doTwFits()
+
+// Calculate the averege of the time-walk fit parameters
+void calcParAvg(UInt_t iplane, UInt_t iside) {
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ // Calculate the weighted average while ignoring fit errors provided by Minuit
+ avgParFit[iplane][iside][ipar] = new TF1("avgParFit", "pol0", 1, nbars[iplane]);
+ avgParFit[iplane][iside][ipar]->SetParName(0, "#bar{"+twFitParNames[ipar]+"}");
+ // Add color to fit lines
+ if (iside == 0) addColorToFitLine(lineStyle, lineWidth, kRed, avgParFit[iplane][iside][ipar]);
+ if (iside == 1) addColorToFitLine(lineStyle, lineWidth, kBlue, avgParFit[iplane][iside][ipar]);
+ // Initialize the parameters
+ if (ipar == 0) avgParFit[iplane][iside][ipar]->SetParameter(0, c0twParInit);
+ if (ipar == 1) avgParFit[iplane][iside][ipar]->SetParameter(0, c1twParInit);
+ if (ipar == 2) avgParFit[iplane][iside][ipar]->SetParameter(0, c2twParInit);
+ // Perform the least squares fit
+ Int_t avgParFitStatus = twFitParGraph[iplane][iside][ipar]->Fit("avgParFit", "REQ");
+ if (avgParFitStatus != 0)
+ cout << "ERROR: Parameter Fit Failed!!! " << "Status = " << avgParFitStatus << " For Plane: " << planeNames[iplane] << " Side: " << sideNames[iside] << endl;
+ // Store the fit parameters
+ avgPar[iplane][iside][ipar] = avgParFit[iplane][iside][ipar]->GetParameter(0);
+ // Add graphs to multi graph
+ twFitParMultiGraph[iplane][ipar]->Add(twFitParGraph[iplane][iside][ipar]);
+ // Calculate the min and max value of each parameter
+ minPar[iplane][iside][ipar] = calcMinOrMax(twFitPar[iplane][iside][ipar], iplane, "min");
+ maxPar[iplane][iside][ipar] = calcMinOrMax(twFitPar[iplane][iside][ipar], iplane, "max");
+ } // Parameter loop
+ return;
+} // calcParAvg()
+
+// Draw the time-walk fit parameter multi-graph
+void drawParams(UInt_t iplane) {
+ // Populate the parameter canvases
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ twFitParCan[ipar]->cd(iplane+1);
+ twFitParMultiGraph[iplane][ipar]->Draw("AP");
+ gPad->Modified(); gPad->Update();
+ // Add aesthetics
+ modifyMultiGraph(twFitParMultiGraph[iplane][ipar], "Paddle Number", "Fit Parameter "+twFitParNames[ipar]);
+ twFitParMultiGraph[iplane][ipar]->GetXaxis()->SetNdivisions(21);
+ twFitParMultiGraph[iplane][ipar]->Draw("AP");
+ // Draw legends
+ twFitParLegend[iplane][ipar] = new TLegend(0.1, 0.7, 0.3, 0.9);
+ avgParLegend[iplane][ipar] = new TLegend(0.6, 0.7, 0.9, 0.9);
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ if (iside == 0) twFitParLegend[iplane][ipar]->AddEntry(twFitParGraph[iplane][iside][ipar], "Positive Side", "p");
+ if (iside == 1) twFitParLegend[iplane][ipar]->AddEntry(twFitParGraph[iplane][iside][ipar], "Negative Side", "p");
+ avgParLegend[iplane][ipar]->AddEntry(avgParFit[iplane][iside][ipar], "#bar{"+twFitParNames[ipar]+"}"+Form(" = %.2f", avgPar[iplane][iside][ipar]),"l");
+ // Define the time-walk fit parameter multigraph y-axis range
+ if (iside == 0) {
+ // Set y-axis range min
+ if (minPar[iplane][iside][ipar] < 0.0 || minPar[iplane][iside+1][ipar] < 0.0) {
+ if (minPar[iplane][iside][ipar] < minPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside][ipar] + minScale*minPar[iplane][iside][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside+1][ipar] + minScale*minPar[iplane][iside+1][ipar]);
+ }
+ if (minPar[iplane][iside][ipar] > 0.0 || minPar[iplane][iside+1][ipar] > 0.0) {
+ if (minPar[iplane][iside][ipar] < minPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside][ipar] - minScale*minPar[iplane][iside][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMinimum(minPar[iplane][iside+1][ipar] - minScale*minPar[iplane][iside+1][ipar]);
+ }
+ // Set y-axis range max
+ if (maxPar[iplane][iside][ipar] < 0.0 || maxPar[iplane][iside+1][ipar] < 0.0) {
+ if (maxPar[iplane][iside][ipar] < maxPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside+1][ipar] - maxScale*maxPar[iplane][iside+1][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside][ipar] - maxScale*maxPar[iplane][iside][ipar]);
+ }
+ if (maxPar[iplane][iside][ipar] > 0.0 || maxPar[iplane][iside+1][ipar] > 0.0) {
+ if (maxPar[iplane][iside][ipar] < maxPar[iplane][iside+1][ipar])
+ twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside+1][ipar] + maxScale*maxPar[iplane][iside+1][ipar]);
+ else twFitParMultiGraph[iplane][ipar]->SetMaximum(maxPar[iplane][iside][ipar] + maxScale*maxPar[iplane][iside][ipar]);
+ }
+ } // Side index = 0 condition
+ } // Side loop
+ twFitParLegend[iplane][ipar]->Draw();
+ avgParLegend[iplane][ipar]->Draw();
+ } // Parameter loop
+ return;
+} // drawParams
+
+//=:=:=:=:=
+//=: Main
+//=:=:=:=:=
+
+void timeWalkCalib() {
+
+ // ROOT settings
+ gStyle->SetTitleFontSize(fontSize);
+ gStyle->SetLabelSize(fontSize, "XY");
+ gStyle->SetTitleSize(fontSize, "XY");
+ gStyle->SetTitleOffset(yTitleOffset, "Y");
+ gStyle->SetStatFormat(".2f");
+ gStyle->SetOptFit(0);
+ gStyle->SetOptStat(0);
+
+ // Read the ROOT file containing the time-walk histos
+ histoFile = new TFile("timeWalkHistos.root", "READ");
+ // Obtain the top level directory
+ dataDir = dynamic_cast <TDirectory*> (histoFile->FindObjectAny("hodoUncalib"));
+ // Create the parameter canvases
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++)
+ twFitParCan[ipar] = makeCan(2, 2, 1600, 800, twFitParCan[ipar], twFitParNames[ipar]+"FitParCan", "Parameter "+twFitParNames[ipar]+" Canvas");
+ // Loop over the planes
+ for(UInt_t iplane = 0; iplane < nPlanes; iplane++) {
+ //for(UInt_t iplane = 0; iplane < 1; iplane++) {
+ // Obtain the plane directory
+ planeDir[iplane] = dynamic_cast <TDirectory*> (dataDir->FindObjectAny(planeNames[iplane]));
+ // Create multigraphs
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++)
+ twFitParMultiGraph[iplane][ipar] = new TMultiGraph(planeNames[iplane]+"_"+twFitParNames[ipar]+"_Multigraph", "Plane "+planeNames[iplane]+" Parameter "+twFitParNames[ipar]);
+ // Loop over the sides
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ // Obtain the side and time walk directories
+ sideDir[iplane][iside] = dynamic_cast <TDirectory*> (planeDir[iplane]->FindObjectAny(sideNames[iside]));
+ twDir[iplane][iside] = dynamic_cast <TDirectory*> (sideDir[iplane][iside]->FindObjectAny("adcTdcTimeDiffWalk"));
+ // Create the time-walk histo and fit canvases
+ if (planeNames[iplane] != "2y") twFitCan[iplane][iside] = makeCan(5, 3, 1600, 800, twFitCan[iplane][iside], planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan", planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan");
+ if (planeNames[iplane] == "2y") twFitCan[iplane][iside] = makeCan(6, 4, 1600, 800, twFitCan[iplane][iside], planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan", planeNames[iplane]+"_"+sideNames[iside]+"_twFitCan");
+ // Loop over the paddles
+ for(UInt_t ipaddle = 0; ipaddle < nbars[iplane]; ipaddle++) {
+ // Populate the paddle index arrays
+ paddleIndex[iplane][iside][ipaddle] = Double_t (ipaddle + 1);
+ // Obtain the time-walk histos
+ h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle] = dynamic_cast <TH2F*> (twDir[iplane][iside]->FindObjectAny(Form("h2_adcTdcTimeDiffWalk_paddle_%d", ipaddle+1)));
+ // Perform the time-walk fits while ignoring the S2Y plane (for now)
+ // if (planeNames[iplane] == "2y") continue;
+ doTwFits(iplane, iside, ipaddle);
+ } // Paddle loop
+ // Produce the time-walk fit parameter graphs
+ for (UInt_t ipar = 0; ipar < nTwFitPars; ipar++) {
+ // Populate graphs and multi-graphs
+ twFitParGraph[iplane][iside][ipar] = new TGraph(nbars[iplane], paddleIndex[iplane][iside], twFitPar[iplane][iside][ipar]);
+ if (iside == 0) addColorToGraph(22, markerSize, kRed, twFitParGraph[iplane][iside][ipar]);
+ if (iside == 1) addColorToGraph(23, markerSize, kBlue, twFitParGraph[iplane][iside][ipar]);
+ // // Add graphs to multi graph
+ // twFitParMultiGraph[iplane][ipar]->Add(twFitParGraph[iplane][iside][ipar]);
+ } // Parameter loop
+ // Calculate the average of the time-walk fit parameters
+ calcParAvg(iplane, iside);
+ } // Side loop
+ // Draw the time-walk parameter graphs
+ drawParams(iplane);
+ } // Plane loop
+ return;
+} // timeWalkFits()
+
+
+
+
+
+
+
diff --git a/CALIBRATION/shms_hodo_calib/timeWalkHistos.C b/CALIBRATION/shms_hodo_calib/timeWalkHistos.C
new file mode 100644
index 0000000000000000000000000000000000000000..0e0dfefdf6ab32db1c02b3928da0e626d73ac13e
--- /dev/null
+++ b/CALIBRATION/shms_hodo_calib/timeWalkHistos.C
@@ -0,0 +1,429 @@
+// Macro to generate time-walk histograms for the hodoscopes
+// Author: Eric Pooser, pooser@jlab.org
+
+#include <time.h>
+
+// Declare replay data file and output file
+TFile *replayFile, *outFile;
+
+// Declare data tree
+TTree *rawDataTree;
+
+// Declare constants
+static const UInt_t nPlanes = 4;
+static const UInt_t nSides = 2;
+static const UInt_t nAdcSignals = 3;
+static const UInt_t nTdcSignals = 1;
+static const UInt_t nSignals = nAdcSignals + nTdcSignals;
+static const UInt_t maxTdcHits = 128;
+static const UInt_t maxAdcHits = 4;
+static const UInt_t nMaxBars = 21;
+
+static const TString spec = "P";
+static const TString detec = "hod";
+static const TString trig = "T";
+static const TString shms = "shms";
+
+static const UInt_t nbars[nPlanes] = {13, 13, 14, 21};
+static const TString planeNames[nPlanes] = {"1x", "1y", "2x", "2y"};
+static const TString sideNames[nSides] = {"pos", "neg"};
+static const TString signalNames[nSignals] = {"Adc", "Tdc"};
+static const TString adcData[nAdcSignals] = {"ErrorFlag", "PulseTimeRaw", "PulseAmp"};
+static const TString tdcData[nTdcSignals] = {"TimeRaw"};
+
+static const Double_t tdcChanToTime = 0.100; // Units of ns
+static const Double_t adcChanToTime = 0.0625; // Units of ns
+static const Double_t adcDynamicRange = 1000.0; // Units of mV
+static const Double_t nAdcChan = 4096.0; // Units of ADC channels
+static const Double_t adcChanTomV = adcDynamicRange/nAdcChan; // Units of mV/ADC Chan
+
+static const Double_t hodoPulseAmpCutLow = 10.0; // Units of mV
+static const Double_t hodoPulseAmpCutHigh = 1000.0; // Units of mV
+static const Double_t refAdcPulseAmpCutLow = 50.0; // Units of mV
+static const Double_t refAdcPulseAmpCutHigh = 60.0; // Units of mV
+static const Double_t refAdcPulseTimeCutLow = 210.0; // Units of ns
+static const Double_t refAdcPulseTimeCutHigh = 225.0; // Units of ns
+static const Double_t adcTdcTimeDiffCutLow = -100.0; // Units of ns
+static const Double_t adcTdcTimeDiffCutHigh = 100.0; // Units of ns
+static const Double_t calEtotCutVal = 0.100; // Units of GeV
+static const Double_t cerNpeSumCutVal = 1.5; // Units of NPE
+// static const Double_t adcTdcTimeDiffCutLow = -6000.0; // Units of ns
+// static const Double_t adcTdcTimeDiffCutHigh = 1000.0; // Units of ns
+
+static const Double_t fontSize = 0.04;
+
+// Declare variables to obtain from data tree
+// Number of ADC & TDC hits
+Int_t adcHits[nPlanes][nSides][nSignals], tdcHits[nPlanes][nSides][nSignals];
+// Paddle number which fired
+Double_t adcPaddle[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+Double_t tdcPaddle[nPlanes][nSides][nSignals][nMaxBars*maxTdcHits];
+// FADC data
+Double_t hodoAdcPulseTimeRaw[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+Double_t hodoAdcPulseAmp[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+Double_t hodoAdcErrorFlag[nPlanes][nSides][nSignals][nMaxBars*maxAdcHits];
+// TDC Data
+Double_t hodoTdcTimeRaw[nPlanes][nSides][nSignals][nMaxBars*maxTdcHits];
+// Trigger apparatus data
+Double_t refAdcPulseTimeRaw, refAdcPulseAmp, refAdcMultiplicity;
+Double_t refT1TdcTimeRaw, refT2TdcTimeRaw, refT3TdcTimeRaw;
+// Declare variables to steer data to arrays
+TString *adcBaseName, *adcNdataName, *adcPaddleName;
+TString *adcErrorFlagName, *adcPulseTimeRawName, *adcPulseAmpName;
+TString *tdcBaseName, *tdcNdataName, *tdcPaddleName;
+TString *tdcTimeRawName;
+// Declare directories for output ROOT file
+TDirectory *trigRawDir, *hodoRawDir, *hodoUncalibDir;
+TDirectory *planeRawDir[nPlanes], *sideRawDir[nPlanes][nSides], *planeUncalibDir[nPlanes], *sideUncalibDir[nPlanes][nSides];
+TDirectory *adcTimeWalkDir[nPlanes][nSides], *tdcTimeWalkDir[nPlanes][nSides], *adcTdcTimeDiffWalkDir[nPlanes][nSides];
+// Declare histos
+// 1D histos
+TH1F *h1_refAdcPulseTimeRaw, *h1_refAdcPulseAmp, *h1_refAdcMultiplicity;
+TH1F *h1_refT1TdcTimeRaw, *h1_refT2TdcTimeRaw, *h1_refT3TdcTimeRaw;
+// 2D histos
+TH2F *h2_adcErrorFlag[nPlanes][nSides];
+TH2F *h2_adcPulseTimeRaw[nPlanes][nSides];
+TH2F *h2_adcPulseAmp[nPlanes][nSides];
+TH2F *h2_tdcTimeRaw[nPlanes][nSides];
+
+TH2F *h2_adcPulseTime[nPlanes][nSides];
+TH2F *h2_tdcTime[nPlanes][nSides];
+TH2F *h2_adcTdcTimeDiff[nPlanes][nSides];
+TH2F *h2_adcPulseAmpCuts[nPlanes][nSides];
+TH2F *h2_adcTimeWalk[nPlanes][nSides][nMaxBars];
+TH2F *h2_tdcTimeWalk[nPlanes][nSides][nMaxBars];
+TH2F *h2_adcTdcTimeDiffWalk[nPlanes][nSides][nMaxBars];
+
+// Declare event data objects
+UInt_t nentries, ievent, adcPaddleNum, tdcPaddleNum;
+
+Int_t numAdcHits, numTdcHits;
+
+Double_t adcErrorFlag, adcPulseTimeRaw, adcPulseTime, adcPulseAmp;
+Double_t tdcTimeRaw, tdcTime, adcTdcTimeDiff;
+Double_t calEtot, cerNpeSum;
+
+Bool_t adcRefMultiplicityCut, adcRefPulseAmpCut, adcRefPulseTimeCut;
+Bool_t edtmCut, adcErrorFlagCut, adcAndTdcHitCut;
+Bool_t adcPulseAmpCut, adcTdcTimeDiffCut;
+Bool_t calEtotCut, cerNpeSumCut;
+
+void generatePlots(UInt_t iplane, UInt_t iside, UInt_t ipaddle) {
+
+ // Create trigger aparatus directory
+ trigRawDir = dynamic_cast <TDirectory*> (outFile->Get("trigAppRaw"));
+ if (!trigRawDir) {trigRawDir = outFile->mkdir("trigAppRaw"); trigRawDir->cd();}
+ else outFile->cd("trigAppRaw");
+ // FADC reference
+ if (!h1_refAdcPulseTimeRaw) h1_refAdcPulseTimeRaw = new TH1F("h1_refAdcPulseTimeRaw", "ROC2 Raw FADC Reference Pulse Time; Raw FADC Pulse Time (ns); Number of Entries / 100 ps", 4000, 0, 400);
+ if (!h1_refAdcPulseAmp) h1_refAdcPulseAmp = new TH1F("h1_refAdcPulseAmp", "ROC2 FADC Reference Pulse Amplitude; FADC Pulse Amplitude (mV); Number of Entries / 1 mV", 1000, 0, 1000);
+ if (!h1_refAdcMultiplicity) h1_refAdcMultiplicity = new TH1F("h1_refAdcMultiplicity", "ROC2 FADC Reference Multiplicity; Raw FADC Multiplicity; Number of Entries", 5, -0.5, 4.5);
+ // TDC reference
+ if (!h1_refT1TdcTimeRaw) h1_refT1TdcTimeRaw = new TH1F("h1_refT1TdcTimeRaw", "ROC2 T1 Raw TDC Reference TDC Time (Slot 20, Channel 15); Raw TDC Time (ns); Number of Entries / 100 ps", 6000, 0, 600);
+ if (!h1_refT2TdcTimeRaw) h1_refT2TdcTimeRaw = new TH1F("h1_refT2TdcTimeRaw", "ROC2 T2 Raw TDC Reference TDC Time (Slot 19, Channel 31); Raw TDC Time (ns); Number of Entries / 100 ps", 6000, 0, 600);
+ if (!h1_refT3TdcTimeRaw) h1_refT3TdcTimeRaw = new TH1F("h1_refT3TdcTimeRaw", "ROC2 T3 Raw TDC Reference TDC Time (Slot 19, Channel 38); Raw TDC Time (ns); Number of Entries / 100 ps", 6000, 0, 600);
+
+ // Create plane directory for raw hodoscope quantities
+ hodoRawDir = dynamic_cast <TDirectory*> (outFile->Get("hodoRaw"));
+ if (!hodoRawDir) {hodoRawDir = outFile->mkdir("hodoRaw"); hodoRawDir->cd();}
+ else outFile->cd("hodoRaw");
+ planeRawDir[iplane] = dynamic_cast <TDirectory*> (hodoRawDir->Get(planeNames[iplane]));
+ if (!planeRawDir[iplane]) {planeRawDir[iplane] = hodoRawDir->mkdir(planeNames[iplane]); planeRawDir[iplane]->cd();}
+ else hodoRawDir->cd(planeNames[iplane]);
+ // Create side directory
+ sideRawDir[iplane][iside] = dynamic_cast <TDirectory*> (planeRawDir[iplane]->Get(sideNames[iside]));
+ if (!sideRawDir[iplane][iside]) {sideRawDir[iplane][iside] = planeRawDir[iplane]->mkdir(sideNames[iside]); sideRawDir[iplane][iside]->cd();}
+ else outFile->cd("hodoRaw/"+planeNames[iplane]+"/"+sideNames[iside]);
+ // Book histos
+ if (!h2_adcErrorFlag[iplane][iside]) h2_adcErrorFlag[iplane][iside] = new TH2F("h2_adcErrorFlag", "FADC Error Flag Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Error Flag", nbars[iplane], 0.5, nbars[iplane]+0.5, 2, -0.5, 1.5);
+ if (!h2_adcPulseTimeRaw[iplane][iside]) h2_adcPulseTimeRaw[iplane][iside] = new TH2F("h2_adcPulseTimeRaw", "Raw FADC Pulse Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; Raw FADC Pulse Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 4000, 0, 400);
+ if (!h2_adcPulseAmp[iplane][iside]) h2_adcPulseAmp[iplane][iside] = new TH2F("h2_adcPulseAmp", "FADC Pulse Amplitude Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Pulse Amplitude (mV) / 1 mV", nbars[iplane], 0.5, nbars[iplane]+0.5, 1000, 0, 1000);
+ if (!h2_tdcTimeRaw[iplane][iside]) h2_tdcTimeRaw[iplane][iside] = new TH2F("h2_tdcTimeRaw", "Raw TDC Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; Raw TDC Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 4000, 0, 400);
+
+ // Create plane directory for uncalibrated hodoscope quantities
+ hodoUncalibDir = dynamic_cast <TDirectory*> (outFile->Get("hodoUncalib"));
+ if (!hodoUncalibDir) {hodoUncalibDir = outFile->mkdir("hodoUncalib"); hodoUncalibDir->cd();}
+ else outFile->cd("hodoUncalib");
+ planeUncalibDir[iplane] = dynamic_cast <TDirectory*> (hodoUncalibDir->Get(planeNames[iplane]));
+ if (!planeUncalibDir[iplane]) {planeUncalibDir[iplane] = hodoUncalibDir->mkdir(planeNames[iplane]); planeUncalibDir[iplane]->cd();}
+ else hodoUncalibDir->cd(planeNames[iplane]);
+ // Create side directory
+ sideUncalibDir[iplane][iside] = dynamic_cast <TDirectory*> (planeUncalibDir[iplane]->Get(sideNames[iside]));
+ if (!sideUncalibDir[iplane][iside]) {sideUncalibDir[iplane][iside] = planeUncalibDir[iplane]->mkdir(sideNames[iside]); sideUncalibDir[iplane][iside]->cd();}
+ else outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]);
+ // Book histos
+ if (!h2_adcPulseTime[iplane][iside]) h2_adcPulseTime[iplane][iside] = new TH2F("h2_adcPulseTime", "FADC Pulse Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Pulse Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 8000, -400, 400);
+ if (!h2_tdcTime[iplane][iside]) h2_tdcTime[iplane][iside] = new TH2F("h2_tdcTime", "TDC Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; TDC Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 8000, -400, 400);
+ if (!h2_adcTdcTimeDiff[iplane][iside]) h2_adcTdcTimeDiff[iplane][iside] = new TH2F("h2_adcTdcTimeDiff", "TDC-ADC Time Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; TDC - ADC Time (ns) / 100 ps", nbars[iplane], 0.5, nbars[iplane]+0.5, 16000, -800, 800);
+ if (!h2_adcPulseAmpCuts[iplane][iside]) h2_adcPulseAmpCuts[iplane][iside] = new TH2F("h2_adcPulseAmpCuts", "FADC Pulse Amplitude Cuts Plane "+planeNames[iplane]+" Side "+sideNames[iside]+"; PMT Number; FADC Pulse Amplitude (mV) / 1 mV", nbars[iplane], 0.5, nbars[iplane]+0.5, 1000, 0, 1000);
+ // Create ADC time-walk directory
+ adcTimeWalkDir[iplane][iside] = dynamic_cast <TDirectory*> (sideUncalibDir[iplane][iside]->Get("adcTimeWalk"));
+ if (!adcTimeWalkDir[iplane][iside]) {adcTimeWalkDir[iplane][iside] = sideUncalibDir[iplane][iside]->mkdir("adcTimeWalk"); adcTimeWalkDir[iplane][iside]->cd();}
+ else (outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]+"/adcTimeWalk"));
+ // Book histos
+ if (!h2_adcTimeWalk[iplane][iside][ipaddle]) h2_adcTimeWalk[iplane][iside][ipaddle] = new TH2F(Form("h2_adcTimeWalk_paddle_%d", ipaddle+1), "ADC Time vs. Pulse Amp Plane "+planeNames[iplane]+" Side "+sideNames[iside]+Form(" Paddle %d", ipaddle+1)+"; Pulse Amplitude (mV) / 1 mV; ADC Time (ns) / 100 ps", 1000, 0, 1000, 4000, -200, 200);
+ // Create TDC time-walk directory
+ tdcTimeWalkDir[iplane][iside] = dynamic_cast <TDirectory*> (sideUncalibDir[iplane][iside]->Get("tdcTimeWalk"));
+ if (!tdcTimeWalkDir[iplane][iside]) {tdcTimeWalkDir[iplane][iside] = sideUncalibDir[iplane][iside]->mkdir("tdcTimeWalk"); tdcTimeWalkDir[iplane][iside]->cd();}
+ else (outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]+"/tdcTimeWalk"));
+ // Book histos
+ if (!h2_tdcTimeWalk[iplane][iside][ipaddle]) h2_tdcTimeWalk[iplane][iside][ipaddle] = new TH2F(Form("h2_tdcTimeWalk_paddle_%d", ipaddle+1), "TDC Time vs. Pulse Amp Plane "+planeNames[iplane]+" Side "+sideNames[iside]+Form(" Paddle %d", ipaddle+1)+"; Pulse Amplitude (mV) / 1 mV; TDC Time (ns) / 100 ps", 1000, 0, 1000, 4000, -300, 100);
+ // Create TDC-ADC time-walk directory
+ adcTdcTimeDiffWalkDir[iplane][iside] = dynamic_cast <TDirectory*> (sideUncalibDir[iplane][iside]->Get("adcTdcTimeDiffWalk"));
+ if (!adcTdcTimeDiffWalkDir[iplane][iside]) {adcTdcTimeDiffWalkDir[iplane][iside] = sideUncalibDir[iplane][iside]->mkdir("adcTdcTimeDiffWalk"); adcTdcTimeDiffWalkDir[iplane][iside]->cd();}
+ else (outFile->cd("hodoUncalib/"+planeNames[iplane]+"/"+sideNames[iside]+"/adcTdcTimeDiffWalk"));
+ // Book histos
+ if (!h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle]) h2_adcTdcTimeDiffWalk[iplane][iside][ipaddle] = new TH2F(Form("h2_adcTdcTimeDiffWalk_paddle_%d", ipaddle+1), "TDC-ADC Time vs. Pulse Amp Plane "+planeNames[iplane]+" Side "+sideNames[iside]+Form(" Paddle %d", ipaddle+1)+"; Pulse Amplitude (mV) / 1 mV; TDC-ADC Time (ns) / 100 ps", 1000, 0, 1000, 500, -20, 30);
+
+} // generatePlots()
+
+void timeWalkHistos(Int_t runNum) {
+
+ // Global ROOT settings
+ gStyle->SetOptFit();
+ gStyle->SetStatFormat(".2f");
+ gStyle->SetTitleFontSize(fontSize);
+ gStyle->SetLabelSize(fontSize, "XY");
+ gStyle->SetTitleSize(fontSize, "XY");
+
+ // Initialize the analysis clock
+ clock_t t;
+ t = clock();
+
+ // Obtain the replay data file and create new output ROOT file
+
+ replayFile = new TFile(Form("ROOTfiles/shms_replay_production_all_%d_-1.root", runNum), "READ");
+ // replayFile = new TFile(Form("ROOTfiles/shms_coin_replay_production_%d_-1.root", runNum), "READ");
+
+ outFile = new TFile("timeWalkHistos_temp.root", "RECREATE");
+ // Obtain the tree
+ rawDataTree = dynamic_cast <TTree*> (replayFile->Get("T"));
+ // Acquire the trigger apparatus data
+ rawDataTree->SetBranchAddress("T.shms.pFADC_TREF_ROC2_adcPulseTimeRaw", &refAdcPulseTimeRaw);
+ rawDataTree->SetBranchAddress("T.shms.pFADC_TREF_ROC2_adcPulseAmp", &refAdcPulseAmp);
+ rawDataTree->SetBranchAddress("T.shms.pFADC_TREF_ROC2_adcMultiplicity", &refAdcMultiplicity);
+ rawDataTree->SetBranchAddress("T.shms.pT1_tdcTimeRaw", &refT1TdcTimeRaw);
+ rawDataTree->SetBranchAddress("T.shms.pT2_tdcTimeRaw", &refT2TdcTimeRaw);
+ rawDataTree->SetBranchAddress("T.shms.pT3_tdcTimeRaw", &refT3TdcTimeRaw);
+ rawDataTree->SetBranchAddress("P.cal.etot", &calEtot);
+ rawDataTree->SetBranchAddress("P.ngcer.npeSum", &cerNpeSum);
+ // Loop over the planes, sides, signals, leafs, and fill data arrays
+ for(UInt_t iplane = 0; iplane < nPlanes; iplane++) {
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ // Generate directory structure and histograms
+ for(UInt_t ipaddle = 0; ipaddle < nbars[iplane]; ipaddle++)
+ generatePlots(iplane, iside, ipaddle);
+ for(UInt_t isignal = 0; isignal < nSignals; isignal++) {
+
+ // Acquire the hodoscope ADC data objects
+ if(signalNames[isignal] == "Adc") {
+ // Define the leaf strings of interest
+ adcBaseName = new TString(spec+"."+detec+"."+planeNames[iplane]+"."+
+ sideNames[iside]+signalNames[isignal]);
+ adcNdataName = new TString("Ndata."+*adcBaseName+"Counter");
+ adcPaddleName = new TString(*adcBaseName+"Counter");
+ // Fill the leafs
+ rawDataTree->SetBranchAddress(*adcNdataName, &adcHits[iplane][iside][isignal]);
+ rawDataTree->SetBranchAddress(*adcPaddleName, &adcPaddle[iplane][iside][isignal]);
+ // Loop over the adc data objects
+ for(UInt_t iadcdata = 0; iadcdata < nAdcSignals; iadcdata++) {
+ if (adcData[iadcdata] == "ErrorFlag") {
+ adcErrorFlagName = new TString(*adcBaseName+"ErrorFlag");
+ rawDataTree->SetBranchAddress(*adcErrorFlagName, &hodoAdcErrorFlag[iplane][iside][isignal]);
+ } // Error flag leaf
+ if (adcData[iadcdata] == "PulseTimeRaw") {
+ adcPulseTimeRawName = new TString(*adcBaseName+"PulseTimeRaw");
+ rawDataTree->SetBranchAddress(*adcPulseTimeRawName, &hodoAdcPulseTimeRaw[iplane][iside][isignal]);
+ } // Raw pulse time leaf
+ if (adcData[iadcdata] == "PulseAmp") {
+ adcPulseAmpName = new TString(*adcBaseName+"PulseAmp");
+ rawDataTree->SetBranchAddress(*adcPulseAmpName, &hodoAdcPulseAmp[iplane][iside][isignal]);
+ } // Pedestal subtracted pulse amplitude leaf
+ } // ADC signal data loop
+ } // ADC signal
+
+ // Acquire the hodoscope TDC data objects
+ if(signalNames[isignal] == "Tdc") {
+ // Define the leaf strings of interest
+ tdcBaseName = new TString(spec+"."+detec+"."+planeNames[iplane]+"."+
+ sideNames[iside]+signalNames[isignal]);
+ tdcNdataName = new TString("Ndata."+*tdcBaseName+"Counter");
+ tdcPaddleName = new TString(*tdcBaseName+"Counter");
+ // Fill the leafs
+ rawDataTree->SetBranchAddress(*tdcNdataName, &tdcHits[iplane][iside][isignal]);
+ rawDataTree->SetBranchAddress(*tdcPaddleName, &tdcPaddle[iplane][iside][isignal]);
+ // Loop over the TDC data objects
+ for(UInt_t itdcdata = 0; itdcdata < nTdcSignals; itdcdata++) {
+ if (tdcData[itdcdata] == "TimeRaw") {
+ tdcTimeRawName = new TString(*tdcBaseName+"TimeRaw");
+ rawDataTree->SetBranchAddress(*tdcTimeRawName, &hodoTdcTimeRaw[iplane][iside][isignal]);
+ } // Raw TDC time leaf
+ } // TDC signal data loop
+ } // TDC signal
+
+ } // Signal loop
+ } // Side loop
+ } // Plane loop
+
+ // Loop over the events and fill histograms
+ // nentries = rawDataTree->GetEntries();
+ nentries = 100000;
+ cout << "\n******************************************" << endl;
+ cout << nentries << " Events Will Be Processed" << endl;
+ cout << "******************************************\n" << endl;
+ for (ievent = 0; ievent < nentries; ievent++) {
+ // Acquire the event from the data tree
+ rawDataTree->GetEntry(ievent);
+ // Fiducial PID cuts
+ calEtotCut = (calEtot < calEtotCutVal);
+ cerNpeSumCut = (cerNpeSum < cerNpeSumCutVal);
+ if (calEtotCut || cerNpeSumCut) continue;
+ // Fill trigger apparatus histos
+ h1_refAdcPulseTimeRaw->Fill(refAdcPulseTimeRaw*adcChanToTime);
+ h1_refAdcPulseAmp->Fill(refAdcPulseAmp);
+ h1_refAdcMultiplicity->Fill(refAdcMultiplicity);
+ h1_refT1TdcTimeRaw->Fill(refT1TdcTimeRaw*tdcChanToTime);
+ h1_refT2TdcTimeRaw->Fill(refT2TdcTimeRaw*tdcChanToTime);
+ h1_refT3TdcTimeRaw->Fill(refT3TdcTimeRaw*tdcChanToTime);
+ // Loop over the planes, sides, signals, data arrays, and fill hodoscope histograms
+ for(UInt_t iplane = 0; iplane < nPlanes; iplane++) {
+ for(UInt_t iside = 0; iside < nSides; iside++) {
+ for(UInt_t isignal = 0; isignal < nSignals; isignal++) {
+
+ // Acquire the hodoscope ADC data objects
+ if(signalNames[isignal] == "Adc") {
+ // Loop over the ADC data objects
+ for(UInt_t iadcdata = 0; iadcdata < nAdcSignals; iadcdata++) {
+ numAdcHits = adcHits[iplane][iside][isignal];
+ for (Int_t iadchit = 0; iadchit < numAdcHits; iadchit++) {
+ // Obtain variables
+ adcPaddleNum = adcPaddle[iplane][iside][isignal][iadchit];
+ adcErrorFlag = hodoAdcErrorFlag[iplane][iside][isignal][iadchit];
+ adcPulseTimeRaw = hodoAdcPulseTimeRaw[iplane][iside][isignal][iadchit]*adcChanToTime;
+ adcPulseTime = adcPulseTimeRaw - refAdcPulseTimeRaw*adcChanToTime;
+ adcPulseAmp = hodoAdcPulseAmp[iplane][iside][isignal][iadchit];
+ // Fill histos
+ if (adcData[iadcdata] == "ErrorFlag")
+ h2_adcErrorFlag[iplane][iside]->Fill(adcPaddleNum, adcErrorFlag);
+ if (adcData[iadcdata] == "PulseTimeRaw") {
+ h2_adcPulseTimeRaw[iplane][iside]->Fill(adcPaddleNum, adcPulseTimeRaw);
+ h2_adcPulseTime[iplane][iside]->Fill(adcPaddleNum, adcPulseTime);
+ }
+ if (adcData[iadcdata] == "PulseAmp")
+ h2_adcPulseAmp[iplane][iside]->Fill(adcPaddleNum, adcPulseAmp);
+ } // ADC hit loop
+ } // ADC signal data loop
+ } // ADC signal
+
+ // Acquire the hodoscope TDC data objects
+ if(signalNames[isignal] == "Tdc") {
+ // Loop over the TDC data objects
+ for(UInt_t itdcdata = 0; itdcdata < nTdcSignals; itdcdata++) {
+ numTdcHits = tdcHits[iplane][iside][isignal];
+ // Define cuts
+ //edtmCut = (numTdcHits == nbars[0] || numTdcHits == nbars[2] || numTdcHits == nbars[3]);
+ // Implement cuts
+ //if (edtmCut) continue;
+ for (Int_t itdchit = 0; itdchit < numTdcHits; itdchit++) {
+ // Obtain variables
+ tdcPaddleNum = tdcPaddle[iplane][iside][isignal][itdchit];
+ tdcTimeRaw = hodoTdcTimeRaw[iplane][iside][isignal][itdchit]*tdcChanToTime;
+ if (iplane == 3 && iside == 0) tdcTime = tdcTimeRaw - refT2TdcTimeRaw*tdcChanToTime;
+ else tdcTime = tdcTimeRaw - refT1TdcTimeRaw*tdcChanToTime;
+ if (tdcData[itdcdata] == "TimeRaw") {
+ h2_tdcTimeRaw[iplane][iside]->Fill(tdcPaddleNum, tdcTimeRaw);
+ h2_tdcTime[iplane][iside]->Fill(tdcPaddleNum, tdcTime);
+ } // RDC raw time signal
+ } // TDC hit loop
+ } // TDC signal data loop
+ } // TDC signal
+
+ // Define cuts
+ adcRefMultiplicityCut = (refAdcMultiplicity != 1.0);
+ adcRefPulseAmpCut = (refAdcPulseAmp < refAdcPulseAmpCutLow || refAdcPulseAmp > refAdcPulseAmpCutHigh);
+ adcRefPulseTimeCut = (refAdcPulseTimeRaw*adcChanToTime < refAdcPulseTimeCutLow || refAdcPulseTimeRaw*adcChanToTime > refAdcPulseTimeCutHigh);
+ // Implement cuts
+ //if (adcRefMultiplicityCut || adcRefPulseAmpCut || adcRefPulseTimeCut) continue;
+ // Acquire the hodoscope ADC data objects
+ if(signalNames[isignal] == "Adc") {
+ // Loop over the signals again
+ for(UInt_t jsignal = 0; jsignal < nSignals; jsignal++) {
+ // Acquire othe hodoscope TDC data objects
+ if(signalNames[jsignal] == "Tdc") {
+ // Loop over the ADC data objects
+ for(UInt_t iadcdata = 0; iadcdata < nAdcSignals; iadcdata++) {
+ if (adcData[iadcdata] != "PulseTimeRaw") continue;
+ // Acquire the number of ADC hits
+ numAdcHits = adcHits[iplane][iside][isignal];
+ for (Int_t iadchit = 0; iadchit < numAdcHits; iadchit++) {
+ // Obtain variables
+ adcPaddleNum = UInt_t (adcPaddle[iplane][iside][isignal][iadchit]);
+ adcErrorFlag = hodoAdcErrorFlag[iplane][iside][isignal][iadchit];
+ adcPulseTimeRaw = hodoAdcPulseTimeRaw[iplane][iside][isignal][iadchit]*adcChanToTime;
+ adcPulseTime = adcPulseTimeRaw - refAdcPulseTimeRaw*adcChanToTime;
+ adcPulseAmp = hodoAdcPulseAmp[iplane][iside][isignal][iadchit];
+ // Define cuts
+ adcErrorFlagCut = (adcErrorFlag != 0);
+ adcPulseAmpCut = (adcPulseAmp < hodoPulseAmpCutLow || adcPulseAmp > hodoPulseAmpCutHigh);
+ // Implement cuts
+ if (adcErrorFlagCut || adcPulseAmpCut) continue;
+ // Loop over the TDC data objects
+ for(UInt_t itdcdata = 0; itdcdata < nTdcSignals; itdcdata++) {
+ if (tdcData[itdcdata] != "TimeRaw") continue;
+ // Acquire the number of TDC hits
+ numTdcHits = tdcHits[iplane][iside][jsignal];
+ // Define cuts
+ //edtmCut = (numTdcHits == nbars[0] || numTdcHits == nbars[2] || numTdcHits == nbars[3]);
+ // Implement cuts
+ //if (edtmCut) continue;
+ for (Int_t itdchit = 0; itdchit < numTdcHits; itdchit++) {
+ // Obtain variables
+ tdcPaddleNum = UInt_t (tdcPaddle[iplane][iside][jsignal][itdchit]);
+ tdcTimeRaw = hodoTdcTimeRaw[iplane][iside][jsignal][itdchit]*tdcChanToTime;
+ if (iplane == 3 && iside == 0) tdcTime = tdcTimeRaw - refT2TdcTimeRaw*tdcChanToTime;
+ else tdcTime = tdcTimeRaw - refT1TdcTimeRaw*tdcChanToTime;
+ adcTdcTimeDiff = tdcTime - adcPulseTime;
+ // Define cuts
+ adcAndTdcHitCut = (adcPaddleNum != tdcPaddleNum);
+ adcTdcTimeDiffCut = (adcTdcTimeDiff < adcTdcTimeDiffCutLow || adcTdcTimeDiff > adcTdcTimeDiffCutHigh);
+ // Implement cuts
+ if (adcAndTdcHitCut || adcTdcTimeDiffCut) continue;
+ h2_adcPulseAmpCuts[iplane][iside]->Fill(tdcPaddleNum, adcPulseAmp);
+ h2_adcTdcTimeDiff[iplane][iside]->Fill(tdcPaddleNum, adcTdcTimeDiff);
+ h2_adcTimeWalk[iplane][iside][tdcPaddleNum-1]->Fill(adcPulseAmp, adcPulseTime);
+ h2_tdcTimeWalk[iplane][iside][tdcPaddleNum-1]->Fill(adcPulseAmp, tdcTime);
+ h2_adcTdcTimeDiffWalk[iplane][iside][tdcPaddleNum-1]->Fill(adcPulseAmp, adcTdcTimeDiff);
+ h2_adcTdcTimeDiffWalk[iplane][iside][tdcPaddleNum-1]->GetXaxis()->CenterTitle();
+ h2_adcTdcTimeDiffWalk[iplane][iside][tdcPaddleNum-1]->GetYaxis()->CenterTitle();
+ } // TDC hit loop
+ } // TDC signal loop
+ } // ADC hit loop
+ } // ADC signal loop
+ } // TDC signal
+ } // Nested signal loop
+ } // ADC signal
+
+ } // Signal loop
+ } // Side loop
+ } // Plane loop
+
+ if (ievent % 100000 == 0 && ievent != 0)
+ cout << ievent << " Events Have Been Processed..." << endl;
+
+ } // rawDataTree event loop
+
+ cout << "\n***************************************" << endl;
+ cout << ievent << " Events Were Processed" << endl;
+ cout << "***************************************\n" << endl;
+
+ // Calculate the analysis rate
+ t = clock() - t;
+ printf ("The Analysis Took %.1f seconds \n", ((float) t) / CLOCKS_PER_SEC);
+ printf ("The Analysis Event Rate Was %.3f kHz \n", (ievent + 1) / (((float) t) / CLOCKS_PER_SEC*1000.));
+
+ outFile->Write();
+ //outFile->Close();
+
+ return 0;
+
+} // time_walk_calib()