Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
TPCPedestalCalibration.C
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file TPCPedestalCalibration.C
1 #include <cstdlib>
2 #include <iostream>
3 #include <map>
4 #include <string>
5 #include <vector>
6 
7 #include "TChain.h"
8 #include "TFile.h"
9 #include "TTree.h"
10 #include "TString.h"
11 #include "TObjString.h"
12 #include "TSystem.h"
13 #include "TROOT.h"
14 #include "TTreeReader.h"
15 #include "TTreeReaderValue.h"
16 #include "TTreeReaderArray.h"
17 
18 /*************************************************************/
19 /* TPC Pedestal Calibration */
20 /* Thomas Marshall,Aditya Dash */
21 /* rosstom@ucla.edu,aditya55@physics.ucla.edu */
22 /*************************************************************/
23 
24 void TPCPedestalCalibration(vector<string> inputFiles = {"/sphenix/user/rosstom/test/testFiles/TPC_ebdc00_pedestal-00010131-0000.prdf_TPCRawDataTree.root"}){
25  for (int fileNum = 0; fileNum < inputFiles.size(); fileNum++){
26  string sectorNum = inputFiles[fileNum];
27  size_t pos = sectorNum.find("ebdc");
28  sectorNum.erase(pos+6);
29  sectorNum.erase(sectorNum.begin(),sectorNum.begin()+pos+4);
30  Int_t sector;
31  if (sectorNum.substr(0,1) == "0")
32  {
33  sectorNum.erase(sectorNum.begin(),sectorNum.end()-1);
34  sector = std::stoi(sectorNum);
35  }
36  else sector = std::stoi(sectorNum);
37 
38  Int_t mod_arr[26]={2,2,1,1,1,3,3,3,3,3,3,2,2,1,2,2,1,1,2,2,3,3,3,3,3,3};
39  Int_t slot_arr[26] = {5,6,1,3,2,12,10,11,9,8,7,1,2,4,8,7,6,5,4,3,1,3,2,4,6,5};
40 
41  TFile *pedestalInput = TFile::Open(inputFiles[fileNum].c_str());
42  TTreeReader myReader("SampleTree", pedestalInput);
43 
44  TTreeReaderValue<Int_t> nSamples_in(myReader, "nSamples");
45  TTreeReaderValue<Int_t> fee_in(myReader, "fee");
46  TTreeReaderArray<UShort_t> adcSamples_in(myReader, "adcSamples");
47  TTreeReaderValue<Int_t> Channel_in(myReader, "Channel");
48 
49 //Adding the output tree
50  TString* outputfilename=new TString("pedestalCalibration_TPC_ebdc"+sectorNum+".root");
51  TFile* outputfile=new TFile(outputfilename->Data(),"recreate");
52  TTree* outputTree=new TTree("outputTree","outputTree");
53  Int_t isAlive=1; // 1 if the channel is working properly, 0 if no signal(number of samples is 0, adc value is 0 or nan), pedestal above 200 or below 10
54  Float_t pedMean; // average pedestal value over all samples
55  Float_t pedStd; // pedestal standard deviation over all samples
56  Int_t chan; // channel number
57  Int_t fee; // fee number
58  Int_t module; // module number (e.g. R1, R2, R3)
59  Int_t slot; // slot number
60  outputTree->Branch("isAlive",&isAlive,"isAlive/I");
61  outputTree->Branch("pedMean",&pedMean,"pedMean/F");
62  outputTree->Branch("pedStd",&pedStd,"pedStd/F");
63  outputTree->Branch("sector",&sector,"sector/I");
64  outputTree->Branch("fee",&fee,"fee/I");
65  outputTree->Branch("chan",&chan,"chan/I");
66  outputTree->Branch("module",&module,"module/I");
67  outputTree->Branch("slot",&slot,"slot/I");
68 
69  Float_t ave_adc_fee_channel[26][256];
70  Float_t std_adc_fee_channel[26][256];
71  Float_t counts_adc_fee_channel[26][256];
72  Int_t alive_array_fee_channel[26][256];
73  for(int fee_no=0;fee_no<26;fee_no++){
74  for(int channel_no=0;channel_no<256;channel_no++)
75  {
76  ave_adc_fee_channel[fee_no][channel_no]=0.0;
77  std_adc_fee_channel[fee_no][channel_no]=0.0;
78  counts_adc_fee_channel[fee_no][channel_no]=0.0;
79  alive_array_fee_channel[fee_no][channel_no]=1;
80  }
81  }
82 
83  while(myReader.Next())
84  {
85  if(*nSamples_in==0)
86  {
87  alive_array_fee_channel[*fee_in][*Channel_in]=0;
88  continue;
89  }
90 
91  bool dead = false;
92  for(int adc_sam_no=0;adc_sam_no<*nSamples_in;adc_sam_no++)
93  {
94  if (adcSamples_in[adc_sam_no] == 0 || TMath::IsNaN(float(adcSamples_in[adc_sam_no])))
95  {
96  alive_array_fee_channel[*fee_in][*Channel_in]=0;
97  }
98  }
99  if (dead) {continue;}
100 
101  for(int adc_sam_no=0;adc_sam_no<*nSamples_in;adc_sam_no++){
102  ave_adc_fee_channel[*fee_in][*Channel_in]+=adcSamples_in[adc_sam_no];
103  std_adc_fee_channel[*fee_in][*Channel_in]+=pow(adcSamples_in[adc_sam_no],2);
104  counts_adc_fee_channel[*fee_in][*Channel_in]+=1.0;
105  }
106 
107  }
108 
109  for(int fee_no=0;fee_no<26;fee_no++)
110  {
111  for(int channel_no=0;channel_no<256;channel_no++)
112  {
113  if(counts_adc_fee_channel[fee_no][channel_no] != 0.0)
114  {
115  float temp1 = ave_adc_fee_channel[fee_no][channel_no]/counts_adc_fee_channel[fee_no][channel_no];
116  float temp2 = std_adc_fee_channel[fee_no][channel_no]/counts_adc_fee_channel[fee_no][channel_no];
117  ave_adc_fee_channel[fee_no][channel_no] = temp1;
118  std_adc_fee_channel[fee_no][channel_no] = temp2;
119  }
120  else
121  {
122  ave_adc_fee_channel[fee_no][channel_no] = 0.0;
123  std_adc_fee_channel[fee_no][channel_no] = 0.0;
124  alive_array_fee_channel[fee_no][channel_no]=0;
125  }
126 
127  if(ave_adc_fee_channel[fee_no][channel_no] > 200 || ave_adc_fee_channel[fee_no][channel_no] < 10)
128  {
129  alive_array_fee_channel[fee_no][channel_no]=0;
130  }
131 
132  pedMean=ave_adc_fee_channel[fee_no][channel_no];
133  pedStd=sqrt(std_adc_fee_channel[fee_no][channel_no] - pow(ave_adc_fee_channel[fee_no][channel_no],2));
134  isAlive=alive_array_fee_channel[fee_no][channel_no];
135  chan=channel_no;
136  fee=fee_no;
137  module=mod_arr[fee_no];
138  slot=slot_arr[fee_no];
139  outputTree->Fill();
140  }
141  }
142  outputfile->cd();
143  outputTree->Write();
144  outputfile->Close();
145  pedestalInput->Close();
146  }
147 }