From 027a58ce18e7e350948329c6930b191c0344071f Mon Sep 17 00:00:00 2001
From: Yero1990 <cyero002@fiu.edu>
Date: Sun, 11 Jun 2017 10:49:41 -0400
Subject: [PATCH] Dc work (#134)
* -added a few lines of code to produce a weighted tzero data file
per wire basis, where the weight is by groups (ribbon cables)
* -started to add lines to read per-wire tzero parameter data file. (work in progress)
* -mofidied plane numbering in 'get_zero_param.C' to match the rest of the
calibration script plane numbering scheme
* -minor fix
* -script that reads a data file containing the per wire tzeros
-the tzero correction is applies on a wire-by-wire basis and the
wire drift times are added to produce the corrected plane drift times
* minor modifications to name of file being read
* -minor modifications to file name being read
---
.../shms_dc_calib/scripts/get_LookUp_Values.C | 8 +-
.../scripts/get_pdc_time_histo.C | 12 +-
.../get_pdc_time_histo_tzero_corrected.C | 246 ++++++++++++------
.../shms_dc_calib/scripts/get_tzero_param.C | 35 ++-
4 files changed, 204 insertions(+), 97 deletions(-)
diff --git a/CALIBRATION/shms_dc_calib/scripts/get_LookUp_Values.C b/CALIBRATION/shms_dc_calib/scripts/get_LookUp_Values.C
index 8e0515e8..aee70556 100644
--- a/CALIBRATION/shms_dc_calib/scripts/get_LookUp_Values.C
+++ b/CALIBRATION/shms_dc_calib/scripts/get_LookUp_Values.C
@@ -15,7 +15,7 @@ void get_LookUp_Values() {
infile >> run_NUM;
//Open root file containing drift time histos
- TFile *f = new TFile(Form("../root_files/run%d/shms_dc_t0_corrected_%d.root", run_NUM, run_NUM),"READ");
+ TFile *f = new TFile(Form("../root_files/run%d/shms_tzero_corr_histos.root", run_NUM),"READ");
//Define histogram array
TH1F *h[NPLANES];
@@ -41,9 +41,9 @@ void get_LookUp_Values() {
//Set headers for subsequent columns of data
ofs << Form("; Lookup Table: RUN %d", run_NUM) << "\n";
- ofs << "; number of bins in Carlos's time to distance lookup table" << "\n";
+ ofs << "; number of bins in time to distance lookup table" << "\n";
ofs << Form("pdriftbins = %d", TOTAL_BINS+1) << "\n";
- ofs << "; number of 1st bin in Carlos's table in ns" << "\n";
+ ofs << "; number of 1st bin in table in ns" << "\n";
ofs << "pdrift1stbin=0" << "\n";
ofs << "; bin size in ns" << "\n";
ofs << "pdriftbinsz=2" << "\n";
@@ -54,7 +54,7 @@ void get_LookUp_Values() {
for (int ip=0; ip<NPLANES; ip++){
- TString drift_time_histo = "pdc"+plane_names[ip]+"_time: t0_corr";
+ TString drift_time_histo = "all_wires_"+plane_names[ip];
//Get drift time histograms from root file
h[ip] = (TH1F*)f->Get(drift_time_histo);
diff --git a/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo.C b/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo.C
index 2c5e4984..9c5cc4fa 100644
--- a/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo.C
+++ b/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo.C
@@ -58,7 +58,7 @@ void get_pdc_time_histo()
TString title = "pdc"+plane_names[ip]+"_drifttime";
//Set Branch Address
- tree->SetBranchAddress(drift_time, &pdc_time[ip][0]);
+ tree->SetBranchAddress(drift_time, pdc_time[ip]);
tree->SetBranchAddress(ndata_name, &Ndata[ip]); /* Ndata represents number of triggers vs number of hits that each trigger produced.
A hit is refer to as when a trigger(traversing particle), ionizes the WC gas and ionized
electrons reach the rearest sense wire, producing a detectable signal in the O'scope */
@@ -67,7 +67,7 @@ void get_pdc_time_histo()
h[ip] = new TH1F(drift_time_histo, title, 200, -50, 350); //set time to 400 ns/200 bins = 2ns/bin
}
-
+
//Declare number of entries in the tree
Long64_t nentries = tree->GetEntries(); //number of triggers (particles that passed through all 4 hodo planes)
@@ -75,10 +75,10 @@ void get_pdc_time_histo()
for(Long64_t i=0; i<nentries; i++)
{
tree->GetEntry(i);
-
+ // cout << i << endl;
//Loop over number of hits for each trigger in each DC plane
- for(ip=0; ip<NPLANES; ip++){
+ for(Int_t ip=0; ip<NPLANES; ip++){
for(Int_t j=0; j<Ndata[ip]; j++){
@@ -89,7 +89,9 @@ void get_pdc_time_histo()
}
}
-
+
//Write histograms to file
+
g->Write();
+
}
diff --git a/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo_tzero_corrected.C b/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo_tzero_corrected.C
index fa028104..07298262 100644
--- a/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo_tzero_corrected.C
+++ b/CALIBRATION/shms_dc_calib/scripts/get_pdc_time_histo_tzero_corrected.C
@@ -1,12 +1,92 @@
-//Script to add t0 correction to SHMS DC drift times
-
#define NPLANES 12
+using namespace std;
void get_pdc_time_histo_tzero_corrected()
{
+
+
+ //this script reads all tzero values inside tzero_group.dat and assigns them to their corresponding planes.
+ //these values will then be used to shift the drift time histos, wire by wire (t0-correction wire-by-wire)
+
+ int ip; //to loop over planes
+ int sw; //to loop over sensewires
+
+
+ int tot_wires[NPLANES] = {107, 107, 79, 79, 107, 107, 107, 107, 79, 79, 107, 107};
+
+ //sum over all wires in both DC
+ static int wire_sum = 0;
+ for (ip=0; ip<NPLANES; ip++) {
+ wire_sum = wire_sum + tot_wires[ip];
+ }
+
+
+ //open and read tzero data file
+ ifstream file;
+ file.open("../data_files/run484/tzero_group.dat");
+
+ string line;
+ int counter;
+ Double_t value;
+ Double_t tzero_offset[wire_sum];
+
-
- //read run number from input file
+ counter = 0;
+
+ while (getline(file, line)) {
+
+ if (line[0]!='#' )
+ {
+ sscanf(line.c_str(), "%lf", &value); //pass data in file to variable 'value'
+ tzero_offset[counter] = value; //write all tzero values for all wires in both DC to an array
+ // cout << tzero_offset[counter] << endl;
+ counter++;
+ }
+ }
+
+
+ //***************************************************************
+
+ //pass all tzero values into their corresponding planes
+
+ //declare a 2d dynamic array to store tzero values to be used for t0-correction
+ Double_t **tzero = new Double_t*[NPLANES];
+
+ for (ip = 0; ip < NPLANES; ip++)
+ {
+ tzero[ip] = new Double_t[tot_wires[ip]];
+ }
+
+ //initialize 2d dynamic array to 0.0
+ for (ip = 0; ip<NPLANES; ip++) {
+ for (sw = 0; sw<tot_wires[ip]; sw++) {
+ tzero[ip][sw] = 0.0;
+
+ }
+ }
+
+ counter = 0;
+ for (ip = 0; ip<NPLANES; ip++) {
+
+ for (sw = 0; sw < tot_wires[ip]; sw++) {
+
+ tzero[ip][sw] = tzero_offset[counter]; //tzero corrections that must be added wire by wire
+ //cout << tzero[ip][sw] << endl;
+ //cout << tzero[ip][sw] << " :: "<<tzero_offset[counter] << endl;
+ counter++;
+ }
+ }
+
+
+
+
+
+ //*****************************************************************************************
+
+ //THIS SCRIPT WILL READ the drift times array tzero[ip][sw] on a wire basis and apply the tzero correction.
+
+
+
int run_NUM;
TString f0 = "input_RUN.txt";
ifstream infile(f0);
@@ -14,94 +94,112 @@ void get_pdc_time_histo_tzero_corrected()
TString run = Form("run%d", run_NUM);
+ int i;
+
+ int nbins = 200; // set number of bins in histos
+ int bin_width = 2;
+ int bin_Content;
+ Double_t shift; //will be the t0 offset
+
+ //Declare plane names to loop over
+
+ TString plane_names[NPLANES]={"1u1", "1u2", "1x1", "1x2", "1v1", "1v2", "2v2", "2v1", "2x2", "2x1", "2u2", "2u1"};
+
+ //Declare a root file to read individual DC cell drift times
+ TString root_file;
+ TFile *f[NPLANES];
+
+ //Declare root file were NEW tzero corrected histograms will be added
+ TString file_name = Form("../root_files/run%d/shms_tzero_corr_histos.root", run_NUM);
+ TFile *g = new TFile(file_name, "RECREATE");
+
- //open file
- TFile *f = new TFile(Form("../../../ROOTfiles/shms_replay_%d.root", run_NUM), "READ");
+ TH1F *h_add[NPLANES]; //t0-corrected plane drift times
- //updates file
- TFile *g = new TFile(Form("../root_files/run%d/shms_dc_t0_corrected_%d.root", run_NUM, run_NUM), "UPDATE"); // create new file to store histo
+
+
- f->cd();
+ //Loop over all planes
+ for (ip = 0; ip < NPLANES; ip++){
- //Get the tree
- TTree *tree = (TTree*)f->Get("T");
+ //READ wire drift time root file per plane
+ root_file = "../root_files/"+run+"/shms_DC_"+plane_names[ip]+Form("_%d.root",run_NUM);
+ f[ip] = new TFile(root_file, "READ");
+
- TString SPECTROMETER="P";
- TString DETECTOR="dc";
- TString plane_names[NPLANES]={"1u1", "1u2", "1x1", "1x2", "1v1", "1v2", "2v2", "2v1", "2x2", "2x1", "2u2", "2u1"};
+ h_add[ip] =new TH1F("all_wires_"+plane_names[ip], "", nbins, -50, 350);
- //Declare Variables to Loop Over
- Int_t Ndata[NPLANES];
- Double_t pdc_time[NPLANES][1000];
-
- //Declare Histogram array to store AVG drift times per plane
- TH1F* h[NPLANES];
-
- g->cd();
-
- //Loop over each plane
- for(Int_t ip=0; ip<NPLANES; ip++){
- TString base_name = SPECTROMETER+"."+DETECTOR+"."+plane_names[ip];
- TString ndata_name = "Ndata."+base_name+".time";
- TString drift_time = base_name+".time";
-
- TString drift_time_histo = "pdc"+plane_names[ip]+"_time: t0_corr";
- TString title = "pdc"+plane_names[ip]+"_drifttime: t0-corrected";
-
- //Set Branch Address
- tree->SetBranchAddress(drift_time, pdc_time[ip]);
- tree->SetBranchAddress(ndata_name, &Ndata[ip]); /* Ndata represents number of triggers vs number of hits that each trigger produced.
- A hit is refer to as when a trigger(traversing particle), ionizes the WC gas and ionized
- electrons reach the rearest sense wire, producing a detectable signal in the O'scope */
-
- //Create Histograms
- h[ip] = new TH1F(drift_time_histo, title, 200, -50, 350); //set time to 400 ns/200 bins = 2ns/bin
-}
-
-
- //open and read tzero data file
- ifstream ifs;
- ifs.open("../data_files/" + run + "/tzero.dat");
-
- double t_zero_offsets[NPLANES];
- for (ip=0; ip < 12; ip++) {
- ifs >> t_zero_offsets[ip]; //add tzero offsets to array
- }
+ TH1F *cell_dt[tot_wires[ip]];
+ TH1F *cell_dt_corr[tot_wires[ip]];
+
+
+ //Get wire histos from root file and loop over each
+ // sense wire of a plane in shms Drift Chambers (DC1 or DC2)
+
+ f[ip]->cd(); //change to file containing uncorrected wire drift times
- //Declare number of entries in the tree
- Long64_t nentries = tree->GetEntries(); //number of triggers (particles that passed through all 4 hodo planes)
+ for (sw=1; sw<=tot_wires[ip]; sw++){
- //Loop over all entries
- for(Long64_t i=0; i<nentries; i++)
- {
- tree->GetEntry(i);
-
- //Loop over number of hits for each trigger in each DC plane
- for(ip=0; ip<NPLANES; ip++){
-
+
+ //set title of histos in root file
+ TString drift_time_histo = Form("wire_%d", sw);
+ //Get drift time histograms from root file
+ cell_dt[sw-1] = (TH1F*)f[ip]->Get(drift_time_histo);
+
+ //Create corrected wire histos
+ cell_dt_corr[sw-1] = new TH1F(plane_names[ip]+Form("_wire_%d_corr", sw), "", nbins, -50, 350);
+
+ shift = tzero[ip][sw-1]; //the shift represents how much the drift time histo needs to be offset
+
+ //cout << "sw: " << sw << " :: " << "offset: " << shift << endl;
- for(Int_t j=0; j<Ndata[ip]; j++){
-
- h[ip]->Fill(pdc_time[ip][j] - t_zero_offsets[ip]); //add t0 offset correction
- }
-
-
+ //************APPLY TZERO OFFSET ON A WIRE-BY-WIRE BASIS TO ALL WIRES IN ALL PLANES***************
-
-
- }
-}
+ //INCLUDE the code 'shift.C ', which shifts a histogram
+ for (i=1; i<=nbins; i++) {
+
+ bin_Content = cell_dt[sw-1]->GetBinContent(i);
+
+
+ cell_dt_corr[sw-1]->SetBinContent(i-shift/bin_width, bin_Content); //apply the t0 correction
+
+ }
-
-//Write histograms to file
-g->Write();
+ //*************************************************************************************************
+
+
+
+ //write wire drift times (after tzero corrections) to file
+ g->cd();
+ cell_dt_corr[sw-1]->Write(plane_names[ip]+Form("_wire_%d", sw), TObject::kWriteDelete);
+
+
+ //add all cell drift times into a single plane
+ h_add[ip]->Add(cell_dt_corr[sw-1]);
+
+
+
+
+ } // end loop over sense wires
+
+
+ //Wire combined wire drift times (t0 -corrected) for each plane to file,
+ g->cd();
+ h_add[ip]->Write();
+
+
+ } // end loop over planes
+
+
+
+
}
diff --git a/CALIBRATION/shms_dc_calib/scripts/get_tzero_param.C b/CALIBRATION/shms_dc_calib/scripts/get_tzero_param.C
index 6f8c06aa..1e1e061f 100644
--- a/CALIBRATION/shms_dc_calib/scripts/get_tzero_param.C
+++ b/CALIBRATION/shms_dc_calib/scripts/get_tzero_param.C
@@ -17,8 +17,15 @@ void get_tzero_param()
TString run = Form("run%d", run_NUM);
+ // 0 1 2 3 4 5 6 7 8 9 10 11
+ //TString planes[NPLANES] = {"1u1","1u2","1x1","1x2","1v1","1v2","2u1","2u2","2x1","2x2","2v1","2v2"};
+
+ TString planes[NPLANES] = {"1u1","1u2","1x1","1x2","1v1","1v2","2v2","2v1","2x2","2x1","2u2","2u1"}; //planes 0, 1, 2, 3, 4, 5 OK ::
+
+ //what needs to be changed
+
+ //2u1 (6) --> 2v2 (11)|| 2u2 (7) --> 2v1 (10)||, 2x1 (8) --> 2x2(9)||,
-TString planes[NPLANES] = {"1u1","1u2","1x1","1x2","1v1","1v2","2u1","2u2","2x1","2x2","2v1","2v2"};
int fNWires[NPLANES] = {107, 107, 79, 79, 107, 107, 107, 107, 79, 79, 107, 107};
int group_size[NPLANES] = {7, 7, 5, 5, 7, 7, 7, 7, 5, 5, 7, 7};
@@ -56,24 +63,24 @@ for (ip=0; ip<NPLANES; ip++){
dc_group_min[ip][0]=1; //set initial array value to 1st sensewire in the group
for (group=1; group < group_size[ip]; group++){
- if (ip==2 || ip==8){
- //dc_group_min[ip][5] = {1, 17, 33, 49, 65}; //plane 2(1x1) or plane 8(2x1)
+ if (ip==2 || ip==9){
+ //dc_group_min[ip][5] = {1, 17, 33, 49, 65}; //plane 2(1x1) or plane9(2x1)
dc_group_min[ip][group] = dc_group_min[ip][group-1] + 16;
}
- else if (ip==3 || ip==9){
- //dc_group_min[ip][5] = {{1, 16, 32, 48, 64}}; //plane 3(1x2) or plane 9(2x2)
+ else if (ip==3 || ip==8){
+ //dc_group_min[ip][5] = {{1, 16, 32, 48, 64}}; //plane 3(1x2) or plane 8(2x2)
dc_group_min[ip][1]=16;
dc_group_min[ip][group+1] = dc_group_min[ip][group] + 16;
}
- else if(ip==0 || ip==4 || ip==6 || ip ==10){
+ else if(ip==0 || ip==4 || ip==11 || ip ==7){
// dc_group_min[ip][7] = {{1, 16, 32, 48, 64, 80, 96}}; //planes: 1u1, 1v1, 2u1, 2v1
dc_group_min[ip][1]=16;
dc_group_min[ip][group+1] = dc_group_min[ip][group] + 16;
}
- else if(ip==1 || ip==5 || ip==7 || ip==11){
+ else if(ip==1 || ip==5 || ip==10 || ip==6){
// dc_group_min[ip][7] = {{1, 13, 29, 45, 61, 77, 93}}; //planes 1u2, 1v2, 2u2, 2v2
dc_group_min[ip][1]=13;
dc_group_min[ip][group+1] = dc_group_min[ip][group] + 16;
@@ -98,25 +105,25 @@ for (ip=0; ip<NPLANES; ip++){
for (group=1; group < group_size[ip]; group++){
- if (ip==2 || ip==8){
- //dc_group_max[ip][5] = {16, 32, 48, 64, 79}; //plane 2(1x1) or plane 8(2x1)
+ if (ip==2 || ip==9){
+ //dc_group_max[ip][5] = {16, 32, 48, 64, 79}; //plane 2(1x1) or plane 9(2x1)
dc_group_max[ip][0]=16;
dc_group_max[ip][group] = dc_group_max[ip][group-1] + 16;
}
- else if (ip==3 || ip==9){
- //dc_group_max[ip][5] = {15, 31, 47, 63, 79}; //plane 3(1x2) or plane 9(2x2)
+ else if (ip==3 || ip==8){
+ //dc_group_max[ip][5] = {15, 31, 47, 63, 79}; //plane 3(1x2) or plane 8(2x2)
dc_group_max[ip][0]=15;
dc_group_max[ip][group] = dc_group_max[ip][group-1] + 16;
}
- else if(ip==0 || ip==4 || ip==6 || ip ==10){
+ else if(ip==0 || ip==4 || ip==11 || ip ==7){
// dc_group_max[ip][7] = {15, 31, 47, 63, 79, 95, 107}; //planes: 1u1, 1v1, 2u1, 2v1
dc_group_max[ip][0]=15;
dc_group_max[ip][group] = dc_group_max[ip][group-1] + 16;
}
- else if(ip==1 || ip==5 || ip==7 || ip==11){
+ else if(ip==1 || ip==5 || ip==10 || ip==6){
// dc_group_max[ip][7] = {12, 28, 44, 60, 76, 92, 107}; //planes 1u2, 1v2, 2u2, 2v2
dc_group_max[ip][0]=12;
dc_group_max[ip][group] = dc_group_max[ip][group-1] + 16;
@@ -180,7 +187,7 @@ for (ip=0; ip<NPLANES; ip++){
if (ifs.is_open()) {
// cout << "File opened!" << endl;
- ifs.ignore(50000, '\n'); //ignore the first line
+ ifs.ignore(50000, '\n'); //ignore the first line
for (sw=0; sw<fNWires[ip]; ++sw) {
ifs >> wire[sw] >> t0[sw] >> t0_err[sw] >> entries[sw];
--
GitLab